spheresld.cpp @ 19aa896

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 19aa896 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

Line
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);
15	thick_inter0 = Parameter(1.0, true);
16	thick_inter0.set_min(0.0);
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);
83	rad_core0 = Parameter(60.0, true);
84	rad_core0.set_min(0.0);
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();
102	dp[6] = 0.0;
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.
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();
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
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;
335
336	}

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

SasView

source: sasview/sansmodels/src/sans/models/c_models/spheresld.cpp @ 19aa896

Download in other formats: