Context Navigation

Func2D.c @ 756

Last change on this file since 756 was 235, checked in by srkline, 15 years ago

Changes to XOP library functions to convert models to take individual SLD's rather than contrast.

Individual SLD's are easier to work with - since they are experimental values, contrast is not.

This, unfortunately, means that what Mathieu has done may need to be tweaked, since our library has been tweaked.

File size: 14.3 KB

Line
1	/*
2	* Func2D.c
3	* SANSAnalysis
4	*
5	* Steve Kline Jan 2008
6	* these are the XOP calls to DANSE 2D functions
7	* from the library supplied by M. Doucet
8	*
9	*/
10
11
12	#include "XOPStandardHeaders.h" // Include ANSI headers, Mac headers, IgorXOP.h, XOP.h and XOPSupport.h
13	#include "SANSAnalysis.h" // structures, two-byte aligned
14	#include "libSANSAnalysis.h" // functions from the libSANS
15	#include "danse.h" // functions from the DANSE part of libSANS
16	#include <math.h>
17	#include "Func2D.h" // declarations for the 2D functions in this file
18
19	int
20	Cylinder_2D(FitParams2DPtr p)
21	{
22	double *dp;
23
24	float *fp; // Pointer to single precision wave data.
25	double qx;
26	double qy;
27	double q, phi;
28	double pars[11];
29	// int i;
30	// char buf[256];
31
32	if (p->waveHandle == NIL) {
33	SetNaN64(&p->result);
34	return NON_EXISTENT_WAVE;
35	}
36
37	qx = p->qx;
38	qy = p->qy;
39
40	// sprintf(buf, "Qx = %g, Qy = %g\r",qx, qy);
41	// XOPNotice(buf);
42
43	q = hypot(qx,qy);
44	phi = atan2(qy,qx);
45
46	switch(WaveType(p->waveHandle)){ // We can handle single and double precision coefficient waves.
47	case NT_FP32:
48	fp= WaveData(p->waveHandle);
49	SetNaN64(&p->result);
50	return REQUIRES_SP_OR_DP_WAVE; //not quite true, but good enough for now AJJ 4/23/07
51	case NT_FP64:
52	dp= WaveData(p->waveHandle);
53
54	// for(i=0; i<11; i++) {
55	// pars[i] = dp[i];
56	// }
57	pars[0] = dp[0];
58	pars[1] = dp[1];
59	pars[2] = dp[2];
60	pars[3] = dp[3] - dp[4];
61	pars[4] = dp[5];
62	pars[5] = dp[6];
63	pars[6] = dp[7];
64	pars[7] = dp[8];
65	pars[8] = dp[9];
66	pars[9] = dp[10];
67	pars[10] = dp[11];
68
69	p->result = disperse_cylinder_analytical_2D( pars, q, phi );
70	return 0;
71	default: // We can't handle this wave data type.
72	SetNaN64(&p->result);
73	return REQUIRES_SP_OR_DP_WAVE;
74	}
75
76	return 0;
77	}
78
79	int
80	Cylinder_2D_Weight2D(FitParams2DWeightPtr p)
81	{
82	double *dp;
83
84	float *fp; // Pointer to single precision wave data.
85	double qx;
86	double qy;
87	double q, phi;
88	double *par_values;
89	double *weight_values;
90	double pars[13];
91	int i, i_theta;
92	double sum;
93	int n_slices;
94
95	if (p->waveHandle == NIL) {
96	SetNaN64(&p->result);
97	return NON_EXISTENT_WAVE;
98	}
99
100	qx = p->qx;
101	qy = p->qy;
102
103	q = hypot(qx,qy);
104	phi = atan2(qy,qx);
105
106	switch(WaveType(p->waveHandle)){ // We can handle single and double precision coefficient waves.
107	case NT_FP32:
108	fp= WaveData(p->waveHandle);
109	SetNaN64(&p->result);
110	return REQUIRES_SP_OR_DP_WAVE; //not quite true, but good enough for now AJJ 4/23/07
111	case NT_FP64:
112	dp= WaveData(p->waveHandle);
113	par_values = WaveData(p->par_values);
114	weight_values = WaveData(p->weight_values);
115
116	for(i=0; i<13; i++) {
117	pars[i] = dp[i];
118	}
119
120	sum = 0.0;
121	n_slices = (int)floor(dp[13]);
122
123	//XOPOKAlert("test","This is a test");
124
125	if(n_slices == 0) {
126	p->result =weight_dispersion( &disperse_cylinder_analytical_2D,
127	par_values, weight_values,
128	(int)floor(pars[11]), (int)floor(pars[12]),
129	pars, q, phi );
130	} else {
131	for(i_theta=0; i_theta<n_slices; i_theta++) {
132	SpinProcess();
133	// average over theta
134
135	// For the cylinder model, theta_cyl=90 degrees
136	// will produce a NAN at phi_cyl=0 and 45 degrees
137	pars[5] = acos(-1.0)/n_slices * i_theta;
138	if( fabs(i_theta / n_slices) - 0.5 < 0.000001 ) {
139	//continue;
140	pars[5] += 0.00001;
141	}
142
143	// Multiply by sin(theta) because we are integrating in
144	// spherical coordinates
145	sum += sin(pars[5])* weight_dispersion( &disperse_cylinder_analytical_2D,
146	par_values, weight_values,
147	(int)floor(pars[11]), (int)floor(pars[12]),
148	pars, q, phi );
149	}
150
151	p->result = sum/n_slices;
152	}
153
154	return 0;
155	default: // We can't handle this wave data type.
156	SetNaN64(&p->result);
157	return REQUIRES_SP_OR_DP_WAVE;
158	}
159
160	return 0;
161	}
162
163	int
164	CoreShellCylinder_2D(FitParams2DPtr p)
165	{
166	double *dp;
167
168	float *fp; // Pointer to single precision wave data.
169	double qx;
170	double qy;
171	double q, phi;
172	double pars[15];
173	int i;
174
175	if (p->waveHandle == NIL) {
176	SetNaN64(&p->result);
177	return NON_EXISTENT_WAVE;
178	}
179
180	qx = p->qx;
181	qy = p->qy;
182
183	q = hypot(qx,qy);
184	phi = atan2(qy,qx);
185
186	switch(WaveType(p->waveHandle)){ // We can handle single and double precision coefficient waves.
187	case NT_FP32:
188	fp= WaveData(p->waveHandle);
189	SetNaN64(&p->result);
190	return REQUIRES_SP_OR_DP_WAVE; //not quite true, but good enough for now AJJ 4/23/07
191	case NT_FP64:
192	dp= WaveData(p->waveHandle);
193
194	for(i=0; i<15; i++) {
195	pars[i] = dp[i];
196	}
197
198	p->result = disperse_core_shell_cylinder_analytical_2D( pars, q, phi );
199	return 0;
200	default: // We can't handle this wave data type.
201	SetNaN64(&p->result);
202	return REQUIRES_SP_OR_DP_WAVE;
203	}
204
205	return 0;
206	}
207
208	int
209	CoreShellCylinder_2D_Weight2D(FitParams2DWeightPtr p)
210	{
211	double *dp;
212
213	float *fp; // Pointer to single precision wave data.
214	double qx;
215	double qy;
216	double q, phi;
217	double *par_values;
218	double *weight_values;
219	double pars[17];
220	int i, i_theta;
221	double sum;
222	int n_slices;
223
224	if (p->waveHandle == NIL) {
225	SetNaN64(&p->result);
226	return NON_EXISTENT_WAVE;
227	}
228
229	qx = p->qx;
230	qy = p->qy;
231
232	q = hypot(qx,qy);
233	phi = atan2(qy,qx);
234
235	switch(WaveType(p->waveHandle)){ // We can handle single and double precision coefficient waves.
236	case NT_FP32:
237	fp= WaveData(p->waveHandle);
238	SetNaN64(&p->result);
239	return REQUIRES_SP_OR_DP_WAVE; //not quite true, but good enough for now AJJ 4/23/07
240	case NT_FP64:
241	dp= WaveData(p->waveHandle);
242	par_values = WaveData(p->par_values);
243	weight_values = WaveData(p->weight_values);
244
245	for(i=0; i<17; i++) {
246	pars[i] = dp[i];
247	}
248
249	sum = 0.0;
250	n_slices = (int)floor(dp[17]);
251
252	if(n_slices == 0) {
253	p->result =weight_dispersion( &disperse_core_shell_cylinder_analytical_2D,
254	par_values, weight_values,
255	(int)floor(pars[15]), (int)floor(pars[16]),
256	pars, q, phi );
257	} else {
258	for(i_theta=0; i_theta<n_slices; i_theta++) {
259	SpinProcess();
260	// average over theta
261
262	// For the cylinder model, theta_cyl=90 degrees
263	// will produce a NAN at phi_cyl=0 and 45 degrees
264	// TODO: integrate from 0 to pi/2 instead of 0 to pi
265	pars[8] = acos(-1.0)/n_slices * i_theta;
266	if( fabs(i_theta / n_slices) - 0.5 < 0.000001 ) {
267	//continue;
268	pars[8] += 0.00001;
269	}
270
271	// Multiply by sin(theta) because we are integrating in
272	// spherical coordinates
273	sum += sin(pars[8])* weight_dispersion( &disperse_core_shell_cylinder_analytical_2D,
274	par_values, weight_values,
275	(int)floor(pars[15]), (int)floor(pars[16]),
276	pars, q, phi );
277	}
278
279	p->result = sum/n_slices;
280	}
281
282	return 0;
283	default: // We can't handle this wave data type.
284	SetNaN64(&p->result);
285	return REQUIRES_SP_OR_DP_WAVE;
286	}
287
288	return 0;
289	}
290
291	int
292	Ellipsoid_2D(FitParams2DPtr p)
293	{
294	double *dp;
295
296	float *fp; // Pointer to single precision wave data.
297	double qx;
298	double qy;
299	double q, phi;
300	double pars[12];
301	// int i;
302
303	if (p->waveHandle == NIL) {
304	SetNaN64(&p->result);
305	return NON_EXISTENT_WAVE;
306	}
307
308	qx = p->qx;
309	qy = p->qy;
310
311	q = hypot(qx,qy);
312	phi = atan2(qy,qx);
313
314	switch(WaveType(p->waveHandle)){ // We can handle single and double precision coefficient waves.
315	case NT_FP32:
316	fp= WaveData(p->waveHandle);
317	SetNaN64(&p->result);
318	return REQUIRES_SP_OR_DP_WAVE; //not quite true, but good enough for now AJJ 4/23/07
319	case NT_FP64:
320	dp= WaveData(p->waveHandle);
321
322	// for(i=0; i<12; i++) {
323	// pars[i] = dp[i];
324	// }
325	pars[0] = dp[0];
326	pars[1] = dp[1];
327	pars[2] = dp[2];
328	pars[3] = dp[3] - dp[4];
329	pars[4] = dp[5];
330	pars[5] = dp[6];
331	pars[6] = dp[7];
332	pars[7] = dp[8];
333	pars[8] = dp[9];
334	pars[9] = dp[10];
335	pars[10] = dp[11];
336	pars[11] = dp[12];
337
338	//p->result = 1.0;
339	p->result = disperse_ellipsoid_analytical_2D( pars, q, phi );
340	return 0;
341	default: // We can't handle this wave data type.
342	SetNaN64(&p->result);
343	return REQUIRES_SP_OR_DP_WAVE;
344	}
345
346	return 0;
347	}
348
349	int
350	Ellipsoid_2D_Weight2D(FitParams2DWeightPtr p)
351	{
352	double *dp;
353
354	float *fp; // Pointer to single precision wave data.
355	double qx;
356	double qy;
357	double q, phi;
358	double *par_values;
359	double *weight_values;
360	double pars[14];
361	int i, i_theta;
362	double sum;
363	int n_slices;
364
365	if (p->waveHandle == NIL) {
366	SetNaN64(&p->result);
367	return NON_EXISTENT_WAVE;
368	}
369
370	qx = p->qx;
371	qy = p->qy;
372
373	q = hypot(qx,qy);
374	phi = atan2(qy,qx);
375
376	switch(WaveType(p->waveHandle)){ // We can handle single and double precision coefficient waves.
377	case NT_FP32:
378	fp= WaveData(p->waveHandle);
379	SetNaN64(&p->result);
380	return REQUIRES_SP_OR_DP_WAVE; //not quite true, but good enough for now AJJ 4/23/07
381	case NT_FP64:
382	dp= WaveData(p->waveHandle);
383	par_values = WaveData(p->par_values);
384	weight_values = WaveData(p->weight_values);
385
386	for(i=0; i<14; i++) {
387	pars[i] = dp[i];
388	}
389
390	sum = 0.0;
391	n_slices = (int)floor(dp[14]);
392
393	if(n_slices == 0) {
394	p->result =weight_dispersion( &disperse_ellipsoid_analytical_2D,
395	par_values, weight_values,
396	(int)floor(pars[12]), (int)floor(pars[13]),
397	pars, q, phi );
398	} else {
399	for(i_theta=0; i_theta<n_slices; i_theta++) {
400	SpinProcess();
401	// average over theta
402
403	// For the cylinder model, theta_cyl=90 degrees
404	// will produce a NAN at phi_cyl=0 and 45 degrees
405	pars[5] = acos(-1.0)/n_slices * i_theta;
406	if( fabs(i_theta / n_slices) - 0.5 < 0.000001 ) {
407	//continue;
408	pars[5] += 0.00001;
409	}
410
411	// Multiply by sin(theta) because we are integrating in
412	// spherical coordinates
413	sum += sin(pars[5])* weight_dispersion( &disperse_ellipsoid_analytical_2D,
414	par_values, weight_values,
415	(int)floor(pars[12]), (int)floor(pars[13]),
416	pars, q, phi );
417	}
418
419	p->result = sum/n_slices;
420	}
421
422	return 0;
423	default: // We can't handle this wave data type.
424	SetNaN64(&p->result);
425	return REQUIRES_SP_OR_DP_WAVE;
426	}
427
428	return 0;
429	}
430
431
432	int
433	EllipticalCylinder_2D(FitParams2DPtr p)
434	{
435	double *dp;
436
437	float *fp; // Pointer to single precision wave data.
438	double qx;
439	double qy;
440	double q, phi;
441	double pars[14];
442	// int i;
443
444	if (p->waveHandle == NIL) {
445	SetNaN64(&p->result);
446	return NON_EXISTENT_WAVE;
447	}
448
449	qx = p->qx;
450	qy = p->qy;
451
452	q = hypot(qx,qy);
453	phi = atan2(qy,qx);
454
455	switch(WaveType(p->waveHandle)){ // We can handle single and double precision coefficient waves.
456	case NT_FP32:
457	fp= WaveData(p->waveHandle);
458	SetNaN64(&p->result);
459	return REQUIRES_SP_OR_DP_WAVE; //not quite true, but good enough for now AJJ 4/23/07
460	case NT_FP64:
461	dp= WaveData(p->waveHandle);
462
463	// for(i=0; i<14; i++) {
464	// pars[i] = dp[i];
465	// }
466	pars[0] = dp[0];
467	pars[1] = dp[1];
468	pars[2] = dp[2];
469	pars[3] = dp[3];
470	pars[4] = dp[4] - dp[5];
471	pars[5] = dp[6];
472	pars[6] = dp[7];
473	pars[7] = dp[8];
474	pars[8] = dp[9];
475	pars[9] = dp[10];
476	pars[10] = dp[11];
477	pars[11] = dp[12];
478	pars[12] = dp[13];
479	pars[13] = dp[14];
480
481
482	p->result = disperse_elliptical_cylinder_analytical_2D( pars, q, phi );
483
484	return 0;
485	default: // We can't handle this wave data type.
486	SetNaN64(&p->result);
487	return REQUIRES_SP_OR_DP_WAVE;
488	}
489
490	return 0;
491	}
492
493	int
494	EllipticalCylinder_2D_Weight2D(FitParams2DWeightPtr p)
495	{
496	double *dp;
497
498	float *fp; // Pointer to single precision wave data.
499	double qx;
500	double qy;
501	double q, phi;
502	double *par_values;
503	double *weight_values;
504	double pars[16];
505	int i, i_theta;
506	double sum;
507	int n_slices;
508
509	if (p->waveHandle == NIL) {
510	SetNaN64(&p->result);
511	return NON_EXISTENT_WAVE;
512	}
513
514	qx = p->qx;
515	qy = p->qy;
516
517	q = hypot(qx,qy);
518	phi = atan2(qy,qx);
519
520	switch(WaveType(p->waveHandle)){ // We can handle single and double precision coefficient waves.
521	case NT_FP32:
522	fp= WaveData(p->waveHandle);
523	SetNaN64(&p->result);
524	return REQUIRES_SP_OR_DP_WAVE; //not quite true, but good enough for now AJJ 4/23/07
525	case NT_FP64:
526	dp= WaveData(p->waveHandle);
527	par_values = WaveData(p->par_values);
528	weight_values = WaveData(p->weight_values);
529
530	for(i=0; i<16; i++) {
531	pars[i] = dp[i];
532	}
533
534	sum = 0.0;
535	n_slices = (int)floor(dp[16]);
536
537	if(n_slices == 0) {
538	p->result =weight_dispersion( &disperse_elliptical_cylinder_analytical_2D,
539	par_values, weight_values,
540	(int)floor(pars[14]), (int)floor(pars[15]),
541	pars, q, phi );
542	} else {
543	for(i_theta=0; i_theta<n_slices; i_theta++) {
544	SpinProcess();
545	// average over theta
546
547	// For the cylinder model, theta_cyl=90 degrees
548	// will produce a NAN at phi_cyl=0 and 45 degrees
549	// TODO: integrate from 0 to pi/2 instead of 0 to pi
550
551	pars[6] = acos(-1.0)/n_slices * i_theta;
552
553	if( fabs(i_theta / n_slices) - 0.5 < 0.000001 ) {
554	//continue;
555	pars[6] += 0.00001;
556	}
557
558	// Multiply by sin(theta) because we are integrating in
559	// spherical coordinates
560	sum += sin(pars[6])* weight_dispersion( &disperse_elliptical_cylinder_analytical_2D,
561	par_values, weight_values,
562	(int)floor(pars[14]), (int)floor(pars[15]),
563	pars, q, phi );
564	}
565
566	p->result = sum/n_slices;
567	}
568
569	return 0;
570	default: // We can't handle this wave data type.
571	SetNaN64(&p->result);
572	return REQUIRES_SP_OR_DP_WAVE;
573	}
574
575	return 0;
576	}
577
578	int
579	Sphere_2D(FitParams2DPtr p)
580	{
581	double *dp;
582
583	float *fp; // Pointer to single precision wave data.
584	double qx;
585	double qy;
586	// double q, phi;
587	double pars[5];
588	int i;
589	// char buf[256];
590
591	if (p->waveHandle == NIL) {
592	SetNaN64(&p->result);
593	return NON_EXISTENT_WAVE;
594	}
595
596	qx = p->qx;
597	qy = p->qy;
598
599	// sprintf(buf, "Qx = %g, Qy = %g\r",qx, qy);
600	// XOPNotice(buf);
601
602
603	//not needed for a symmetric scattering function like this
604	// q = hypot(qx,qy);
605	// phi = atan2(qy,qx);
606
607	switch(WaveType(p->waveHandle)){ // We can handle single and double precision coefficient waves.
608	case NT_FP32:
609	fp= WaveData(p->waveHandle);
610	SetNaN64(&p->result);
611	return REQUIRES_SP_OR_DP_WAVE; //not quite true, but good enough for now AJJ 4/23/07
612	case NT_FP64:
613	dp= WaveData(p->waveHandle);
614
615	for(i=0; i<5; i++) {
616	pars[i] = dp[i];
617	}
618
619	p->result = SphereForm(pars, sqrt(qxqx+qyqy)); //kind of a trivial example...
620	return 0;
621	default: // We can't handle this wave data type.
622	SetNaN64(&p->result);
623	return REQUIRES_SP_OR_DP_WAVE;
624	}
625
626	return 0;
627	}

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

Context Navigation

source: sans/XOP_Dev/SANSAnalysis/XOP/Func2D.c @ 756

Download in other formats: