Ignore:
Timestamp:
Dec 6, 2010 12:51:24 PM (12 years ago)
Author:
srkline
Message:

Adding new model functions to the picker list

Adding resolution-smeared functions to the 2D analysis models

Constraints are now functional during a 2D fit

vesib.cor example data is now 6-column

Location:
sans/Dev/trunk/NCNR_User_Procedures/Analysis
Files:
9 edited

Legend:

Unmodified
Added
Removed
  • sans/Dev/trunk/NCNR_User_Procedures/Analysis/Models/Models_2D/CoreShellCyl2D_v40.ipf

    r753 r771  
    100100End 
    101101 
    102 //AAO version, uses XOP if available 
    103 // simply calls the original single point calculation with 
    104 // a wave assignment (this will behave nicely if given point ranges) 
    105 // 
    106 // NON-THREADED IMPLEMENTATION 
    107 // 
    108 //Function CoreShellCylinder2D(cw,zw,xw,yw) : FitFunc 
    109 //      Wave cw,zw,xw,yw 
    110 //       
    111 //#if exists("CoreShellCylinderModel_D") 
    112 //      Make/O/D/N=15 CSCyl2D_tmp 
    113 //      CSCyl2D_tmp = cw 
    114 //      CSCyl2D_tmp[14] = 25 
    115 // 
    116 //      zw = CoreShellCylinderModel_D(CSCyl2D_tmp,xw,yw) 
    117 //       
    118 ////    zw = CoreShellCylinderModel_D(cw,xw,yw) 
    119 //#else 
    120 //      Abort "You do not have the SANS Analysis XOP installed" 
    121 //#endif 
    122 //      return(0) 
    123 //End 
    124 // 
    125  
    126 ////threaded version of the function 
    127 //ThreadSafe Function CoreShellCylinder2D_T(cw,zw,xw,yw,p1,p2) 
    128 //      WAVE cw,zw,xw,yw 
    129 //      Variable p1,p2 
    130 //       
    131 //#if exists("CoreShellCylinder_2DX")                   //to hide the function if XOP not installed 
    132 // 
    133 //      Make/O/D/N=15 CSCyl2D_tmp 
    134 //      CSCyl2D_tmp = cw 
    135 //      CSCyl2D_tmp[14] = 25 
    136 //      CSCyl2D_tmp[7] = 0              //send a zero background to the calculation, add it in later 
    137 // 
    138 //      zw[p1,p2]= CoreShellCylinder_2DX(CSCyl2D_tmp,xw,yw) + cw[7] 
    139 //       
    140 //#endif 
    141 // 
    142 //      return 0 
    143 //End 
    144 // 
    145 //// 
    146 ////  Fit function that is actually a wrapper to dispatch the calculation to N threads 
    147 //// 
    148 //// nthreads is 1 or an even number, typically 2 
    149 //// it doesn't matter if npt is odd. In this case, fractional point numbers are passed 
    150 //// and the wave indexing works just fine - I tested this with test waves of 7 and 8 points 
    151 //// and the points "2.5" and "3.5" evaluate correctly as 2 and 3 
    152 //// 
    153 //Function CoreShellCylinder2D(cw,zw,xw,yw) : FitFunc 
    154 //      Wave cw,zw,xw,yw 
    155 //       
    156 //      Variable npt=numpnts(yw) 
    157 //      Variable i,nthreads= ThreadProcessorCount 
    158 //      variable mt= ThreadGroupCreate(nthreads) 
    159 // 
    160 //      for(i=0;i<nthreads;i+=1) 
    161 //      //      Print (i*npt/nthreads),((i+1)*npt/nthreads-1) 
    162 //              ThreadStart mt,i,CoreShellCylinder2D_T(cw,zw,xw,yw,(i*npt/nthreads),((i+1)*npt/nthreads-1)) 
    163 //      endfor 
    164 // 
    165 //      do 
    166 //              variable tgs= ThreadGroupWait(mt,100) 
    167 //      while( tgs != 0 ) 
    168 // 
    169 //      variable dummy= ThreadGroupRelease(mt) 
    170 //       
    171 //      return(0) 
    172 //End 
     102 
     103// 
     104Proc PlotSmearedCoreShellCylinder2D(str)                                                 
     105        String str 
     106        Prompt str,"Pick the data folder containing the 2D data",popup,getAList(4) 
     107 
     108        if (!exists("CoreShellCylinder_2DX")) 
     109                Abort "You must have the SANSAnalysis XOP installed to use 2D models" 
     110        endif 
     111         
     112        SetDataFolder $("root:"+str) 
     113         
     114        // Setup parameter table for model function 
     115//      make/O/T/N=15 parameters_CSCyl2D 
     116//      Make/O/D/N=15 coef_CSCyl2D 
     117        make/O/T/N=14 smear_parameters_CSCyl2D 
     118        Make/O/D/N=14 smear_coef_CSCyl2D 
     119         
     120        smear_coef_CSCyl2D[0] = 1.0 
     121        smear_coef_CSCyl2D[1] = 20.0 
     122        smear_coef_CSCyl2D[2] = 10.0 
     123        smear_coef_CSCyl2D[3] = 200.0 
     124        smear_coef_CSCyl2D[4] = 1.0e-6 
     125        smear_coef_CSCyl2D[5] = 4.0e-6 
     126        smear_coef_CSCyl2D[6] = 1.0e-6 
     127        smear_coef_CSCyl2D[7] = 0.0 
     128        smear_coef_CSCyl2D[8] = 1.57 
     129        smear_coef_CSCyl2D[9] = 0.0 
     130        smear_coef_CSCyl2D[10] = 0.0 
     131        smear_coef_CSCyl2D[11] = 0.0 
     132        smear_coef_CSCyl2D[12] = 0.0 
     133        smear_coef_CSCyl2D[13] = 0.0 
     134        //hard-wire the number of integration points 
     135//      smear_coef_CSCyl2D[14] = 10 
     136         
     137        smear_parameters_CSCyl2D[0] = "Scale" 
     138        smear_parameters_CSCyl2D[1] = "Radius" 
     139        smear_parameters_CSCyl2D[2] = "Thickness" 
     140        smear_parameters_CSCyl2D[3] = "Length" 
     141        smear_parameters_CSCyl2D[4] = "Core SLD" 
     142        smear_parameters_CSCyl2D[5] = "Shell SLD" 
     143        smear_parameters_CSCyl2D[6] = "Solvent SLD" 
     144        smear_parameters_CSCyl2D[7] = "Background" 
     145        smear_parameters_CSCyl2D[8] = "Axis Theta" 
     146        smear_parameters_CSCyl2D[9] = "Axis Phi" 
     147        smear_parameters_CSCyl2D[10] = "Sigma of polydisp in Radius [Angstrom]" 
     148        smear_parameters_CSCyl2D[11] = "Sigma of polydisp in Thickness [Angstrom]" 
     149        smear_parameters_CSCyl2D[12] = "Sigma of polydisp in Theta [rad]" 
     150        smear_parameters_CSCyl2D[13] = "Sigma of polydisp in Phi [rad]" 
     151//      smear_parameters_CSCyl2D[14] = "Num of polydisp points" 
     152         
     153        Edit smear_parameters_CSCyl2D,smear_coef_CSCyl2D                                         
     154         
     155        // generate the triplet representation 
     156        Duplicate/O $(str+"_qx") smeared_CSCyl2D 
     157        SetScale d,0,0,"1/cm",smeared_CSCyl2D                                    
     158                 
     159        Variable/G gs_CSCyl2D=0 
     160        gs_CSCyl2D := fSmearedCoreShellCylinder2D(smear_coef_CSCyl2D,smeared_CSCyl2D)           //wrapper to fill the STRUCT 
     161         
     162        Display $(str+"_qy") vs $(str+"_qx") 
     163        modifygraph log=0 
     164        ModifyGraph mode=3,marker=16,zColor($(str+"_qy"))={smeared_CSCyl2D,*,*,YellowHot,0} 
     165        ModifyGraph standoff=0 
     166        ModifyGraph width={Aspect,1} 
     167        ModifyGraph lowTrip=0.001 
     168        Label bottom "qx (A\\S-1\\M)" 
     169        Label left "qy (A\\S-1\\M)" 
     170        AutoPositionWindow/M=1/R=$(WinName(0,1)) $WinName(0,2) 
     171         
     172        // generate the matrix representation 
     173        Duplicate/O $(str+"_qx"), sm_qx 
     174        Duplicate/O $(str+"_qy"), sm_qy         // I can't use local variables in dependencies, so I need the name (that I can't get) 
     175         
     176        // generate the matrix representation 
     177        ConvertQxQy2Mat(sm_qx,sm_qy,smeared_CSCyl2D,"sm_CSCyl2D_mat") 
     178        Duplicate/O $"sm_CSCyl2D_mat",$"sm_CSCyl2D_lin"                 //keep a linear-scaled version of the data 
     179        // _mat is for display, _lin is the real calculation 
     180 
     181        // not a function evaluation - this simply keeps the matrix for display in sync with the triplet calculation 
     182        Variable/G gs_CSCyl2Dmat=0 
     183        gs_CSCyl2Dmat := UpdateQxQy2Mat(sm_qx,sm_qy,smeared_CSCyl2D,sm_CSCyl2D_lin,sm_CSCyl2D_mat) 
     184         
     185         
     186        SetDataFolder root: 
     187        AddModelToStrings("SmearedCoreShellCylinder2D","smear_coef_CSCyl2D","smear_parameters_CSCyl2D","CSCyl2D") 
     188End 
     189 
    173190 
    174191// 
     
    186203 
    187204        Make/O/D/N=15 CSCyl2D_tmp 
    188         CSCyl2D_tmp = cw 
     205        CSCyl2D_tmp[0,13] = cw 
    189206        CSCyl2D_tmp[14] = 25 
    190207        CSCyl2D_tmp[7] = 0              //send a zero background to the calculation, add it in later 
     
    196213        return(0) 
    197214End 
     215 
     216 
     217 
     218 
     219/////////////////////smeared functions ////////////////////// 
     220 
     221Function SmearedCoreShellCylinder2D(s) 
     222        Struct ResSmear_2D_AAOStruct &s 
     223         
     224//// non-threaded, but generic calculation 
     225//// the last param is nord      
     226//      Smear_2DModel_PP(CoreShellCylinder2D_noThread,s,10) 
     227 
     228 
     229//// the last param is nord 
     230        SmearedCoreShellCylinder2D_THR(s,10)             
     231 
     232        return(0) 
     233end 
     234 
     235// for the plot dependency only 
     236Function fSmearedCoreShellCylinder2D(coefW,resultW) 
     237        Wave coefW,resultW 
     238         
     239        String str = getWavesDataFolder(resultW,0) 
     240        String DF="root:"+str+":" 
     241         
     242        WAVE qx = $(DF+str+"_qx") 
     243        WAVE qy = $(DF+str+"_qy") 
     244        WAVE qz = $(DF+str+"_qz") 
     245        WAVE sQpl = $(DF+str+"_sQpl") 
     246        WAVE sQpp = $(DF+str+"_sQpp") 
     247        WAVE shad = $(DF+str+"_fs") 
     248         
     249        STRUCT ResSmear_2D_AAOStruct s 
     250        WAVE s.coefW = coefW     
     251        WAVE s.zw = resultW      
     252        WAVE s.xw[0] = qx 
     253        WAVE s.xw[1] = qy 
     254        WAVE s.qz = qz 
     255        WAVE s.sQpl = sQpl 
     256        WAVE s.sQpp = sQpp 
     257        WAVE s.fs = shad 
     258         
     259        Variable err 
     260        err = SmearedCoreShellCylinder2D(s) 
     261         
     262        return (0) 
     263End 
     264 
     265// 
     266// NON-THREADED IMPLEMENTATION 
     267// -- same as threaded, but no MultiThread KW 
     268// 
     269ThreadSafe Function CoreShellCylinder2D_noThread(cw,zw,xw,yw) 
     270        Wave cw,zw,xw,yw 
     271         
     272        //      Variable t1=StopMSTimer(-2) 
     273 
     274#if exists("CoreShellCylinder_2DX")                     //to hide the function if XOP not installed 
     275 
     276        Make/O/D/N=15 CSCyl2D_tmp 
     277        CSCyl2D_tmp[0,13] = cw 
     278        CSCyl2D_tmp[14] = 25 
     279        CSCyl2D_tmp[7] = 0              //send a zero background to the calculation, add it in later 
     280 
     281        zw= CoreShellCylinder_2DX(CSCyl2D_tmp,xw,yw) + cw[7] 
     282         
     283#endif 
     284 
     285//      Print "elapsed time = ",(StopMSTimer(-2) - t1)/1e6 
     286         
     287        return(0) 
     288End 
     289 
     290// 
     291// this is the threaded version, that dispatches the calculation out to threads 
     292// 
     293// must be written specific to each 2D function 
     294//  
     295Function SmearedCoreShellCylinder2D_THR(s,nord) 
     296        Struct ResSmear_2D_AAOStruct &s 
     297        Variable nord 
     298         
     299        String weightStr,zStr 
     300         
     301// create all of the necessary quadrature waves here - rather than inside a threadsafe function 
     302        switch(nord)     
     303                case 5:          
     304                        weightStr="gauss5wt" 
     305                        zStr="gauss5z" 
     306                        if (WaveExists($weightStr) == 0) 
     307                                Make/O/D/N=(nord) $weightStr,$zStr 
     308                                Make5GaussPoints($weightStr,$zStr)       
     309                        endif 
     310                        break                            
     311                case 10:                 
     312                        weightStr="gauss10wt" 
     313                        zStr="gauss10z" 
     314                        if (WaveExists($weightStr) == 0) 
     315                                Make/O/D/N=(nord) $weightStr,$zStr 
     316                                Make10GaussPoints($weightStr,$zStr)      
     317                        endif 
     318                        break                            
     319                case 20:                 
     320                        weightStr="gauss20wt" 
     321                        zStr="gauss20z" 
     322                        if (WaveExists($weightStr) == 0) 
     323                                Make/O/D/N=(nord) $weightStr,$zStr 
     324                                Make20GaussPoints($weightStr,$zStr)      
     325                        endif 
     326                        break 
     327                default:                                                         
     328                        Abort "Smear_2DModel_PP_Threaded called with invalid nord value"                                         
     329        endswitch 
     330         
     331        Wave/Z wt = $weightStr 
     332        Wave/Z xi = $zStr               // wave references to pass 
     333 
     334        Variable npt=numpnts(s.xw[0]) 
     335        Variable i,nthreads= ThreadProcessorCount 
     336        variable mt= ThreadGroupCreate(nthreads) 
     337 
     338        Variable t1=StopMSTimer(-2) 
     339         
     340        for(i=0;i<nthreads;i+=1) 
     341//              Print (i*npt/nthreads),((i+1)*npt/nthreads-1) 
     342                ThreadStart mt,i,SmearedCoreShellCylinder2D_T(s.coefW,s.xw[0],s.xw[1],s.qz,s.sQpl,s.sQpp,s.fs,s.zw,wt,xi,(i*npt/nthreads),((i+1)*npt/nthreads-1),nord) 
     343        endfor 
     344 
     345        do 
     346                variable tgs= ThreadGroupWait(mt,100) 
     347        while( tgs != 0 ) 
     348 
     349        variable dummy= ThreadGroupRelease(mt) 
     350         
     351// comment out the threading + uncomment this for testing to make sure that the single thread works 
     352//      nThreads=1 
     353//      SmearedCoreShellCylinder2D_T(s.coefW,s.xw[0],s.xw[1],s.qz,s.sQpl,s.sQpp,s.fs,s.zw,wt,xi,(i*npt/nthreads),((i+1)*npt/nthreads-1),nord) 
     354 
     355        Print "elapsed time = ",(StopMSTimer(-2) - t1)/1e6 
     356         
     357        return(0) 
     358end 
     359 
     360// 
     361// - worker function for threads of Sphere2D 
     362// 
     363ThreadSafe Function SmearedCoreShellCylinder2D_T(coef,qxw,qyw,qzw,sxw,syw,fsw,zw,wt,xi,pt1,pt2,nord) 
     364        WAVE coef,qxw,qyw,qzw,sxw,syw,fsw,zw,wt,xi 
     365        Variable pt1,pt2,nord 
     366         
     367// now passed in.... 
     368//      Wave wt = $weightStr 
     369//      Wave xi = $zStr          
     370 
     371        Variable ii,jj,kk,num 
     372        Variable qx,qy,qz,qval,sx,sy,fs 
     373        Variable qy_pt,qx_pt,res_x,res_y,answer,sumIn,sumOut 
     374         
     375        Variable normFactor,phi,theta,maxSig,numStdDev=3 
     376         
     377/// keep these waves local 
     378        Make/O/D/N=(nord) fcnRet,xptW,res_tot,yptW 
     379         
     380        // now just loop over the points as specified 
     381         
     382        answer=0 
     383         
     384        Variable spl,spp,apl,app,bpl,bpp,phi_pt,qpl_pt 
     385 
     386        //loop over q-values 
     387        for(ii=pt1;ii<(pt2+1);ii+=1) 
     388                 
     389                qx = qxw[ii] 
     390                qy = qyw[ii] 
     391                qz = qzw[ii] 
     392                qval = sqrt(qx^2+qy^2+qz^2) 
     393                spl = sxw[ii] 
     394                spp = syw[ii] 
     395                fs = fsw[ii] 
     396                 
     397                normFactor = 2*pi*spl*spp 
     398                 
     399                phi = FindPhi(qx,qy) 
     400                 
     401                apl = -numStdDev*spl + qval             //parallel = q integration limits 
     402                bpl = numStdDev*spl + qval 
     403                app = -numStdDev*spp + phi              //perpendicular = phi integration limits 
     404                bpp = numStdDev*spp + phi 
     405                 
     406                //make sure the limits are reasonable. 
     407                if(apl < 0) 
     408                        apl = 0 
     409                endif 
     410                // do I need to specially handle limits when phi ~ 0? 
     411         
     412                 
     413                sumOut = 0 
     414                for(jj=0;jj<nord;jj+=1)         // call phi the "outer' 
     415                        phi_pt = (xi[jj]*(bpp-app)+app+bpp)/2 
     416 
     417                        sumIn=0 
     418                        for(kk=0;kk<nord;kk+=1)         //at phi, integrate over Qpl 
     419 
     420                                qpl_pt = (xi[kk]*(bpl-apl)+apl+bpl)/2 
     421                                 
     422                                FindQxQy(qpl_pt,phi_pt,qx_pt,qy_pt)             //find the corresponding QxQy to the Q,phi 
     423                                yPtw[kk] = qy_pt                                        //phi is the same in this loop, but qy is not 
     424                                xPtW[kk] = qx_pt                                        //qx is different here too, as we're varying Qpl 
     425                                 
     426                                res_tot[kk] = exp(-0.5*( (qpl_pt-qval)^2/spl/spl + (phi_pt-phi)^2/spp/spp ) ) 
     427                                res_tot[kk] /= normFactor 
     428//                              res_tot[kk] *= fs 
     429 
     430                        endfor 
     431                         
     432                        CoreShellCylinder2D_noThread(coef,fcnRet,xptw,yptw)                     //fcn passed in is an AAO 
     433                         
     434                        //sumIn += wt[jj]*wt[kk]*res_tot*fcnRet[0] 
     435                        fcnRet *= wt[jj]*wt*res_tot 
     436                        // 
     437                        answer += (bpl-apl)/2.0*sum(fcnRet)             // 
     438                endfor 
     439 
     440                answer *= (bpp-app)/2.0 
     441                zw[ii] = answer 
     442        endfor 
     443         
     444        return(0) 
     445end 
  • sans/Dev/trunk/NCNR_User_Procedures/Analysis/Models/Models_2D/Cylinder_2D_v40.ipf

    r753 r771  
    1414///  REQUIRES DANSE XOP for 2D FUNCTIONS 
    1515 
     16// 
     17// 4.8x speedup for threading (N=4 + Nvirt=4) 
     18// 
    1619 
    1720// 
     
    3033        // NOTE THAT THE COEFFICIENTS [N] ARE IN A DIFFERENT ORDER !!! 
    3134        // Setup parameter table for model function 
    32 //      make/O/T/N=11 parameters_Cyl2D 
    33 //      Make/O/D/N=11 coef_Cyl2D 
    3435        make/O/T/N=11 parameters_Cyl2D 
    3536        Make/O/D/N=11 coef_Cyl2D 
     
    9798End 
    9899 
    99 //AAO version, uses XOP if available 
    100 // simply calls the original single point calculation with 
    101 // a wave assignment (this will behave nicely if given point ranges) 
    102 // 
    103 // NON-THREADED IMPLEMENTATION 
    104 // 
    105 //Function Cylinder2D(cw,zw,xw,yw) : FitFunc 
    106 //      Wave cw,zw,xw,yw 
    107 //       
    108 //#if exists("CylinderModel_D") 
    109 // 
    110 //      Make/O/D/N=12 Cyl2D_tmp                         // there seems to be no speed penalty for doing this... 
    111 //      Cyl2D_tmp = cw 
    112 //      Cyl2D_tmp[11] = 25                                      // hard-wire the number of integration points 
    113 //       
    114 //      zw= CylinderModel_D(Cyl2D_tmp,xw,yw) 
    115 // 
    116 //      //zw = CylinderModel_D(cw,xw,yw) 
    117 //#else 
    118 //      Abort "You do not have the SANS Analysis XOP installed" 
    119 //#endif 
    120 //      return(0) 
    121 //End 
    122 // 
    123  
    124 ////threaded version of the function 
    125 //ThreadSafe Function Cylinder2D_T(cw,zw,xw,yw,p1,p2) 
    126 //      WAVE cw,zw,xw,yw 
    127 //      Variable p1,p2 
    128 //       
    129 //#if exists("Cylinder_2DX")                    //to hide the function if XOP not installed 
    130 // 
    131 //      Make/O/D/N=12 Cyl2D_tmp                         // there seems to be no speed penalty for doing this... 
    132 //      Cyl2D_tmp = cw 
    133 //      Cyl2D_tmp[11] = 25                                      // hard-wire the number of integration points 
    134 //      Cyl2D_tmp[5] = 0                                                // send a background of zero 
    135 //       
    136 //      zw[p1,p2]= Cylinder_2DX(Cyl2D_tmp,xw,yw) + cw[5]                //add in the proper background here 
    137 // 
    138 //#endif 
    139 // 
    140 //      return 0 
    141 //End 
    142 // 
     100 
     101// - sets up a dependency to a wrapper, not the actual SmearedModelFunction 
     102Proc PlotSmearedCylinder2D(str)                                          
     103        String str 
     104        Prompt str,"Pick the data folder containing the 2D data",popup,getAList(4) 
     105         
     106        if (!exists("Cylinder_2DX") ) 
     107                Abort "You must have the SANSAnalysis XOP installed to use 2D models" 
     108        endif 
     109        SetDataFolder $("root:"+str) 
     110         
     111        // NOTE THAT THE COEFFICIENTS [N] ARE IN A DIFFERENT ORDER !!! 
     112        // Setup parameter table for model function 
     113        make/O/T/N=11 smear_parameters_Cyl2D 
     114        Make/O/D/N=11 smear_coef_Cyl2D 
     115        smear_coef_Cyl2D[0] = 1.0 
     116        smear_coef_Cyl2D[1] = 20.0 
     117        smear_coef_Cyl2D[2] = 60.0 
     118        smear_coef_Cyl2D[3] = 1e-6 
     119        smear_coef_Cyl2D[4] = 6.3e-6 
     120        smear_coef_Cyl2D[5] = 0.0 
     121        smear_coef_Cyl2D[6] = 1.57 
     122        smear_coef_Cyl2D[7] = 0.0 
     123        smear_coef_Cyl2D[8] = 0.0 
     124        smear_coef_Cyl2D[9] = 0.0 
     125        smear_coef_Cyl2D[10] = 0.0 
     126         
     127        // currently, the number of integration points is hard-wired to be 25 in Cylinder2D_T 
     128        //smear_coef_Cyl2D[11] = 25 
     129        // 
     130        smear_parameters_Cyl2D[0] = "Scale" 
     131        smear_parameters_Cyl2D[1] = "Radius" 
     132        smear_parameters_Cyl2D[2] = "Length" 
     133        smear_parameters_Cyl2D[3] = "SLD cylinder (A^-2)" 
     134        smear_parameters_Cyl2D[4] = "SLD solvent" 
     135        smear_parameters_Cyl2D[5] = "Background" 
     136        smear_parameters_Cyl2D[6] = "Axis Theta" 
     137        smear_parameters_Cyl2D[7] = "Axis Phi" 
     138         
     139        smear_parameters_Cyl2D[9] = "Sigma of polydisp in Theta [rad]"          //***** 
     140        smear_parameters_Cyl2D[10] = "Sigma of polydisp in Phi [rad]"                   //***** 
     141        smear_parameters_Cyl2D[8] = "Sigma of polydisp in Radius [A]"           //***** 
     142         
     143//      smear_parameters_Cyl2D[11] = "number of integration points" 
     144 
     145        Edit smear_parameters_Cyl2D,smear_coef_Cyl2D                                     
     146         
     147        // generate the triplet representation 
     148        Duplicate/O $(str+"_qx") smeared_Cyl2D 
     149        SetScale d,0,0,"1/cm",smeared_Cyl2D                                      
     150                 
     151        Variable/G gs_Cyl2D=0 
     152        gs_Cyl2D := fSmearedCylinder2D(smear_coef_Cyl2D,smeared_Cyl2D)          //wrapper to fill the STRUCT 
     153         
     154        Display $(str+"_qy") vs $(str+"_qx") 
     155        modifygraph log=0 
     156        ModifyGraph mode=3,marker=16,zColor($(str+"_qy"))={smeared_Cyl2D,*,*,YellowHot,0} 
     157        ModifyGraph standoff=0 
     158        ModifyGraph width={Aspect,1} 
     159        ModifyGraph lowTrip=0.001 
     160        Label bottom "qx (A\\S-1\\M)" 
     161        Label left "qy (A\\S-1\\M)" 
     162        AutoPositionWindow/M=1/R=$(WinName(0,1)) $WinName(0,2) 
     163         
     164        // generate the matrix representation 
     165        Duplicate/O $(str+"_qx"), sm_qx 
     166        Duplicate/O $(str+"_qy"), sm_qy         // I can't use local variables in dependencies, so I need the name (that I can't get) 
     167         
     168        // generate the matrix representation 
     169        ConvertQxQy2Mat(sm_qx,sm_qy,smeared_Cyl2D,"sm_Cyl2D_mat") 
     170        Duplicate/O $"sm_Cyl2D_mat",$"sm_Cyl2D_lin"             //keep a linear-scaled version of the data 
     171        // _mat is for display, _lin is the real calculation 
     172 
     173        // not a function evaluation - this simply keeps the matrix for display in sync with the triplet calculation 
     174        Variable/G gs_Cyl2Dmat=0 
     175        gs_Cyl2Dmat := UpdateQxQy2Mat(sm_qx,sm_qy,smeared_Cyl2D,sm_Cyl2D_lin,sm_Cyl2D_mat) 
     176         
     177         
     178        SetDataFolder root: 
     179        AddModelToStrings("SmearedCylinder2D","smear_coef_Cyl2D","smear_parameters_Cyl2D","Cyl2D") 
     180End 
     181 
     182 
     183 
     184 
     185 
     186 
    143187//// 
    144188////  Fit function that is actually a wrapper to dispatch the calculation to N threads 
     
    185229 
    186230        Make/O/D/N=12 Cyl2D_tmp                         // there seems to be no speed penalty for doing this... 
    187         Cyl2D_tmp = cw 
     231        Cyl2D_tmp[0,10] = cw 
    188232        Cyl2D_tmp[11] = 25                                      // hard-wire the number of integration points 
    189233        Cyl2D_tmp[5] = 0                                                // send a background of zero 
     
    197241        return(0) 
    198242End 
     243 
     244 
     245/////////////////////smeared functions ////////////////////// 
     246 
     247Function SmearedCylinder2D(s) 
     248        Struct ResSmear_2D_AAOStruct &s 
     249         
     250//// non-threaded, but generic calculation 
     251//// the last param is nord      
     252//      Smear_2DModel_PP(Cylinder2D_noThread,s,10) 
     253 
     254 
     255//// the last param is nord 
     256        SmearedCylinder2D_THR(s,10)              
     257 
     258        return(0) 
     259end 
     260 
     261// for the plot dependency only 
     262Function fSmearedCylinder2D(coefW,resultW) 
     263        Wave coefW,resultW 
     264         
     265        String str = getWavesDataFolder(resultW,0) 
     266        String DF="root:"+str+":" 
     267         
     268        WAVE qx = $(DF+str+"_qx") 
     269        WAVE qy = $(DF+str+"_qy") 
     270        WAVE qz = $(DF+str+"_qz") 
     271        WAVE sQpl = $(DF+str+"_sQpl") 
     272        WAVE sQpp = $(DF+str+"_sQpp") 
     273        WAVE shad = $(DF+str+"_fs") 
     274         
     275        STRUCT ResSmear_2D_AAOStruct s 
     276        WAVE s.coefW = coefW     
     277        WAVE s.zw = resultW      
     278        WAVE s.xw[0] = qx 
     279        WAVE s.xw[1] = qy 
     280        WAVE s.qz = qz 
     281        WAVE s.sQpl = sQpl 
     282        WAVE s.sQpp = sQpp 
     283        WAVE s.fs = shad 
     284         
     285        Variable err 
     286        err = SmearedCylinder2D(s) 
     287         
     288        return (0) 
     289End 
     290 
     291// 
     292// NON-THREADED IMPLEMENTATION 
     293// -- same as threaded, but no MultiThread KW 
     294// 
     295ThreadSafe Function Cylinder2D_noThread(cw,zw,xw,yw) 
     296        Wave cw,zw,xw,yw 
     297         
     298        //      Variable t1=StopMSTimer(-2) 
     299 
     300#if exists("Cylinder_2DX")                      //to hide the function if XOP not installed 
     301 
     302        Make/O/D/N=12 Cyl2D_tmp                         // there seems to be no speed penalty for doing this... 
     303        Cyl2D_tmp[0,10] = cw 
     304        Cyl2D_tmp[11] = 25                                      // hard-wire the number of integration points 
     305        Cyl2D_tmp[5] = 0                                                // send a background of zero 
     306         
     307        zw = Cylinder_2DX(Cyl2D_tmp,xw,yw) + cw[5]              //add in the proper background here 
     308 
     309#endif 
     310 
     311//      Print "elapsed time = ",(StopMSTimer(-2) - t1)/1e6 
     312         
     313        return(0) 
     314End 
     315 
     316// 
     317// this is the threaded version, that dispatches the calculation out to threads 
     318// 
     319// must be written specific to each 2D function 
     320//  
     321Function SmearedCylinder2D_THR(s,nord) 
     322        Struct ResSmear_2D_AAOStruct &s 
     323        Variable nord 
     324         
     325        String weightStr,zStr 
     326         
     327// create all of the necessary quadrature waves here - rather than inside a threadsafe function 
     328        switch(nord)     
     329                case 5:          
     330                        weightStr="gauss5wt" 
     331                        zStr="gauss5z" 
     332                        if (WaveExists($weightStr) == 0) 
     333                                Make/O/D/N=(nord) $weightStr,$zStr 
     334                                Make5GaussPoints($weightStr,$zStr)       
     335                        endif 
     336                        break                            
     337                case 10:                 
     338                        weightStr="gauss10wt" 
     339                        zStr="gauss10z" 
     340                        if (WaveExists($weightStr) == 0) 
     341                                Make/O/D/N=(nord) $weightStr,$zStr 
     342                                Make10GaussPoints($weightStr,$zStr)      
     343                        endif 
     344                        break                            
     345                case 20:                 
     346                        weightStr="gauss20wt" 
     347                        zStr="gauss20z" 
     348                        if (WaveExists($weightStr) == 0) 
     349                                Make/O/D/N=(nord) $weightStr,$zStr 
     350                                Make20GaussPoints($weightStr,$zStr)      
     351                        endif 
     352                        break 
     353                default:                                                         
     354                        Abort "Smear_2DModel_PP_Threaded called with invalid nord value"                                         
     355        endswitch 
     356         
     357        Wave/Z wt = $weightStr 
     358        Wave/Z xi = $zStr               // wave references to pass 
     359 
     360        Variable npt=numpnts(s.xw[0]) 
     361        Variable i,nthreads= ThreadProcessorCount 
     362        variable mt= ThreadGroupCreate(nthreads) 
     363 
     364        Variable t1=StopMSTimer(-2) 
     365         
     366        for(i=0;i<nthreads;i+=1) 
     367//              Print (i*npt/nthreads),((i+1)*npt/nthreads-1) 
     368                ThreadStart mt,i,SmearedCylinder2D_T(s.coefW,s.xw[0],s.xw[1],s.qz,s.sQpl,s.sQpp,s.fs,s.zw,wt,xi,(i*npt/nthreads),((i+1)*npt/nthreads-1),nord) 
     369        endfor 
     370 
     371        do 
     372                variable tgs= ThreadGroupWait(mt,100) 
     373        while( tgs != 0 ) 
     374 
     375        variable dummy= ThreadGroupRelease(mt) 
     376         
     377// comment out the threading + uncomment this for testing to make sure that the single thread works 
     378//      nThreads=1 
     379//      SmearedCylinder2D_T(s.coefW,s.xw[0],s.xw[1],s.qz,s.sQpl,s.sQpp,s.fs,s.zw,wt,xi,(i*npt/nthreads),((i+1)*npt/nthreads-1),nord) 
     380 
     381        Print "elapsed time = ",(StopMSTimer(-2) - t1)/1e6 
     382         
     383        return(0) 
     384end 
     385 
     386// 
     387// - worker function for threads of Sphere2D 
     388// 
     389ThreadSafe Function SmearedCylinder2D_T(coef,qxw,qyw,qzw,sxw,syw,fsw,zw,wt,xi,pt1,pt2,nord) 
     390        WAVE coef,qxw,qyw,qzw,sxw,syw,fsw,zw,wt,xi 
     391        Variable pt1,pt2,nord 
     392         
     393// now passed in.... 
     394//      Wave wt = $weightStr 
     395//      Wave xi = $zStr          
     396 
     397        Variable ii,jj,kk,num 
     398        Variable qx,qy,qz,qval,sx,sy,fs 
     399        Variable qy_pt,qx_pt,res_x,res_y,answer,sumIn,sumOut 
     400         
     401        Variable normFactor,phi,theta,maxSig,numStdDev=3 
     402         
     403/// keep these waves local 
     404        Make/O/D/N=(nord) fcnRet,xptW,res_tot,yptW 
     405         
     406        // now just loop over the points as specified 
     407         
     408        answer=0 
     409         
     410        Variable spl,spp,apl,app,bpl,bpp,phi_pt,qpl_pt 
     411 
     412        //loop over q-values 
     413        for(ii=pt1;ii<(pt2+1);ii+=1) 
     414                 
     415                qx = qxw[ii] 
     416                qy = qyw[ii] 
     417                qz = qzw[ii] 
     418                qval = sqrt(qx^2+qy^2+qz^2) 
     419                spl = sxw[ii] 
     420                spp = syw[ii] 
     421                fs = fsw[ii] 
     422                 
     423                normFactor = 2*pi*spl*spp 
     424                 
     425                phi = FindPhi(qx,qy) 
     426                 
     427                apl = -numStdDev*spl + qval             //parallel = q integration limits 
     428                bpl = numStdDev*spl + qval 
     429                app = -numStdDev*spp + phi              //perpendicular = phi integration limits 
     430                bpp = numStdDev*spp + phi 
     431                 
     432                //make sure the limits are reasonable. 
     433                if(apl < 0) 
     434                        apl = 0 
     435                endif 
     436                // do I need to specially handle limits when phi ~ 0? 
     437         
     438                 
     439                sumOut = 0 
     440                for(jj=0;jj<nord;jj+=1)         // call phi the "outer' 
     441                        phi_pt = (xi[jj]*(bpp-app)+app+bpp)/2 
     442 
     443                        sumIn=0 
     444                        for(kk=0;kk<nord;kk+=1)         //at phi, integrate over Qpl 
     445 
     446                                qpl_pt = (xi[kk]*(bpl-apl)+apl+bpl)/2 
     447                                 
     448                                FindQxQy(qpl_pt,phi_pt,qx_pt,qy_pt)             //find the corresponding QxQy to the Q,phi 
     449                                yPtw[kk] = qy_pt                                        //phi is the same in this loop, but qy is not 
     450                                xPtW[kk] = qx_pt                                        //qx is different here too, as we're varying Qpl 
     451                                 
     452                                res_tot[kk] = exp(-0.5*( (qpl_pt-qval)^2/spl/spl + (phi_pt-phi)^2/spp/spp ) ) 
     453                                res_tot[kk] /= normFactor 
     454//                              res_tot[kk] *= fs 
     455 
     456                        endfor 
     457                         
     458                        Cylinder2D_noThread(coef,fcnRet,xptw,yptw)                      //fcn passed in is an AAO 
     459                         
     460                        //sumIn += wt[jj]*wt[kk]*res_tot*fcnRet[0] 
     461                        fcnRet *= wt[jj]*wt*res_tot 
     462                        // 
     463                        answer += (bpl-apl)/2.0*sum(fcnRet)             // 
     464                endfor 
     465 
     466                answer *= (bpp-app)/2.0 
     467                zw[ii] = answer 
     468        endfor 
     469         
     470        return(0) 
     471end 
     472 
     473 
  • sans/Dev/trunk/NCNR_User_Procedures/Analysis/Models/Models_2D/Ellipsoid2D_v40.ipf

    r753 r771  
    100100End 
    101101 
    102 //AAO version, uses XOP if available 
    103 // simply calls the original single point calculation with 
    104 // a wave assignment (this will behave nicely if given point ranges) 
    105 // 
    106 // NON-THREADED IMPLEMENTATION 
    107 // 
    108 //Function Ellipsoid2D(cw,zw,xw,yw) : FitFunc 
    109 //      Wave cw,zw,xw,yw 
    110 //       
    111 //#if exists("EllipsoidModel_D") 
    112 // 
    113 //      Make/O/D/N=13 Ellip2D_tmp 
    114 //      Ellip2D_tmp = cw 
    115 //      Ellip2D_tmp[12] = 25 
    116 //       
    117 //      zw = EllipsoidModel_D(Ellip2D_tmp,xw,yw) 
    118 //       
    119 ////    zw = EllipsoidModel_D(cw,xw,yw) 
    120 //#else 
    121 //      Abort "You do not have the SANS Analysis XOP installed" 
    122 //#endif 
    123 //      return(0) 
    124 //End 
    125 //// 
    126 // 
    127 ////threaded version of the function 
    128 //ThreadSafe Function Ellipsoid2D_T(cw,zw,xw,yw,p1,p2) 
    129 //      WAVE cw,zw,xw,yw 
    130 //      Variable p1,p2 
    131 //       
    132 //#if exists("Ellipsoid_2DX")                   //to hide the function if XOP not installed 
    133 // 
    134 //      Make/O/D/N=13 Ellip2D_tmp 
    135 //      Ellip2D_tmp = cw 
    136 //      Ellip2D_tmp[12] = 25 
    137 //      Ellip2D_tmp[5] = 0              //pass in a zero background and add it in later 
    138 //       
    139 //      zw[p1,p2]= Ellipsoid_2DX(Ellip2D_tmp,xw,yw) + cw[5] 
    140 //       
    141 //#endif 
    142 // 
    143 //      return 0 
    144 //End 
    145 // 
    146 //// 
    147 ////  Fit function that is actually a wrapper to dispatch the calculation to N threads 
    148 //// 
    149 //// nthreads is 1 or an even number, typically 2 
    150 //// it doesn't matter if npt is odd. In this case, fractional point numbers are passed 
    151 //// and the wave indexing works just fine - I tested this with test waves of 7 and 8 points 
    152 //// and the points "2.5" and "3.5" evaluate correctly as 2 and 3 
    153 //// 
    154 //Function Ellipsoid2D(cw,zw,xw,yw) : FitFunc 
    155 //      Wave cw,zw,xw,yw 
    156 //       
    157 //      Variable npt=numpnts(yw) 
    158 //      Variable i,nthreads= ThreadProcessorCount 
    159 //      variable mt= ThreadGroupCreate(nthreads) 
    160 // 
    161 //      for(i=0;i<nthreads;i+=1) 
    162 //      //      Print (i*npt/nthreads),((i+1)*npt/nthreads-1) 
    163 //              ThreadStart mt,i,Ellipsoid2D_T(cw,zw,xw,yw,(i*npt/nthreads),((i+1)*npt/nthreads-1)) 
    164 //      endfor 
    165 // 
    166 //      do 
    167 //              variable tgs= ThreadGroupWait(mt,100) 
    168 //      while( tgs != 0 ) 
    169 // 
    170 //      variable dummy= ThreadGroupRelease(mt) 
    171 //       
    172 //      return(0) 
    173 //End 
     102 
     103// - sets up a dependency to a wrapper, not the actual SmearedModelFunction 
     104Proc PlotSmearedEllipsoid2D(str)                                                 
     105        String str 
     106        Prompt str,"Pick the data folder containing the 2D data",popup,getAList(4) 
     107 
     108        if (!exists("Ellipsoid_2DX")) 
     109                Abort "You must have the SANSAnalysis XOP installed to use 2D models" 
     110        endif 
     111         
     112        SetDataFolder $("root:"+str) 
     113         
     114        // Setup parameter table for model function 
     115        make/O/T/N=12 smear_parameters_Ellip2D 
     116        Make/O/D/N=12 smear_coef_Ellip2D 
     117         
     118        smear_coef_Ellip2D[0] = 1.0 
     119        smear_coef_Ellip2D[1] = 20.0 
     120        smear_coef_Ellip2D[2] = 60.0 
     121        smear_coef_Ellip2D[3] = 1.0e-6 
     122        smear_coef_Ellip2D[4] = 6.3e-6 
     123        smear_coef_Ellip2D[5] = 0.0 
     124        smear_coef_Ellip2D[6] = 1.57 
     125        smear_coef_Ellip2D[7] = 0.0 
     126        smear_coef_Ellip2D[8] = 0.0 
     127        smear_coef_Ellip2D[9] = 0.0 
     128        smear_coef_Ellip2D[10] = 0.0 
     129        smear_coef_Ellip2D[11] = 0.0 
     130        // hard-wire the number of integration points 
     131//      smear_coef_Ellip2D[12] = 10 
     132         
     133        smear_parameters_Ellip2D[0] = "Scale" 
     134        smear_parameters_Ellip2D[1] = "Radius_a (rotation axis)" 
     135        smear_parameters_Ellip2D[2] = "Radius_b" 
     136        smear_parameters_Ellip2D[3] = "SLD cylinder (A^-2)" 
     137        smear_parameters_Ellip2D[4] = "SLD solvent" 
     138        smear_parameters_Ellip2D[5] = "Background" 
     139        smear_parameters_Ellip2D[6] = "Axis Theta" 
     140        smear_parameters_Ellip2D[7] = "Axis Phi" 
     141        smear_parameters_Ellip2D[8] = "Sigma of polydisp in R_a [Angstrom]" 
     142        smear_parameters_Ellip2D[9] = "Sigma of polydisp in R_b [Angstrom]" 
     143        smear_parameters_Ellip2D[10] = "Sigma of polydisp in Theta [rad]" 
     144        smear_parameters_Ellip2D[11] = "Sigma of polydisp in Phi [rad]" 
     145         
     146//      smear_parameters_Ellip2D[12] = "Num of polydisp points" 
     147 
     148         
     149        Edit smear_parameters_Ellip2D,smear_coef_Ellip2D                                         
     150         
     151        // generate the triplet representation 
     152        Duplicate/O $(str+"_qx") smeared_Ellip2D 
     153        SetScale d,0,0,"1/cm",smeared_Ellip2D                                    
     154                 
     155        Variable/G gs_Ellip2D=0 
     156        gs_Ellip2D := fSmearedEllipsoid2D(smear_coef_Ellip2D,smeared_Ellip2D)           //wrapper to fill the STRUCT 
     157         
     158        Display $(str+"_qy") vs $(str+"_qx") 
     159        modifygraph log=0 
     160        ModifyGraph mode=3,marker=16,zColor($(str+"_qy"))={smeared_Ellip2D,*,*,YellowHot,0} 
     161        ModifyGraph standoff=0 
     162        ModifyGraph width={Aspect,1} 
     163        ModifyGraph lowTrip=0.001 
     164        Label bottom "qx (A\\S-1\\M)" 
     165        Label left "qy (A\\S-1\\M)" 
     166        AutoPositionWindow/M=1/R=$(WinName(0,1)) $WinName(0,2) 
     167         
     168        // generate the matrix representation 
     169        Duplicate/O $(str+"_qx"), sm_qx 
     170        Duplicate/O $(str+"_qy"), sm_qy         // I can't use local variables in dependencies, so I need the name (that I can't get) 
     171         
     172        // generate the matrix representation 
     173        ConvertQxQy2Mat(sm_qx,sm_qy,smeared_Ellip2D,"sm_Ellip2D_mat") 
     174        Duplicate/O $"sm_Ellip2D_mat",$"sm_Ellip2D_lin"                 //keep a linear-scaled version of the data 
     175        // _mat is for display, _lin is the real calculation 
     176 
     177        // not a function evaluation - this simply keeps the matrix for display in sync with the triplet calculation 
     178        Variable/G gs_Ellip2Dmat=0 
     179        gs_Ellip2Dmat := UpdateQxQy2Mat(sm_qx,sm_qy,smeared_Ellip2D,sm_Ellip2D_lin,sm_Ellip2D_mat) 
     180         
     181        SetDataFolder root: 
     182        AddModelToStrings("SmearedEllipsoid2D","smear_coef_Ellip2D","smear_parameters_Ellip2D","Ellip2D") 
     183End 
     184 
     185 
    174186 
    175187// 
     
    187199 
    188200        Make/O/D/N=13 Ellip2D_tmp 
    189         Ellip2D_tmp = cw 
     201        Ellip2D_tmp[0,11] = cw 
    190202        Ellip2D_tmp[12] = 25 
    191203        Ellip2D_tmp[5] = 0              //pass in a zero background and add it in later 
     
    197209        return(0) 
    198210End 
     211 
     212 
     213/////////////////////smeared functions ////////////////////// 
     214 
     215Function SmearedEllipsoid2D(s) 
     216        Struct ResSmear_2D_AAOStruct &s 
     217         
     218//// non-threaded, but generic calculation 
     219//// the last param is nord      
     220//      Smear_2DModel_PP(Ellipsoid2D_noThread,s,10) 
     221 
     222 
     223//// the last param is nord 
     224        SmearedEllipsoid2D_THR(s,10)             
     225 
     226        return(0) 
     227end 
     228 
     229// for the plot dependency only 
     230Function fSmearedEllipsoid2D(coefW,resultW) 
     231        Wave coefW,resultW 
     232         
     233        String str = getWavesDataFolder(resultW,0) 
     234        String DF="root:"+str+":" 
     235         
     236        WAVE qx = $(DF+str+"_qx") 
     237        WAVE qy = $(DF+str+"_qy") 
     238        WAVE qz = $(DF+str+"_qz") 
     239        WAVE sQpl = $(DF+str+"_sQpl") 
     240        WAVE sQpp = $(DF+str+"_sQpp") 
     241        WAVE shad = $(DF+str+"_fs") 
     242         
     243        STRUCT ResSmear_2D_AAOStruct s 
     244        WAVE s.coefW = coefW     
     245        WAVE s.zw = resultW      
     246        WAVE s.xw[0] = qx 
     247        WAVE s.xw[1] = qy 
     248        WAVE s.qz = qz 
     249        WAVE s.sQpl = sQpl 
     250        WAVE s.sQpp = sQpp 
     251        WAVE s.fs = shad 
     252         
     253        Variable err 
     254        err = SmearedEllipsoid2D(s) 
     255         
     256        return (0) 
     257End 
     258 
     259// 
     260// NON-THREADED IMPLEMENTATION 
     261// -- same as threaded, but no MultiThread KW 
     262// 
     263ThreadSafe Function Ellipsoid2D_noThread(cw,zw,xw,yw) 
     264        Wave cw,zw,xw,yw 
     265         
     266        //      Variable t1=StopMSTimer(-2) 
     267 
     268#if exists("Ellipsoid_2DX")                     //to hide the function if XOP not installed 
     269 
     270        Make/O/D/N=13 Ellip2D_tmp 
     271        Ellip2D_tmp[0,11] = cw 
     272        Ellip2D_tmp[12] = 25 
     273        Ellip2D_tmp[5] = 0              //pass in a zero background and add it in later 
     274         
     275        zw= Ellipsoid_2DX(Ellip2D_tmp,xw,yw) + cw[5] 
     276         
     277#endif 
     278 
     279//      Print "elapsed time = ",(StopMSTimer(-2) - t1)/1e6 
     280         
     281        return(0) 
     282End 
     283 
     284// 
     285// this is the threaded version, that dispatches the calculation out to threads 
     286// 
     287// must be written specific to each 2D function 
     288//  
     289Function SmearedEllipsoid2D_THR(s,nord) 
     290        Struct ResSmear_2D_AAOStruct &s 
     291        Variable nord 
     292         
     293        String weightStr,zStr 
     294         
     295// create all of the necessary quadrature waves here - rather than inside a threadsafe function 
     296        switch(nord)     
     297                case 5:          
     298                        weightStr="gauss5wt" 
     299                        zStr="gauss5z" 
     300                        if (WaveExists($weightStr) == 0) 
     301                                Make/O/D/N=(nord) $weightStr,$zStr 
     302                                Make5GaussPoints($weightStr,$zStr)       
     303                        endif 
     304                        break                            
     305                case 10:                 
     306                        weightStr="gauss10wt" 
     307                        zStr="gauss10z" 
     308                        if (WaveExists($weightStr) == 0) 
     309                                Make/O/D/N=(nord) $weightStr,$zStr 
     310                                Make10GaussPoints($weightStr,$zStr)      
     311                        endif 
     312                        break                            
     313                case 20:                 
     314                        weightStr="gauss20wt" 
     315                        zStr="gauss20z" 
     316                        if (WaveExists($weightStr) == 0) 
     317                                Make/O/D/N=(nord) $weightStr,$zStr 
     318                                Make20GaussPoints($weightStr,$zStr)      
     319                        endif 
     320                        break 
     321                default:                                                         
     322                        Abort "Smear_2DModel_PP_Threaded called with invalid nord value"                                         
     323        endswitch 
     324         
     325        Wave/Z wt = $weightStr 
     326        Wave/Z xi = $zStr               // wave references to pass 
     327 
     328        Variable npt=numpnts(s.xw[0]) 
     329        Variable i,nthreads= ThreadProcessorCount 
     330        variable mt= ThreadGroupCreate(nthreads) 
     331 
     332        Variable t1=StopMSTimer(-2) 
     333         
     334        for(i=0;i<nthreads;i+=1) 
     335//              Print (i*npt/nthreads),((i+1)*npt/nthreads-1) 
     336                ThreadStart mt,i,SmearedEllipsoid2D_T(s.coefW,s.xw[0],s.xw[1],s.qz,s.sQpl,s.sQpp,s.fs,s.zw,wt,xi,(i*npt/nthreads),((i+1)*npt/nthreads-1),nord) 
     337        endfor 
     338 
     339        do 
     340                variable tgs= ThreadGroupWait(mt,100) 
     341        while( tgs != 0 ) 
     342 
     343        variable dummy= ThreadGroupRelease(mt) 
     344         
     345// comment out the threading + uncomment this for testing to make sure that the single thread works 
     346//      nThreads=1 
     347//      SmearedEllipsoid2D_T(s.coefW,s.xw[0],s.xw[1],s.qz,s.sQpl,s.sQpp,s.fs,s.zw,wt,xi,(i*npt/nthreads),((i+1)*npt/nthreads-1),nord) 
     348 
     349        Print "elapsed time = ",(StopMSTimer(-2) - t1)/1e6 
     350         
     351        return(0) 
     352end 
     353 
     354// 
     355// - worker function for threads of Sphere2D 
     356// 
     357ThreadSafe Function SmearedEllipsoid2D_T(coef,qxw,qyw,qzw,sxw,syw,fsw,zw,wt,xi,pt1,pt2,nord) 
     358        WAVE coef,qxw,qyw,qzw,sxw,syw,fsw,zw,wt,xi 
     359        Variable pt1,pt2,nord 
     360         
     361// now passed in.... 
     362//      Wave wt = $weightStr 
     363//      Wave xi = $zStr          
     364 
     365        Variable ii,jj,kk,num 
     366        Variable qx,qy,qz,qval,sx,sy,fs 
     367        Variable qy_pt,qx_pt,res_x,res_y,answer,sumIn,sumOut 
     368         
     369        Variable normFactor,phi,theta,maxSig,numStdDev=3 
     370         
     371/// keep these waves local 
     372        Make/O/D/N=(nord) fcnRet,xptW,res_tot,yptW 
     373         
     374        // now just loop over the points as specified 
     375         
     376        answer=0 
     377         
     378        Variable spl,spp,apl,app,bpl,bpp,phi_pt,qpl_pt 
     379 
     380        //loop over q-values 
     381        for(ii=pt1;ii<(pt2+1);ii+=1) 
     382                 
     383                qx = qxw[ii] 
     384                qy = qyw[ii] 
     385                qz = qzw[ii] 
     386                qval = sqrt(qx^2+qy^2+qz^2) 
     387                spl = sxw[ii] 
     388                spp = syw[ii] 
     389                fs = fsw[ii] 
     390                 
     391                normFactor = 2*pi*spl*spp 
     392                 
     393                phi = FindPhi(qx,qy) 
     394                 
     395                apl = -numStdDev*spl + qval             //parallel = q integration limits 
     396                bpl = numStdDev*spl + qval 
     397                app = -numStdDev*spp + phi              //perpendicular = phi integration limits 
     398                bpp = numStdDev*spp + phi 
     399                 
     400                //make sure the limits are reasonable. 
     401                if(apl < 0) 
     402                        apl = 0 
     403                endif 
     404                // do I need to specially handle limits when phi ~ 0? 
     405         
     406                 
     407                sumOut = 0 
     408                for(jj=0;jj<nord;jj+=1)         // call phi the "outer' 
     409                        phi_pt = (xi[jj]*(bpp-app)+app+bpp)/2 
     410 
     411                        sumIn=0 
     412                        for(kk=0;kk<nord;kk+=1)         //at phi, integrate over Qpl 
     413 
     414                                qpl_pt = (xi[kk]*(bpl-apl)+apl+bpl)/2 
     415                                 
     416                                FindQxQy(qpl_pt,phi_pt,qx_pt,qy_pt)             //find the corresponding QxQy to the Q,phi 
     417                                yPtw[kk] = qy_pt                                        //phi is the same in this loop, but qy is not 
     418                                xPtW[kk] = qx_pt                                        //qx is different here too, as we're varying Qpl 
     419                                 
     420                                res_tot[kk] = exp(-0.5*( (qpl_pt-qval)^2/spl/spl + (phi_pt-phi)^2/spp/spp ) ) 
     421                                res_tot[kk] /= normFactor 
     422//                              res_tot[kk] *= fs 
     423 
     424                        endfor 
     425                         
     426                        Ellipsoid2D_noThread(coef,fcnRet,xptw,yptw)                     //fcn passed in is an AAO 
     427                         
     428                        //sumIn += wt[jj]*wt[kk]*res_tot*fcnRet[0] 
     429                        fcnRet *= wt[jj]*wt*res_tot 
     430                        // 
     431                        answer += (bpl-apl)/2.0*sum(fcnRet)             // 
     432                endfor 
     433 
     434                answer *= (bpp-app)/2.0 
     435                zw[ii] = answer 
     436        endfor 
     437         
     438        return(0) 
     439end 
     440 
     441 
  • sans/Dev/trunk/NCNR_User_Procedures/Analysis/Models/Models_2D/EllipticalCylinder2D_v40.ipf

    r753 r771  
    101101End 
    102102 
    103 //AAO version, uses XOP if available 
    104 // simply calls the original single point calculation with 
    105 // a wave assignment (this will behave nicely if given point ranges) 
    106 // 
    107 // NON-THREADED IMPLEMENTATION 
    108 // 
    109 //Function EllipticalCylinder2D(cw,zw,xw,yw) : FitFunc 
    110 //      Wave cw,zw,xw,yw 
    111 //       
    112 //#if exists("EllipticalCylinderModel_D") 
    113 // 
    114 //      Make/O/D/N=15 EllCyl2D_tmp 
    115 //      EllCyl2D_tmp = cw 
    116 //      EllCyl2D_tmp[14] = 25 
    117 //       
    118 //      zw = EllipticalCylinderModel_D(EllCyl2D_tmp,xw,yw) 
    119 //       
    120 ////    zw = EllipticalCylinderModel_D(cw,xw,yw) 
    121 //#else 
    122 //      Abort "You do not have the SANS Analysis XOP installed" 
    123 //#endif 
    124 //      return(0) 
    125 //End 
    126 // 
    127  
    128 ////threaded version of the function 
    129 //ThreadSafe Function EllipticalCylinder2D_T(cw,zw,xw,yw,p1,p2) 
    130 //      WAVE cw,zw,xw,yw 
    131 //      Variable p1,p2 
    132 //       
    133 //#if exists("EllipticalCylinder_2DX")                  //to hide the function if XOP not installed 
    134 // 
    135 //      Make/O/D/N=15 EllCyl2D_tmp 
    136 //      EllCyl2D_tmp = cw 
    137 //      EllCyl2D_tmp[14] = 25 
    138 //      EllCyl2D_tmp[6] = 0             //pass in a zero background and add it in later 
    139 //       
    140 //      zw[p1,p2]= EllipticalCylinder_2DX(EllCyl2D_tmp,xw,yw) + cw[6] 
    141 //       
    142 //#endif 
    143 // 
    144 //      return 0 
    145 //End 
    146 // 
    147 //// 
    148 ////  Fit function that is actually a wrapper to dispatch the calculation to N threads 
    149 //// 
    150 //// nthreads is 1 or an even number, typically 2 
    151 //// it doesn't matter if npt is odd. In this case, fractional point numbers are passed 
    152 //// and the wave indexing works just fine - I tested this with test waves of 7 and 8 points 
    153 //// and the points "2.5" and "3.5" evaluate correctly as 2 and 3 
    154 //// 
    155 //Function EllipticalCylinder2D(cw,zw,xw,yw) : FitFunc 
    156 //      Wave cw,zw,xw,yw 
    157 //       
    158 //      Variable npt=numpnts(yw) 
    159 //      Variable i,nthreads= ThreadProcessorCount 
    160 //      variable mt= ThreadGroupCreate(nthreads) 
    161 // 
    162 //      for(i=0;i<nthreads;i+=1) 
    163 //      //      Print (i*npt/nthreads),((i+1)*npt/nthreads-1) 
    164 //              ThreadStart mt,i,EllipticalCylinder2D_T(cw,zw,xw,yw,(i*npt/nthreads),((i+1)*npt/nthreads-1)) 
    165 //      endfor 
    166 // 
    167 //      do 
    168 //              variable tgs= ThreadGroupWait(mt,100) 
    169 //      while( tgs != 0 ) 
    170 // 
    171 //      variable dummy= ThreadGroupRelease(mt) 
    172 //       
    173 //      return(0) 
    174 //End 
     103 
     104Proc PlotSmearedEllipticalCyl2D(str)                                             
     105        String str 
     106        Prompt str,"Pick the data folder containing the 2D data",popup,getAList(4) 
     107 
     108        if (!exists("EllipticalCylinder_2DX")) 
     109                Abort "You must have the SANSAnalysis XOP installed to use 2D models" 
     110        endif 
     111         
     112        SetDataFolder $("root:"+str) 
     113         
     114        // Setup parameter table for model function 
     115        make/O/T/N=14 smear_parameters_EllCyl2D 
     116        Make/O/D/N=14 smear_coef_EllCyl2D 
     117         
     118        smear_coef_EllCyl2D[0] = 1.0 
     119        smear_coef_EllCyl2D[1] = 20.0 
     120        smear_coef_EllCyl2D[2] = 1.5 
     121        smear_coef_EllCyl2D[3] = 400.0 
     122        smear_coef_EllCyl2D[4] = 3e-6 
     123        smear_coef_EllCyl2D[5] = 6.3e-6 
     124        smear_coef_EllCyl2D[6] = 0.0 
     125        smear_coef_EllCyl2D[7] = 1.57 
     126        smear_coef_EllCyl2D[8] = 0.0 
     127        smear_coef_EllCyl2D[9] = 0.0 
     128        smear_coef_EllCyl2D[10] = 0.0 
     129        smear_coef_EllCyl2D[11] = 0.0 
     130        smear_coef_EllCyl2D[12] = 0.0 
     131        smear_coef_EllCyl2D[13] = 0.0 
     132         
     133        // now hard-wire the # of integration points 
     134        //smear_coef_EllCyl2D[14] = 25 
     135                 
     136        smear_parameters_EllCyl2D[0] = "Scale" 
     137        smear_parameters_EllCyl2D[1] = "R_minor" 
     138        smear_parameters_EllCyl2D[2] = "R_ratio (major/minor)" 
     139        smear_parameters_EllCyl2D[3] = "Length" 
     140        smear_parameters_EllCyl2D[4] = "SLD cylinder (A^-2)" 
     141        smear_parameters_EllCyl2D[5] = "SLD solvent" 
     142        smear_parameters_EllCyl2D[6] = "Background" 
     143        smear_parameters_EllCyl2D[7] = "Axis Theta" 
     144        smear_parameters_EllCyl2D[8] = "Axis Phi" 
     145        smear_parameters_EllCyl2D[9] = "Ellipse Psi" 
     146        smear_parameters_EllCyl2D[10] = "Sigma of polydisp in R_minor [Angstrom]" 
     147        smear_parameters_EllCyl2D[11] = "Sigma of polydisp in R_ratio" 
     148        smear_parameters_EllCyl2D[12] = "Sigma of polydisp in Theta [rad]" 
     149        smear_parameters_EllCyl2D[13] = "Sigma of polydisp in Phi [rad]" 
     150        //smear_parameters_EllCyl2D[14] = "Num of polydisp points" 
     151         
     152        Edit smear_parameters_EllCyl2D,smear_coef_EllCyl2D                                       
     153         
     154        // generate the triplet representation 
     155        Duplicate/O $(str+"_qx") smeared_EllCyl2D 
     156        SetScale d,0,0,"1/cm",smeared_EllCyl2D                                   
     157                 
     158        Variable/G gs_EllCyl2D=0 
     159        gs_EllCyl2D := fSmearedEllipticalCyl2D(smear_coef_EllCyl2D,smeared_EllCyl2D)            //wrapper to fill the STRUCT 
     160         
     161        Display $(str+"_qy") vs $(str+"_qx") 
     162        modifygraph log=0 
     163        ModifyGraph mode=3,marker=16,zColor($(str+"_qy"))={smeared_EllCyl2D,*,*,YellowHot,0} 
     164        ModifyGraph standoff=0 
     165        ModifyGraph width={Aspect,1} 
     166        ModifyGraph lowTrip=0.001 
     167        Label bottom "qx (A\\S-1\\M)" 
     168        Label left "qy (A\\S-1\\M)" 
     169        AutoPositionWindow/M=1/R=$(WinName(0,1)) $WinName(0,2) 
     170         
     171        // generate the matrix representation 
     172        Duplicate/O $(str+"_qx"), sm_qx 
     173        Duplicate/O $(str+"_qy"), sm_qy         // I can't use local variables in dependencies, so I need the name (that I can't get) 
     174         
     175        // generate the matrix representation 
     176        ConvertQxQy2Mat(sm_qx,sm_qy,smeared_EllCyl2D,"sm_EllCyl2D_mat") 
     177        Duplicate/O $"sm_EllCyl2D_mat",$"sm_EllCyl2D_lin"               //keep a linear-scaled version of the data 
     178        // _mat is for display, _lin is the real calculation 
     179 
     180        // not a function evaluation - this simply keeps the matrix for display in sync with the triplet calculation 
     181        Variable/G gs_EllCyl2Dmat=0 
     182        gs_EllCyl2Dmat := UpdateQxQy2Mat(sm_qx,sm_qy,smeared_EllCyl2D,sm_EllCyl2D_lin,sm_EllCyl2D_mat) 
     183         
     184         
     185        SetDataFolder root: 
     186        AddModelToStrings("SmearedEllipticalCyl2D","smear_coef_EllCyl2D","smear_parameters_EllCyl2D","EllCyl2D") 
     187End 
     188 
    175189 
    176190// 
     
    188202 
    189203        Make/O/D/N=15 EllCyl2D_tmp 
    190         EllCyl2D_tmp = cw 
     204        EllCyl2D_tmp[0,13] = cw 
    191205        EllCyl2D_tmp[14] = 25 
    192206        EllCyl2D_tmp[6] = 0             //pass in a zero background and add it in later 
     
    198212        return(0) 
    199213End 
     214 
     215 
     216/////////////////////smeared functions ////////////////////// 
     217 
     218Function SmearedEllipticalCyl2D(s) 
     219        Struct ResSmear_2D_AAOStruct &s 
     220         
     221//// non-threaded, but generic calculation 
     222//// the last param is nord      
     223//      Smear_2DModel_PP(EllipticalCylinder2D_noThread,s,10) 
     224 
     225 
     226//// the last param is nord 
     227        SmearedEllipticalCyl2D_THR(s,10)                 
     228 
     229        return(0) 
     230end 
     231 
     232// for the plot dependency only 
     233Function fSmearedEllipticalCyl2D(coefW,resultW) 
     234        Wave coefW,resultW 
     235         
     236        String str = getWavesDataFolder(resultW,0) 
     237        String DF="root:"+str+":" 
     238         
     239        WAVE qx = $(DF+str+"_qx") 
     240        WAVE qy = $(DF+str+"_qy") 
     241        WAVE qz = $(DF+str+"_qz") 
     242        WAVE sQpl = $(DF+str+"_sQpl") 
     243        WAVE sQpp = $(DF+str+"_sQpp") 
     244        WAVE shad = $(DF+str+"_fs") 
     245         
     246        STRUCT ResSmear_2D_AAOStruct s 
     247        WAVE s.coefW = coefW     
     248        WAVE s.zw = resultW      
     249        WAVE s.xw[0] = qx 
     250        WAVE s.xw[1] = qy 
     251        WAVE s.qz = qz 
     252        WAVE s.sQpl = sQpl 
     253        WAVE s.sQpp = sQpp 
     254        WAVE s.fs = shad 
     255         
     256        Variable err 
     257        err = SmearedEllipticalCyl2D(s) 
     258         
     259        return (0) 
     260End 
     261 
     262// 
     263// NON-THREADED IMPLEMENTATION 
     264// -- same as threaded, but no MultiThread KW 
     265// 
     266ThreadSafe Function EllipticalCylinder2D_noThread(cw,zw,xw,yw) 
     267        Wave cw,zw,xw,yw 
     268         
     269        //      Variable t1=StopMSTimer(-2) 
     270 
     271#if exists("EllipticalCylinder_2DX")                    //to hide the function if XOP not installed 
     272 
     273        Make/O/D/N=15 EllCyl2D_tmp 
     274        EllCyl2D_tmp[0,13] = cw 
     275        EllCyl2D_tmp[14] = 25 
     276        EllCyl2D_tmp[6] = 0             //pass in a zero background and add it in later 
     277         
     278        zw= EllipticalCylinder_2DX(EllCyl2D_tmp,xw,yw) + cw[6] 
     279         
     280#endif 
     281 
     282//      Print "elapsed time = ",(StopMSTimer(-2) - t1)/1e6 
     283         
     284        return(0) 
     285End 
     286 
     287// 
     288// this is the threaded version, that dispatches the calculation out to threads 
     289// 
     290// must be written specific to each 2D function 
     291//  
     292Function SmearedEllipticalCyl2D_THR(s,nord) 
     293        Struct ResSmear_2D_AAOStruct &s 
     294        Variable nord 
     295         
     296        String weightStr,zStr 
     297         
     298// create all of the necessary quadrature waves here - rather than inside a threadsafe function 
     299        switch(nord)     
     300                case 5:          
     301                        weightStr="gauss5wt" 
     302                        zStr="gauss5z" 
     303                        if (WaveExists($weightStr) == 0) 
     304                                Make/O/D/N=(nord) $weightStr,$zStr 
     305                                Make5GaussPoints($weightStr,$zStr)       
     306                        endif 
     307                        break                            
     308                case 10:                 
     309                        weightStr="gauss10wt" 
     310                        zStr="gauss10z" 
     311                        if (WaveExists($weightStr) == 0) 
     312                                Make/O/D/N=(nord) $weightStr,$zStr 
     313                                Make10GaussPoints($weightStr,$zStr)      
     314                        endif 
     315                        break                            
     316                case 20:                 
     317                        weightStr="gauss20wt" 
     318                        zStr="gauss20z" 
     319                        if (WaveExists($weightStr) == 0) 
     320                                Make/O/D/N=(nord) $weightStr,$zStr 
     321                                Make20GaussPoints($weightStr,$zStr)      
     322                        endif 
     323                        break 
     324                default:                                                         
     325                        Abort "Smear_2DModel_PP_Threaded called with invalid nord value"                                         
     326        endswitch 
     327         
     328        Wave/Z wt = $weightStr 
     329        Wave/Z xi = $zStr               // wave references to pass 
     330 
     331        Variable npt=numpnts(s.xw[0]) 
     332        Variable i,nthreads= ThreadProcessorCount 
     333        variable mt= ThreadGroupCreate(nthreads) 
     334 
     335        Variable t1=StopMSTimer(-2) 
     336         
     337        for(i=0;i<nthreads;i+=1) 
     338//              Print (i*npt/nthreads),((i+1)*npt/nthreads-1) 
     339                ThreadStart mt,i,SmearedEllipticalCyl2D_T(s.coefW,s.xw[0],s.xw[1],s.qz,s.sQpl,s.sQpp,s.fs,s.zw,wt,xi,(i*npt/nthreads),((i+1)*npt/nthreads-1),nord) 
     340        endfor 
     341 
     342        do 
     343                variable tgs= ThreadGroupWait(mt,100) 
     344        while( tgs != 0 ) 
     345 
     346        variable dummy= ThreadGroupRelease(mt) 
     347         
     348// comment out the threading + uncomment this for testing to make sure that the single thread works 
     349//      nThreads=1 
     350//      SmearedEllipticalCyl2D_T(s.coefW,s.xw[0],s.xw[1],s.qz,s.sQpl,s.sQpp,s.fs,s.zw,wt,xi,(i*npt/nthreads),((i+1)*npt/nthreads-1),nord) 
     351 
     352        Print "elapsed time = ",(StopMSTimer(-2) - t1)/1e6 
     353         
     354        return(0) 
     355end 
     356 
     357// 
     358// - worker function for threads of Sphere2D 
     359// 
     360ThreadSafe Function SmearedEllipticalCyl2D_T(coef,qxw,qyw,qzw,sxw,syw,fsw,zw,wt,xi,pt1,pt2,nord) 
     361        WAVE coef,qxw,qyw,qzw,sxw,syw,fsw,zw,wt,xi 
     362        Variable pt1,pt2,nord 
     363         
     364// now passed in.... 
     365//      Wave wt = $weightStr 
     366//      Wave xi = $zStr          
     367 
     368        Variable ii,jj,kk,num 
     369        Variable qx,qy,qz,qval,sx,sy,fs 
     370        Variable qy_pt,qx_pt,res_x,res_y,answer,sumIn,sumOut 
     371         
     372        Variable normFactor,phi,theta,maxSig,numStdDev=3 
     373         
     374/// keep these waves local 
     375        Make/O/D/N=(nord) fcnRet,xptW,res_tot,yptW 
     376         
     377        // now just loop over the points as specified 
     378         
     379        answer=0 
     380         
     381        Variable spl,spp,apl,app,bpl,bpp,phi_pt,qpl_pt 
     382 
     383        //loop over q-values 
     384        for(ii=pt1;ii<(pt2+1);ii+=1) 
     385                 
     386                qx = qxw[ii] 
     387                qy = qyw[ii] 
     388                qz = qzw[ii] 
     389                qval = sqrt(qx^2+qy^2+qz^2) 
     390                spl = sxw[ii] 
     391                spp = syw[ii] 
     392                fs = fsw[ii] 
     393                 
     394                normFactor = 2*pi*spl*spp 
     395                 
     396                phi = FindPhi(qx,qy) 
     397                 
     398                apl = -numStdDev*spl + qval             //parallel = q integration limits 
     399                bpl = numStdDev*spl + qval 
     400                app = -numStdDev*spp + phi              //perpendicular = phi integration limits 
     401                bpp = numStdDev*spp + phi 
     402                 
     403                //make sure the limits are reasonable. 
     404                if(apl < 0) 
     405                        apl = 0 
     406                endif 
     407                // do I need to specially handle limits when phi ~ 0? 
     408         
     409                 
     410                sumOut = 0 
     411                for(jj=0;jj<nord;jj+=1)         // call phi the "outer' 
     412                        phi_pt = (xi[jj]*(bpp-app)+app+bpp)/2 
     413 
     414                        sumIn=0 
     415                        for(kk=0;kk<nord;kk+=1)         //at phi, integrate over Qpl 
     416 
     417                                qpl_pt = (xi[kk]*(bpl-apl)+apl+bpl)/2 
     418                                 
     419                                FindQxQy(qpl_pt,phi_pt,qx_pt,qy_pt)             //find the corresponding QxQy to the Q,phi 
     420                                yPtw[kk] = qy_pt                                        //phi is the same in this loop, but qy is not 
     421                                xPtW[kk] = qx_pt                                        //qx is different here too, as we're varying Qpl 
     422                                 
     423                                res_tot[kk] = exp(-0.5*( (qpl_pt-qval)^2/spl/spl + (phi_pt-phi)^2/spp/spp ) ) 
     424                                res_tot[kk] /= normFactor 
     425//                              res_tot[kk] *= fs 
     426 
     427                        endfor 
     428                         
     429                        EllipticalCylinder2D_noThread(coef,fcnRet,xptw,yptw)                    //fcn passed in is an AAO 
     430                         
     431                        //sumIn += wt[jj]*wt[kk]*res_tot*fcnRet[0] 
     432                        fcnRet *= wt[jj]*wt*res_tot 
     433                        // 
     434                        answer += (bpl-apl)/2.0*sum(fcnRet)             // 
     435                endfor 
     436 
     437                answer *= (bpp-app)/2.0 
     438                zw[ii] = answer 
     439        endfor 
     440         
     441        return(0) 
     442end 
     443 
     444 
  • sans/Dev/trunk/NCNR_User_Procedures/Analysis/Models/Models_2D/Sphere_2D_v40.ipf

    r754 r771  
    119119        Wave cw,zw,xw,yw 
    120120 
    121  
    122121#if exists("Sphere_2DX")                        //to hide the function if XOP not installed 
    123          
    124122        MultiThread     zw= Sphere_2DX(cw,xw,yw) 
    125  
    126123#endif 
    127124 
     
    217214 
    218215//// the last param is nord 
    219         SmearSphere2D_THR(s,10)          
     216        SmearedSphere2D_THR(s,10)                
    220217 
    221218        return(0) 
     
    260257// must be written specific to each 2D function 
    261258//  
    262 Function SmearSphere2D_THR(s,nord) 
     259Function SmearedSphere2D_THR(s,nord) 
    263260        Struct ResSmear_2D_AAOStruct &s 
    264261        Variable nord 
     
    307304        for(i=0;i<nthreads;i+=1) 
    308305//              Print (i*npt/nthreads),((i+1)*npt/nthreads-1) 
    309                 ThreadStart mt,i,SmearSphere2D_T(s.coefW,s.xw[0],s.xw[1],s.qz,s.sQpl,s.sQpp,s.fs,s.zw,wt,xi,(i*npt/nthreads),((i+1)*npt/nthreads-1),nord) 
     306                ThreadStart mt,i,SmearedSphere2D_T(s.coefW,s.xw[0],s.xw[1],s.qz,s.sQpl,s.sQpp,s.fs,s.zw,wt,xi,(i*npt/nthreads),((i+1)*npt/nthreads-1),nord) 
    310307        endfor 
    311308 
     
    328325// - worker function for threads of Sphere2D 
    329326// 
    330 ThreadSafe Function SmearSphere2D_T(coef,qxw,qyw,qzw,sxw,syw,fsw,zw,wt,xi,pt1,pt2,nord) 
     327ThreadSafe Function SmearedSphere2D_T(coef,qxw,qyw,qzw,sxw,syw,fsw,zw,wt,xi,pt1,pt2,nord) 
    331328        WAVE coef,qxw,qyw,qzw,sxw,syw,fsw,zw,wt,xi 
    332329        Variable pt1,pt2,nord 
  • sans/Dev/trunk/NCNR_User_Procedures/Analysis/Models/NewModels_2010/PolyRaspberry_v40.ipf

    r770 r771  
    11#pragma rtGlobals=1             // Use modern global access method. 
    22 
     3 
     4//////////////////////////////////////////////////// 
     5// Raspberry model 
     6// 
     7// Larson-Smith et al. Small angle scattering model for Pickering emulsions and raspberry particles. 
     8//   Journal of colloid and interface science (2010) vol. 343 (1) pp. 36-41 
     9// 
     10//////////////////////////////////////////////////// 
    311 
    412// Raspberry particles with polydisperse large sphere 
  • sans/Dev/trunk/NCNR_User_Procedures/Analysis/Models/NewModels_2010/Raspberry_v40.ipf

    r770 r771  
    33 
    44//////////////////////////////////////////////////// 
    5 // Raspberry model - Pozzo & Larson 
     5// Raspberry model 
     6// 
     7// Larson-Smith et al. Small angle scattering model for Pickering emulsions and raspberry particles. 
     8//   Journal of colloid and interface science (2010) vol. 343 (1) pp. 36-41 
     9// 
    610// Default parameters are for a 5000A hexadecane drop stabilized by 100A silica particles 
    711// in D2O. The particles are 50% inserted into the interface (delta = 0) and surface coverage is 50% 
  • sans/Dev/trunk/NCNR_User_Procedures/Analysis/Packages/ModelPicker/SANSModelPicker_v40.ipf

    r749 r771  
    139139////paste here... after deleting the old make statement and list 
    140140         
    141   Make/O/T/N=92  SANS_Model_List 
     141  Make/O/T/N=95  SANS_Model_List 
    142142 
    143143  SANS_Model_List[0] = "BE_Polyelectrolyte.ipf" 
     
    236236  SANS_Model_List[90] = "Guinier_Porod.ipf" 
    237237  SANS_Model_List[91] = "PolymerExcludVol.ipf" 
     238  SANS_Model_List[92] = "Raspberry.ipf" 
     239  SANS_Model_List[93] = "PolyRaspberry.ipf" 
     240  SANS_Model_List[94] = "Yukawa_SQ.ipf" 
    238241 
    239242  ///end paste here 
  • sans/Dev/trunk/NCNR_User_Procedures/Analysis/Packages/ModelPicker/SANS_Model_List.itx

    r515 r771  
    7070        "Core_and_NShells.ipf" 
    7171        "PolyCore_and_NShells.ipf" 
    72         "Fractal_Polysphere.ipf" 
     72        "Fractal_PolySphere.ipf" 
    7373        "GaussLorentzGel.ipf" 
    7474        "PolyGaussCoil.ipf" 
     
    9090        "CSParallelepiped.ipf" 
    9191        "Fractal_PolyCore.ipf" 
     92        "FuzzySpheres.ipf" 
     93        "FuzzySpheres_Sq.ipf" 
     94        "Guinier_Porod.ipf" 
     95        "PolymerExcludVol.ipf" 
     96        "Raspberry.ipf" 
     97        "PolyRaspberry.ipf" 
     98        "Yukawa_SQ.ipf" 
    9299END 
    93100X SetScale/P x 0,1,"", SANS_Model_List; SetScale y 0,0,"", SANS_Model_List 
Note: See TracChangeset for help on using the changeset viewer.