spheresld.cpp @ 3a084629

ESS_GUIESS_GUI_DocsESS_GUI_batch_fittingESS_GUI_bumps_abstractionESS_GUI_iss1116ESS_GUI_iss879ESS_GUI_iss959ESS_GUI_openclESS_GUI_orderingESS_GUI_sync_sascalccostrafo411magnetic_scattrelease-4.1.1release-4.1.2release-4.2.2release_4.0.1ticket-1009ticket-1094-headlessticket-1242-2d-resolutionticket-1243ticket-1249ticket885unittest-saveload

Last change on this file since 3a084629 was 7a8faf8, checked in by Jae Cho <jhjcho@…>, 14 years ago
Added polydispersion in a couple of params in spheresldmodel
Property mode set to `100644`
File size: 9.5 KB

Rev	Line
[0164899a]	1
	2	#include <math.h>
	3	#include "models.hh"
	4	#include "parameters.hh"
	5	#include <stdio.h>
	6	using namespace std;
	7
	8	extern "C" {
	9	#include "spheresld.h"
	10	}
	11
	12	SphereSLDModel :: SphereSLDModel() {
	13	n_shells = Parameter(1.0);
	14	scale = Parameter(1.0);
[7a8faf8]	15	thick_inter0 = Parameter(1.0, true);
	16	thick_inter0.set_min(0.0);
[0164899a]	17	func_inter0 = Parameter(0);
	18	sld_core0 = Parameter(2.07e-06);
	19	sld_solv = Parameter(1.0e-06);
	20	background = Parameter(0.0);
	21
	22
	23	sld_flat1 = Parameter(2.7e-06);
	24	sld_flat2 = Parameter(3.5e-06);
	25	sld_flat3 = Parameter(4.0e-06);
	26	sld_flat4 = Parameter(3.5e-06);
	27	sld_flat5 = Parameter(4.0e-06);
	28	sld_flat6 = Parameter(3.5e-06);
	29	sld_flat7 = Parameter(4.0e-06);
	30	sld_flat8 = Parameter(3.5e-06);
	31	sld_flat9 = Parameter(4.0e-06);
	32	sld_flat10 = Parameter(3.5e-06);
	33
	34
	35	thick_inter1 = Parameter(1.0);
	36	thick_inter2 = Parameter(1.0);
	37	thick_inter3 = Parameter(1.0);
	38	thick_inter4 = Parameter(1.0);
	39	thick_inter5 = Parameter(1.0);
	40	thick_inter6 = Parameter(1.0);
	41	thick_inter7 = Parameter(1.0);
	42	thick_inter8 = Parameter(1.0);
	43	thick_inter9 = Parameter(1.0);
	44	thick_inter10 = Parameter(1.0);
	45
	46
	47	thick_flat1 = Parameter(100.0);
	48	thick_flat2 = Parameter(100.0);
	49	thick_flat3 = Parameter(100.0);
	50	thick_flat4 = Parameter(100.0);
	51	thick_flat5 = Parameter(100.0);
	52	thick_flat6 = Parameter(100.0);
	53	thick_flat7 = Parameter(100.0);
	54	thick_flat8 = Parameter(100.0);
	55	thick_flat9 = Parameter(100.0);
	56	thick_flat10 = Parameter(100.0);
	57
	58
	59	func_inter1 = Parameter(0);
	60	func_inter2 = Parameter(0);
	61	func_inter3 = Parameter(0);
	62	func_inter4 = Parameter(0);
	63	func_inter5 = Parameter(0);
	64	func_inter6 = Parameter(0);
	65	func_inter7 = Parameter(0);
	66	func_inter8 = Parameter(0);
	67	func_inter9 = Parameter(0);
	68	func_inter10 = Parameter(0);
	69
	70	nu_inter1 = Parameter(2.5);
	71	nu_inter2 = Parameter(2.5);
	72	nu_inter3 = Parameter(2.5);
	73	nu_inter4 = Parameter(2.5);
	74	nu_inter5 = Parameter(2.5);
	75	nu_inter6 = Parameter(2.5);
	76	nu_inter7 = Parameter(2.5);
	77	nu_inter8 = Parameter(2.5);
	78	nu_inter9 = Parameter(2.5);
	79	nu_inter10 = Parameter(2.5);
	80
	81	npts_inter = Parameter(35.0);
	82	nu_inter0 = Parameter(2.5);
[7a8faf8]	83	rad_core0 = Parameter(60.0, true);
	84	rad_core0.set_min(0.0);
[0164899a]	85	}
	86
	87	/**
	88	* Function to evaluate 1D SphereSLD function
	89	* @param q: q-value
	90	* @return: function value
	91	*/
	92	double SphereSLDModel :: operator()(double q) {
	93	double dp[60];
	94	// Fill parameter array for IGOR library
	95	// Add the background after averaging
	96	dp[0] = n_shells();
	97	dp[1] = scale();
	98	dp[2] = thick_inter0();
	99	dp[3] = func_inter0();
	100	dp[4] = sld_core0();
	101	dp[5] = sld_solv();
[7a8faf8]	102	dp[6] = 0.0;
[0164899a]	103
	104	dp[7] = sld_flat1();
	105	dp[8] = sld_flat2();
	106	dp[9] = sld_flat3();
	107	dp[10] = sld_flat4();
	108	dp[11] = sld_flat5();
	109	dp[12] = sld_flat6();
	110	dp[13] = sld_flat7();
	111	dp[14] = sld_flat8();
	112	dp[15] = sld_flat9();
	113	dp[16] = sld_flat10();
	114
	115	dp[17] = thick_inter1();
	116	dp[18] = thick_inter2();
	117	dp[19] = thick_inter3();
	118	dp[20] = thick_inter4();
	119	dp[21] = thick_inter5();
	120	dp[22] = thick_inter6();
	121	dp[23] = thick_inter7();
	122	dp[24] = thick_inter8();
	123	dp[25] = thick_inter9();
	124	dp[26] = thick_inter10();
	125
	126	dp[27] = thick_flat1();
	127	dp[28] = thick_flat2();
	128	dp[29] = thick_flat3();
	129	dp[30] = thick_flat4();
	130	dp[31] = thick_flat5();
	131	dp[32] = thick_flat6();
	132	dp[33] = thick_flat7();
	133	dp[34] = thick_flat8();
	134	dp[35] = thick_flat9();
	135	dp[36] = thick_flat10();
	136
	137	dp[37] = func_inter1();
	138	dp[38] = func_inter2();
	139	dp[39] = func_inter3();
	140	dp[40] = func_inter4();
	141	dp[41] = func_inter5();
	142	dp[42] = func_inter6();
	143	dp[43] = func_inter7();
	144	dp[44] = func_inter8();
	145	dp[45] = func_inter9();
	146	dp[46] = func_inter10();
	147
	148	dp[47] = nu_inter1();
	149	dp[48] = nu_inter2();
	150	dp[49] = nu_inter3();
	151	dp[50] = nu_inter4();
	152	dp[51] = nu_inter5();
	153	dp[52] = nu_inter6();
	154	dp[53] = nu_inter7();
	155	dp[54] = nu_inter8();
	156	dp[55] = nu_inter9();
	157	dp[56] = nu_inter10();
	158
	159
	160	dp[57] = npts_inter();
	161	dp[58] = nu_inter0();
	162	dp[59] = rad_core0();
	163
	164	// No polydispersion supported in this model.
[7a8faf8]	165	// Get the dispersion points for the radius
	166	vector<WeightPoint> weights_rad_core0;
	167	rad_core0.get_weights(weights_rad_core0);
	168	vector<WeightPoint> weights_thick_inter0;
	169	thick_inter0.get_weights(weights_thick_inter0);
	170	// Perform the computation, with all weight points
	171	double sum = 0.0;
	172	double norm = 0.0;
	173	double vol = 0.0;
	174
	175	// Loop over core weight points
	176	for(int i=0; i<weights_rad_core0.size(); i++) {
	177	dp[59] = weights_rad_core0[i].value;
	178	// Loop over thick_inter0 weight points
	179	for(int j=0; j<weights_thick_inter0.size(); j++) {
	180	dp[2] = weights_thick_inter0[j].value;
	181
	182	//Un-normalize Sphere by volume
	183	sum += weights_rad_core0[i].weight * weights_thick_inter0[j].weight
	184	* sphere_sld_kernel(dp,q) * pow((weights_rad_core0[i].value +
	185	weights_thick_inter0[j].value),3.0);
	186	//Find average volume
	187	vol += weights_rad_core0[i].weight * weights_thick_inter0[j].weight
	188	* pow((weights_rad_core0[i].value+weights_thick_inter0[j].value),3.0);
	189
	190	norm += weights_rad_core0[i].weight * weights_thick_inter0[j].weight;
	191	}
	192	}
	193
	194	if (vol != 0.0 && norm != 0.0) {
	195	//Re-normalize by avg volume
	196	sum = sum/(vol/norm);}
	197
	198	return sum/norm + background();
[0164899a]	199	}
	200
	201	/**
	202	* Function to evaluate 2D SphereSLD function
	203	* @param q_x: value of Q along x
	204	* @param q_y: value of Q along y
	205	* @return: function value
	206	*/
	207	double SphereSLDModel :: operator()(double qx, double qy) {
	208	double q = sqrt(qxqx + qyqy);
	209	return (*this).operator()(q);
	210	}
	211
	212	/**
	213	* Function to evaluate SphereSLD function
	214	* @param pars: parameters of the SphereSLD
	215	* @param q: q-value
	216	* @param phi: angle phi
	217	* @return: function value
	218	*/
	219	double SphereSLDModel :: evaluate_rphi(double q, double phi) {
	220	return (*this).operator()(q);
	221	}
	222
	223	/**
	224	* Function to calculate TOTAL radius
	225	* ToDo: Find What is the effective radius for this model.
	226	* @return: effective radius value
	227	*/
	228	// No polydispersion supported in this model.
	229	// Calculate max radius assumming max_radius = effective radius
	230	// Note that this max radius is not affected by sld of layer, sld of interface, or
	231	// sld of solvent.
	232	double SphereSLDModel :: calculate_ER() {
	233	SphereSLDParameters dp;
	234
	235	dp.n_shells = n_shells();
	236
	237	dp.rad_core0 = rad_core0();
	238	dp.thick_flat1 = thick_flat1();
	239	dp.thick_flat2 = thick_flat2();
	240	dp.thick_flat3 = thick_flat3();
	241	dp.thick_flat4 = thick_flat4();
	242	dp.thick_flat5 = thick_flat5();
	243	dp.thick_flat6 = thick_flat6();
	244	dp.thick_flat7 = thick_flat7();
	245	dp.thick_flat8 = thick_flat8();
	246	dp.thick_flat9 = thick_flat9();
	247	dp.thick_flat10 = thick_flat10();
	248
	249	dp.thick_inter0 = thick_inter0();
	250	dp.thick_inter1 = thick_inter1();
	251	dp.thick_inter2 = thick_inter2();
	252	dp.thick_inter3 = thick_inter3();
	253	dp.thick_inter4 = thick_inter4();
	254	dp.thick_inter5 = thick_inter5();
	255	dp.thick_inter6 = thick_inter6();
	256	dp.thick_inter7 = thick_inter7();
	257	dp.thick_inter8 = thick_inter8();
	258	dp.thick_inter9 = thick_inter9();
	259	dp.thick_inter10 = thick_inter10();
	260
[7a8faf8]	261	double rad_out = 0.0;
	262	double out = 0.0;
	263	// Perform the computation, with all weight points
	264	double sum = 0.0;
	265	double norm = 0.0;
	266
	267	// Get the dispersion points for the radius
	268	vector<WeightPoint> weights_rad_core0;
	269	rad_core0.get_weights(weights_rad_core0);
	270
	271	// Get the dispersion points for the thick 1
	272	vector<WeightPoint> weights_thick_inter0;
	273	thick_inter0.get_weights(weights_thick_inter0);
	274	// Loop over radius weight points
	275	for(int i=0; i<weights_rad_core0.size(); i++) {
	276	dp.rad_core0 = weights_rad_core0[i].value;
	277	// Loop over radius weight points
	278	for(int j=0; j<weights_thick_inter0.size(); j++) {
	279	dp.thick_inter0 = weights_thick_inter0[j].value;
	280	rad_out = dp.rad_core0 + dp.thick_inter0;
	281	if (dp.n_shells > 0)
	282	rad_out += dp.thick_flat1 + dp.thick_inter1;
	283	if (dp.n_shells > 1)
	284	rad_out += dp.thick_flat2 + dp.thick_inter2;
	285	if (dp.n_shells > 2)
	286	rad_out += dp.thick_flat3 + dp.thick_inter3;
	287	if (dp.n_shells > 3)
	288	rad_out += dp.thick_flat4 + dp.thick_inter4;
	289	if (dp.n_shells > 4)
	290	rad_out += dp.thick_flat5 + dp.thick_inter5;
	291	if (dp.n_shells > 5)
	292	rad_out += dp.thick_flat6 + dp.thick_inter6;
	293	if (dp.n_shells > 6)
	294	rad_out += dp.thick_flat7 + dp.thick_inter7;
	295	if (dp.n_shells > 7)
	296	rad_out += dp.thick_flat8 + dp.thick_inter8;
	297	if (dp.n_shells > 8)
	298	rad_out += dp.thick_flat9 + dp.thick_inter9;
	299	if (dp.n_shells > 9)
	300	rad_out += dp.thick_flat10 + dp.thick_inter10;
	301	sum += weights_rad_core0[i].weight*weights_thick_inter0[j].weight
	302	* (rad_out);
	303	norm += weights_rad_core0[i].weight*weights_thick_inter0[j].weight;
	304	}
	305	}
	306	if (norm != 0){
	307	//return the averaged value
	308	out = sum/norm;}
	309	else{
	310	//return normal value
	311	out = dp.rad_core0 + dp.thick_inter0;
	312	if (dp.n_shells > 0)
	313	out += dp.thick_flat1 + dp.thick_inter1;
	314	if (dp.n_shells > 1)
	315	out += dp.thick_flat2 + dp.thick_inter2;
	316	if (dp.n_shells > 2)
	317	out += dp.thick_flat3 + dp.thick_inter3;
	318	if (dp.n_shells > 3)
	319	out += dp.thick_flat4 + dp.thick_inter4;
	320	if (dp.n_shells > 4)
	321	out += dp.thick_flat5 + dp.thick_inter5;
	322	if (dp.n_shells > 5)
	323	out += dp.thick_flat6 + dp.thick_inter6;
	324	if (dp.n_shells > 6)
	325	out += dp.thick_flat7 + dp.thick_inter7;
	326	if (dp.n_shells > 7)
	327	out += dp.thick_flat8 + dp.thick_inter8;
	328	if (dp.n_shells > 8)
	329	out += dp.thick_flat9 + dp.thick_inter9;
	330	if (dp.n_shells > 9)
	331	out += dp.thick_flat10 + dp.thick_inter10;
	332	}
	333
	334	return out;
[0164899a]	335
	336	}

Note: See TracBrowser for help on using the repository browser.

SasView

source: sasview/sansmodels/src/sans/models/c_models/spheresld.cpp @ 3a084629

Download in other formats: