source: sans/Analysis/branches/ajj_23APR07/IGOR_Package_Files/Put in User Procedures/SANS_Models_v3.00/Models_2D/Ellipsoid2D.ipf @ 198

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

Adding three new 2D functions, the easy ones that use a Gaussian distribution of parameters.

Changed the parsing of constraints to only build a constraint wave for those parameters that are not being held. failing to do so results in an "unknown error" from curve fitting. This fix was done in both the 1D and 2D wrappers.

Changed the default mask that was two pixels all around to now mask out the beamstop for 2D data sets. 8 pixel radius is currently hard-wired in PlotUtils2D - and is generated as the fit is dispatched, automatically determining the beam center from the data set.

File size: 4.4 KB
Line 
1#pragma rtGlobals=1             // Use modern global access method.
2#pragma IgorVersion=6.0
3
4//
5// The plotting macro sets up TWO dependencies
6// - one for the triplet calculation
7// - one for a matrix to display, a copy of the triplet
8//
9// For display, there are two copies of the matrix. One matrix is linear, and is a copy of the
10// triplet (which is ALWAYS linear). The other matrix is toggled log/lin for display
11// in the same way the 2D SANS data matrix is handled.
12//
13
14///  REQUIRES DANSE XOP for 2D FUNCTIONS
15
16//
17// the calculation is done as for the QxQy data set:
18// three waves XYZ, then converted to a matrix
19//
20Proc PlotEllipsoid2D(str)                                               
21        String str
22        Prompt str,"Pick the data folder containing the 2D data",popup,getAList(4)
23       
24        SetDataFolder $("root:"+str)
25       
26        // Setup parameter table for model function
27        make/O/T/N=12 parameters_Ellip2D
28        Make/O/D/N=12 coef_Ellip2D
29       
30        coef_Ellip2D[0] = 1.0
31        coef_Ellip2D[1] = 20.0
32        coef_Ellip2D[2] = 60.0
33        coef_Ellip2D[3] = 1.0e-6
34        coef_Ellip2D[4] = 0.0
35        coef_Ellip2D[5] = 1.57
36        coef_Ellip2D[6] = 0.0
37        coef_Ellip2D[7] = 0.0
38        coef_Ellip2D[8] = 0.0
39        coef_Ellip2D[9] = 0.0
40        coef_Ellip2D[10] = 0.0
41        coef_Ellip2D[11] = 10
42       
43        parameters_Ellip2D[0] = "Scale"
44        parameters_Ellip2D[1] = "Radius_a"
45        parameters_Ellip2D[2] = "Radius_b (rotation axis)"
46        parameters_Ellip2D[3] = "Contrast"
47        parameters_Ellip2D[4] = "Background"
48        parameters_Ellip2D[5] = "Axis Theta"
49        parameters_Ellip2D[6] = "Axis Phi"
50        parameters_Ellip2D[7] = "Sigma of polydisp in R_a [Angstrom]"
51        parameters_Ellip2D[8] = "Sigma of polydisp in R_b [Angstrom]"
52        parameters_Ellip2D[9] = "Sigma of polydisp in Theta [rad]"
53        parameters_Ellip2D[10] = "Sigma of polydisp in Phi [rad]"
54        parameters_Ellip2D[11] = "Num of polydisp points"
55
56       
57        Edit parameters_Ellip2D,coef_Ellip2D                                   
58       
59        // generate the triplet representation
60        Duplicate/O $(str+"_qx") xwave_Ellip2D
61        Duplicate/O $(str+"_qy") ywave_Ellip2D,zwave_Ellip2D                   
62               
63        Variable/G gs_Ellip2D=0
64        gs_Ellip2D := Ellipsoid2D(coef_Ellip2D,zwave_Ellip2D,xwave_Ellip2D,ywave_Ellip2D)       //AAO 2D calculation
65       
66        Display ywave_Ellip2D vs xwave_Ellip2D
67        modifygraph log=0
68        ModifyGraph mode=3,marker=16,zColor(ywave_Ellip2D)={zwave_Ellip2D,*,*,YellowHot,0}
69        ModifyGraph standoff=0
70        ModifyGraph width={Aspect,1}
71        ModifyGraph lowTrip=0.001
72        Label bottom "qx (A\\S-1\\M)"
73        Label left "qy (A\\S-1\\M)"
74        AutoPositionWindow/M=1/R=$(WinName(0,1)) $WinName(0,2)
75       
76        // generate the matrix representation
77        ConvertQxQy2Mat(xwave_Ellip2D,ywave_Ellip2D,zwave_Ellip2D,"Ellip2D_mat")
78        Duplicate/O $"Ellip2D_mat",$"Ellip2D_lin"               //keep a linear-scaled version of the data
79        // _mat is for display, _lin is the real calculation
80
81        // not a function evaluation - this simply keeps the matrix for display in sync with the triplet calculation
82        Variable/G gs_Ellip2Dmat=0
83        gs_Ellip2Dmat := UpdateQxQy2Mat(xwave_Ellip2D,ywave_Ellip2D,zwave_Ellip2D,Ellip2D_lin,Ellip2D_mat)
84       
85       
86        SetDataFolder root:
87        AddModelToStrings("Ellipsoid2D","coef_Ellip2D","Ellip2D")
88End
89
90//AAO version, uses XOP if available
91// simply calls the original single point calculation with
92// a wave assignment (this will behave nicely if given point ranges)
93//
94// NON-THREADED IMPLEMENTATION
95//
96//Function Ellipsoid2D(cw,zw,xw,yw) : FitFunc
97//      Wave cw,zw,xw,yw
98//     
99//#if exists("EllipsoidModel_D")
100//      zw = EllipsoidModel_D(cw,xw,yw)
101//#else
102//      Abort "You do not have the SANS Analysis XOP installed"
103//#endif
104//      return(0)
105//End
106//
107
108//threaded version of the function
109ThreadSafe Function Ellipsoid2D_T(cw,zw,xw,yw,p1,p2)
110        WAVE cw,zw,xw,yw
111        Variable p1,p2
112       
113#if exists("EllipsoidModel_D")                  //to hide the function if XOP not installed
114        zw[p1,p2]= EllipsoidModel_D(cw,xw,yw)
115#endif
116
117        return 0
118End
119
120//
121//  Fit function that is actually a wrapper to dispatch the calculation to N threads
122//
123// nthreads is 1 or an even number, typically 2
124// it doesn't matter if npt is odd. In this case, fractional point numbers are passed
125// and the wave indexing works just fine - I tested this with test waves of 7 and 8 points
126// and the points "2.5" and "3.5" evaluate correctly as 2 and 3
127//
128Function Ellipsoid2D(cw,zw,xw,yw) : FitFunc
129        Wave cw,zw,xw,yw
130       
131        Variable npt=numpnts(yw)
132        Variable i,nthreads= ThreadProcessorCount
133        variable mt= ThreadGroupCreate(nthreads)
134
135        for(i=0;i<nthreads;i+=1)
136        //      Print (i*npt/nthreads),((i+1)*npt/nthreads-1)
137                ThreadStart mt,i,Ellipsoid2D_T(cw,zw,xw,yw,(i*npt/nthreads),((i+1)*npt/nthreads-1))
138        endfor
139
140        do
141                variable tgs= ThreadGroupWait(mt,100)
142        while( tgs != 0 )
143
144        variable dummy= ThreadGroupRelease(mt)
145       
146        return(0)
147End
Note: See TracBrowser for help on using the repository browser.