- Timestamp:
- May 18, 2009 5:15:08 PM (14 years ago)
- Location:
- sans/Dev/trunk/NCNR_User_Procedures
- Files:
-
- 6 edited
-
Analysis/Models/Models_2D/Sphere_2D_v40.ipf (modified) (3 diffs)
-
Common/GaussUtils_v40.ipf (modified) (2 diffs)
-
Common/Packages/PlotManager/PlotUtils2D_v40.ipf (modified) (8 diffs)
-
Reduction/SANS/NCNR_Utils.ipf (modified) (1 diff)
-
Reduction/SANS/ProtocolAsPanel.ipf (modified) (1 diff)
-
Reduction/SANS/WriteQIS.ipf (modified) (5 diffs)
Legend:
- Unmodified
- Added
- Removed
-
sans/Dev/trunk/NCNR_User_Procedures/Analysis/Models/Models_2D/Sphere_2D_v40.ipf
r253 r509 32 32 Duplicate/O $(str+"_qy") ywave_sf2D,zwave_sf2D 33 33 34 Variable/G g s_sf2D=035 g s_sf2D := Sphere2D(coef_sf2D,zwave_sf2D,xwave_sf2D,ywave_sf2D) //AAO 2D calculation34 Variable/G g_sf2D=0 35 g_sf2D := Sphere2D(coef_sf2D,zwave_sf2D,xwave_sf2D,ywave_sf2D) //AAO 2D calculation 36 36 37 37 Display ywave_sf2D vs xwave_sf2D … … 51 51 52 52 // not a function evaluation - this simply keeps the matrix for display in sync with the triplet calculation 53 Variable/G g s_sf2Dmat=054 g s_sf2Dmat := UpdateQxQy2Mat(xwave_sf2D,ywave_sf2D,zwave_sf2D,sf2D_lin,sf2D_mat)53 Variable/G g_sf2Dmat=0 54 g_sf2Dmat := UpdateQxQy2Mat(xwave_sf2D,ywave_sf2D,zwave_sf2D,sf2D_lin,sf2D_mat) 55 55 56 56 57 57 SetDataFolder root: 58 58 AddModelToStrings("Sphere2D","coef_sf2D","sf2D") 59 End 60 61 // - sets up a dependency to a wrapper, not the actual SmearedModelFunction 62 Proc PlotSmearedSphere2D(str) 63 String str 64 Prompt str,"Pick the data folder containing the 2D data",popup,getAList(4) 65 66 // if any of the resolution waves are missing => abort 67 // if(ResolutionWavesMissingDF(str)) //updated to NOT use global strings (in GaussUtils) 68 // Abort 69 // endif 70 71 SetDataFolder $("root:"+str) 72 73 // Setup parameter table for model function 74 Make/O/D smear_coef_sf2D = {1.,60,1e-6,6.3e-6,0.01} 75 make/o/t smear_parameters_sf2D = {"scale","Radius (A)","SLD sphere (A-2)","SLD solvent (A-2)","bkgd (cm-1)"} 76 Edit smear_parameters_sf2D,smear_coef_sf2D 77 78 Duplicate/O $(str+"_qx") smeared_sf2D //1d place for the smeared model 79 SetScale d,0,0,"1/cm",smeared_sf2D 80 81 Variable/G gs_sf2D=0 82 gs_sf2D := fSmearedSphere2D(smear_coef_sf2D,smeared_sf2D) //this wrapper fills the STRUCT 83 84 Display $(str+"_qy") vs $(str+"_qx") 85 modifygraph log=0 86 ModifyGraph mode=3,marker=16,zColor($(str+"_qy"))={smeared_sf2D,*,*,YellowHot,0} 87 ModifyGraph standoff=0 88 ModifyGraph width={Aspect,1} 89 ModifyGraph lowTrip=0.001 90 Label bottom "qx (A\\S-1\\M)" 91 Label left "qy (A\\S-1\\M)" 92 AutoPositionWindow/M=1/R=$(WinName(0,1)) $WinName(0,2) 93 94 SetDataFolder root: 95 AddModelToStrings("SmearedSphere2D","smear_coef_sf2D","sf2D") 59 96 End 60 97 … … 106 143 return(0) 107 144 End 145 146 147 //threaded version of the function 148 Function Sphere2D_noThread(cw,zw,xw,yw) 149 WAVE cw,zw, xw,yw 150 151 #if exists("Sphere_2DX") //to hide the function if XOP not installed 152 153 zw= Sphere_2DX(cw,xw,yw) 154 155 #endif 156 157 return 0 158 End 159 160 161 Function SmearedSphere2D(s) 162 Struct ResSmear_2D_AAOStruct &s 163 164 Smear_2DModel_5(Sphere2D_noThread,s) 165 return(0) 166 end 167 168 169 Function fSmearedSphere2D(coefW,resultW) 170 Wave coefW,resultW 171 172 String str = getWavesDataFolder(resultW,0) 173 String DF="root:"+str+":" 174 175 WAVE qx = $(DF+str+"_qx") 176 WAVE qy = $(DF+str+"_qy") 177 WAVE qz = $(DF+str+"_qz") 178 WAVE sigQx = $(DF+str+"_sigQx") 179 WAVE sigQy = $(DF+str+"_sigQy") 180 WAVE shad = $(DF+str+"_fs") 181 182 STRUCT ResSmear_2D_AAOStruct s 183 WAVE s.coefW = coefW 184 WAVE s.zw = resultW 185 WAVE s.qx = qx 186 WAVE s.qy = qy 187 WAVE s.qz = qz 188 WAVE s.sigQx = sigQx 189 WAVE s.sigQy = sigQy 190 WAVE s.fs = shad 191 192 Variable err 193 err = SmearedSphere2D(s) 194 195 return (0) 196 End -
sans/Dev/trunk/NCNR_User_Procedures/Common/GaussUtils_v40.ipf
r480 r509 61 61 EndStructure 62 62 63 64 // tentative pass at 2D resolution smearing 65 // 66 Structure ResSmear_2D_AAOStruct 67 Wave coefW 68 Wave zw //answer 69 Wave qy // q-value 70 Wave qx 71 Wave qz 72 Wave sigQx //resolution 73 Wave sigQy 74 Wave fs 75 String info 76 EndStructure 63 77 64 78 … … 1062 1076 end 1063 1077 1078 // prototype function for 2D smearing routine 1079 Function SANS_2D_ModelAAO_proto(w,zw,xw,yw) 1080 Wave w,zw,xw,yw 1081 1082 Print "in SANSModelAAO_proto function" 1083 return(1) 1084 end 1085 1086 // prototype function for fit wrapper using 2D smearing 1087 // not used (yet) 1088 Function SANS_2D_ModelSTRUCT_proto(s) 1089 Struct ResSmear_2D_AAOStruct &s 1090 1091 Print "in SANSModelSTRUCT_proto function" 1092 return(1) 1093 end 1094 1064 1095 //Numerical Recipes routine to calculate the nn(th) stage 1065 1096 //refinement of a trapezoid integration -
sans/Dev/trunk/NCNR_User_Procedures/Common/Packages/PlotManager/PlotUtils2D_v40.ipf
r490 r509 15 15 16 16 // 17 // changed May 2009 to read in the resolution information (now 7 columns) 18 // -- subject to change -- 17 19 // 18 20 Proc LoadQxQy() … … 24 26 25 27 String w0,w1,w2,n0,n1,n2 28 String w3,w4,w5,w6,n3,n4,n5,n6 26 29 27 30 // put the names of the three loaded waves into local names … … 29 32 n1 = StringFromList(1, S_waveNames ,";" ) 30 33 n2 = StringFromList(2, S_waveNames ,";" ) 34 n3 = StringFromList(3, S_waveNames ,";" ) 35 n4 = StringFromList(4, S_waveNames ,";" ) 36 n5 = StringFromList(5, S_waveNames ,";" ) 37 n6 = StringFromList(6, S_waveNames ,";" ) 31 38 32 39 //remove the semicolon AND period from files from the VAX … … 34 41 w1 = CleanupName((S_fileName + "_qy"),0) 35 42 w2 = CleanupName((S_fileName + "_i"),0) 36 43 w3 = CleanupName((S_fileName + "_qz"),0) 44 w4 = CleanupName((S_fileName + "_sigQx"),0) 45 w5 = CleanupName((S_fileName + "_sigQy"),0) 46 w6 = CleanupName((S_fileName + "_fs"),0) 47 37 48 String baseStr=w1[0,strlen(w1)-4] 38 49 if(DataFolderExists("root:"+baseStr)) … … 40 51 if(V_flag==2) //user selected No, don't load the data 41 52 SetDataFolder root: 42 KillWaves $n0,$n1,$n2 // kill the default waveX that were loaded53 KillWaves $n0,$n1,$n2,$n3,$n4,$n5,$n6 // kill the default waveX that were loaded 43 54 return //quits the macro 44 55 endif … … 69 80 Duplicate/O $("root:"+n1), $w1 70 81 Duplicate/O $("root:"+n2), $w2 82 Duplicate/O $("root:"+n3), $w3 83 Duplicate/O $("root:"+n4), $w4 84 Duplicate/O $("root:"+n5), $w5 85 Duplicate/O $("root:"+n6), $w6 71 86 72 87 Variable/G gIsLogScale = 0 … … 88 103 //clean up 89 104 SetDataFolder root: 90 KillWaves/Z $n0,$n1,$n2 105 KillWaves/Z $n0,$n1,$n2,$n3,$n4,$n5,$n6 91 106 EndMacro 92 107 … … 413 428 //make log scale 414 429 w=log(linW) 430 w = w[p][q] == -inf ? NaN : w[p][q] //remove the -inf for display 415 431 endif 416 432 endfor -
sans/Dev/trunk/NCNR_User_Procedures/Reduction/SANS/NCNR_Utils.ipf
r412 r509 152 152 Return results 153 153 End 154 155 // 156 // 157 // this should end up in NCNR_Utils once it's fully functional 158 // 159 // lots of unnecessary junk... 160 // 161 //********************** 162 // 2D resolution function calculation - in terms of X and Y 163 // 164 // based on notes from David Mildner, 2008 165 // 166 // - 21 MAR 07 uses projected BS diameter on the detector 167 // - APR 07 still need to add resolution with lenses. currently there is no flag in the 168 // raw data header to indicate the presence of lenses. 169 // 170 // - Aug 07 - added input to switch calculation based on lenses (==1 if in) 171 // 172 // passed values are read from RealsRead 173 // except DDet and apOff, which are set from globals before passing 174 // 175 // phi is the azimuthal angle, CCW from +x axis 176 // r_dist is the real-space distance from ctr of detector to QxQy pixel location 177 // 178 Function/S get2DResolution(inQ,phi,lambda,lambdaWidth,DDet,apOff,S1,S2,L1,L2,BS,del_r,usingLenses,r_dist,SigmaQX,SigmaQY,fSubS) 179 Variable inQ, phi,lambda, lambdaWidth, DDet, apOff, S1, S2, L1, L2, BS, del_r,usingLenses,r_dist 180 Variable &SigmaQX,&SigmaQY,&fSubS //these are the output quantities at the input Q value 181 182 //lots of calculation variables 183 Variable a2, lp, v_lambda, v_b, v_d, vz, yg, v_g 184 Variable r0, delta, inc_gamma, fr, fv, rmd, v_r1, rm, v_r 185 186 //Constants 187 Variable vz_1 = 3.956e5 //velocity [cm/s] of 1 A neutron 188 Variable g = 981.0 //gravity acceleration [cm/s^2] 189 Variable h_m = 3995 // h/m [=] A*m/s 190 191 String results 192 results ="Failure" 193 194 S1 *= 0.5*0.1 //convert to radius and [cm] 195 S2 *= 0.5*0.1 196 197 L1 *= 100.0 // [cm] 198 L1 -= apOff //correct the distance 199 200 L2 *= 100.0 201 L2 += apOff 202 del_r *= 0.1 //width of annulus, convert mm to [cm] 203 204 BS *= 0.5*0.1 //nominal BS diameter passed in, convert to radius and [cm] 205 // 21 MAR 07 SRK - use the projected BS diameter, based on a point sample aperture 206 Variable LB 207 LB = 20.1 + 1.61*BS //distance in cm from beamstop to anode plane (empirical) 208 BS = bs + bs*lb/(l2-lb) //adjusted diameter of shadow from parallax 209 210 //Start resolution calculation 211 a2 = S1*L2/L1 + S2*(L1+L2)/L1 212 lp = 1.0/( 1.0/L1 + 1.0/L2) 213 214 v_lambda = lambdaWidth^2/6.0 215 216 // if(usingLenses==1) //SRK 2007 217 if(usingLenses != 0) //SRK 2008 allows for the possibility of different numbers of lenses in header 218 v_b = 0.25*(S1*L2/L1)^2 +0.25*(2/3)*(lambdaWidth/lambda)^2*(S2*L2/lp)^2 //correction to 2nd term 219 else 220 v_b = 0.25*(S1*L2/L1)^2 +0.25*(S2*L2/lp)^2 //original form 221 endif 222 223 v_d = (DDet/2.3548)^2 + del_r^2/12.0 224 vz = vz_1 / lambda 225 yg = 0.5*g*L2*(L1+L2)/vz^2 226 v_g = 2.0*(2.0*yg^2*v_lambda) //factor of 2 correction, B. Hammouda, 2007 227 228 r0 = L2*tan(2.0*asin(lambda*inQ/(4.0*Pi) )) 229 delta = 0.5*(BS - r0)^2/v_d 230 231 if (r0 < BS) 232 inc_gamma=exp(gammln(1.5))*(1-gammp(1.5,delta)) 233 else 234 inc_gamma=exp(gammln(1.5))*(1+gammp(1.5,delta)) 235 endif 236 237 fSubS = 0.5*(1.0+erf( (r0-BS)/sqrt(2.0*v_d) ) ) 238 if (fSubS <= 0.0) 239 fSubS = 1.e-10 240 endif 241 // fr = 1.0 + sqrt(v_d)*exp(-1.0*delta) /(r0*fSubS*sqrt(2.0*Pi)) 242 // fv = inc_gamma/(fSubS*sqrt(Pi)) - r0^2*(fr-1.0)^2/v_d 243 // 244 // rmd = fr*r0 245 // v_r1 = v_b + fv*v_d +v_g 246 // 247 // rm = rmd + 0.5*v_r1/rmd 248 // v_r = v_r1 - 0.5*(v_r1/rmd)^2 249 // if (v_r < 0.0) 250 // v_r = 0.0 251 // endif 252 253 Variable kap,a_val 254 255 kap = 2*pi/lambda 256 a_val = (L1+L2)*g/2/(h_m)^2 257 258 SigmaQX = 3*(S1/L1)^2 + 3*(S2/LP)^2 + 2*(DDet/L2)^2 + 2*(r_dist/L2)^2*(lambdaWidth)^2*(cos(phi))^2 259 260 SigmaQY = 3*(S1/L1)^2 + 3*(S2/LP)^2 + 2*(DDet/L2)^2 + 2*(r_dist/L2)^2*(lambdaWidth)^2*(sin(phi))^2 + 8*(a_val/L2)^2*lambda^4*lambdaWidth^2 261 262 SigmaQX = sqrt(kap*kap/12*SigmaQX) 263 SigmaQy = sqrt(kap*kap/12*SigmaQY) 264 265 results = "success" 266 Return results 267 End 268 269 154 270 155 271 -
sans/Dev/trunk/NCNR_User_Procedures/Reduction/SANS/ProtocolAsPanel.ipf
r418 r509 2073 2073 // Print "box is now ",x1,x2,y1,y2 2074 2074 detCnt = SumCountsInBox(x1,x2,y1,y2,"RAW") 2075 if(cmpstr(tw[9],"ILL ")==0) 2076 detCnt /= 4 // for cerca detector, header is right, sum(data) is 4x too large 2077 // this is usually corrected in the Add step 2078 endif 2075 2079 // 2076 2080 kappa = detCnt/countTime/attenTrans*1.0e8/(monCnt/countTime)*(pixel/sdd)^2 -
sans/Dev/trunk/NCNR_User_Procedures/Reduction/SANS/WriteQIS.ipf
r474 r509 630 630 631 631 632 //ASCII export of data as 3-columns qx-qy-Intensity 632 // NEW additions - May 2009 633 //ASCII export of data as 7-columns qx-qy-Intensity-qz-sigmaQx-sigmaQy-fShad 633 634 //limited header information? 634 635 // … … 718 719 labelWave[13] = "" 719 720 labelWave[14] = "*** Data written from "+type+" folder and may not be a fully corrected data file ***" 720 labelWave[15] = "Data columns are Qx - Qy - I(Qx,Qy) "721 labelWave[15] = "Data columns are Qx - Qy - I(Qx,Qy) - Qz - SigmaQx - SigmaQy - fSubS(beam stop shadow)" 721 722 labelWave[16] = "" 722 723 labelWave[17] = "ASCII data created " +date()+" "+time() … … 736 737 // Endif 737 738 738 Duplicate/O data,qx_val,qy_val,z_val,qval,qz_val 739 Duplicate/O data,qx_val,qy_val,z_val,qval,qz_val,phi,r_dist 739 740 740 741 // Redimension/N=(pixelsX*pixelsY) qx_val,qy_val,z_val 741 742 // MyMat2XYZ(data,qx_val,qy_val,z_val) //x and y are [p][q] indexes, not q-vals yet 742 743 744 Variable xctr,yctr,sdd,lambda,pixSize 745 xctr = rw[16] 746 yctr = rw[17] 747 sdd = rw[18] 748 lambda = rw[26] 749 pixSize = rw[13]/10 //convert mm to cm (x and y are the same size pixels) 750 743 751 qx_val = CalcQx(p+1,q+1,rw[16],rw[17],rw[18],rw[26],rw[13]/10) //+1 converts to detector coordinate system 744 752 qy_val = CalcQy(p+1,q+1,rw[16],rw[17],rw[18],rw[26],rw[13]/10) 745 753 746 754 Redimension/N=(pixelsX*pixelsY) qx_val,qy_val,z_val 747 755 756 ///************ 757 // do everything to write out the resolution too 748 758 // un-comment these if you want to write out qz_val and qval too, then use the proper save command 749 // qval = CalcQval(p+1,q+1,rw[16],rw[17],rw[18],rw[26],rw[13]/10) 750 // qz_val = CalcQz(p+1,q+1,rw[16],rw[17],rw[18],rw[26],rw[13]/10) 751 // Redimension/N=(pixelsX*pixelsY) qz_val,qval 759 qval = CalcQval(p+1,q+1,rw[16],rw[17],rw[18],rw[26],rw[13]/10) 760 qz_val = CalcQz(p+1,q+1,rw[16],rw[17],rw[18],rw[26],rw[13]/10) 761 phi = FindPhi( pixSize*((p+1)-xctr) , pixSize*((q+1)-yctr)) //(dx,dy) 762 r_dist = sqrt( (pixSize*((p+1)-xctr))^2 + (pixSize*((q+1)-yctr))^2 ) //radial distance from ctr to pt 763 Redimension/N=(pixelsX*pixelsY) qz_val,qval,phi,r_dist 764 //everything in 1D now 765 Duplicate/O qval SigmaQX,SigmaQY,fsubS 766 767 Variable L2 = rw[18] 768 Variable BS = rw[21] 769 Variable S1 = rw[23] 770 Variable S2 = rw[24] 771 Variable L1 = rw[25] 772 Variable lambdaWidth = rw[27] 773 Variable usingLenses = rw[28] //new 2007 774 775 //Two parameters DDET and APOFF are instrument dependent. Determine 776 //these from the instrument name in the header. 777 //From conversation with JB on 01.06.99 these are the current good values 778 Variable DDet 779 NVAR apOff = root:myGlobals:apOff //in cm 780 DDet = rw[10]/10 // header value (X) is in mm, want cm here 781 782 Variable ret1,ret2,ret3,nq 783 nq = pixelsX*pixelsY 784 ii=0 785 786 do 787 get2DResolution(qval[ii],phi[ii],lambda,lambdaWidth,DDet,apOff,S1,S2,L1,L2,BS,pixSize,usingLenses,r_dist[ii],ret1,ret2,ret3) 788 SigmaQX[ii] = ret1 789 SigmaQY[ii] = ret2 790 fsubs[ii] = ret3 791 ii+=1 792 while(ii<nq) 793 794 //********************* 752 795 753 796 //not demo-compatible, but approx 8x faster!! 754 797 #if(cmpstr(stringbykey("IGORKIND",IgorInfo(0),":",";"),"pro") == 0) 755 Save/G/M="\r\n" labelWave,qx_val,qy_val,z_val as fullpath // /M=termStr specifies terminator756 // Save/G/M="\r\n" labelWave,qx_val,qy_val,qz_val,qval,z_valas fullpath // for debugging, write out everything798 // Save/G/M="\r\n" labelWave,qx_val,qy_val,z_val as fullpath // /M=termStr specifies terminator 799 Save/G/M="\r\n" labelWave,qx_val,qy_val,z_val,qz_val,SigmaQx,SigmaQy,fSubS as fullpath // for debugging, write out everything 757 800 #else 758 801 Open refNum as fullpath … … 763 806 #endif 764 807 765 Killwaves/Z spWave,labelWave,qx_val,qy_val,z_val,qval,qz_val 808 Killwaves/Z spWave,labelWave,qx_val,qy_val,z_val,qval,qz_val,sigmaQx,SigmaQy,fSubS,phi,r_dist 766 809 767 810 Print "QxQy_Export File written: ", GetFileNameFromPathNoSemi(fullPath) … … 812 855 return(0) 813 856 End 814
Note: See TracChangeset
for help on using the changeset viewer.
