Changeset 467 for sans/Dev

Timestamp:

Jan 7, 2009 1:17:09 PM (14 years ago)

Author:

srkline

Message:

Final fixes, w/Charles Dewhurst, to finally calculate Qx and Qy correctly as vector components of the vector Q, as it is projected onto the XY planar detector. An additional function has been added to calculate Qz as well, but is not written out anywhere. This calculation only has an effect at the highest q-values, for QxQy? export. At a Q ~ 0.45, the correct Qx and Qy are about 2% different than the old calculation. |Q| has always been correct.

Location:

sans/Dev/trunk

Files:

• NCNR_Help_Files/USANSReductionHelp_v2.21.ihf (modified) (previous)
• NCNR_User_Procedures/Reduction/SANS/RawWindowHook.ipf (modified) (6 diffs)
• NCNR_User_Procedures/Reduction/SANS/WriteQIS.ipf (modified) (3 diffs)

sans/Dev/trunk/NCNR_User_Procedures/Reduction/SANS/RawWindowHook.ipf

@@ -271 +271 @@
 //returned magnitude of Q is in 1/Angstroms
 //
-// repaired the dumb error of incorrect qx and qy calculation 3 dec 08 SRK (Lionel...)
-//
 Function CalcQval(xaxval,yaxval,xctr,yctr,sdd,lam,pixSize)
         Variable xaxval,yaxval,xctr,yctr,sdd,lam,pixSize
         
-        Variable dx,dy,qval,theta,dist
+        Variable dx,dy,qval,two_theta,dist
         
         Variable pixSizeX=pixSize
@@ -286 +284 @@
         dist = sqrt(dx^2 + dy^2)
         
-        theta = atan(dist/sdd)
-
-        qval = 4*Pi/lam*sin(theta/2)
+        two_theta = atan(dist/sdd)
+
+        qval = 4*Pi/lam*sin(two_theta/2)
         
         return qval
@@ -301 +299 @@
 //wavelength is in Angstroms
 //
-// repaired the dumb error of incorrect qx and qy calculation 3 dec 08 SRK (Lionel...)
+// repaired incorrect qx and qy calculation 3 dec 08 SRK (Lionel and C. Dewhurst)
+// now properly accounts for qz
 //
 Function CalcQX(xaxval,yaxval,xctr,yctr,sdd,lam,pixSize)
         Variable xaxval,yaxval,xctr,yctr,sdd,lam,pixSize
 
-        Variable qx,qval,phi,dx,dy
+        Variable qx,qval,phi,dx,dy,dist,two_theta
         
         qval = CalcQval(xaxval,yaxval,xctr,yctr,sdd,lam,pixSize)
@@ -315 +314 @@
         phi = FindPhi(dx,dy)
         
-        qx = qval*cos(phi)
+        //get scattering angle to project onto flat detector => Qr = qval*cos(theta)
+        dist = sqrt(dx^2 + dy^2)
+        two_theta = atan(dist/sdd)
+
+        qx = qval*cos(two_theta/2)*cos(phi)
         
         return qx
@@ -327 +330 @@
 //wavelength is in Angstroms
 //
-// repaired the dumb error of incorrect qx and qy calculation 3 dec 08 SRK (Lionel...)
+// repaired incorrect qx and qy calculation 3 dec 08 SRK (Lionel and C. Dewhurst)
+// now properly accounts for qz
 //
 Function CalcQY(xaxval,yaxval,xctr,yctr,sdd,lam,pixSize)
         Variable xaxval,yaxval,xctr,yctr,sdd,lam,pixSize
         
-        Variable dy,qval,dx,phi,qy
+        Variable dy,qval,dx,phi,qy,dist,two_theta
         
         qval = CalcQval(xaxval,yaxval,xctr,yctr,sdd,lam,pixSize)
@@ -341 +345 @@
         phi = FindPhi(dx,dy)
         
-        qy = qval*sin(phi)
+        //get scattering angle to project onto flat detector => Qr = qval*cos(theta)
+        dist = sqrt(dx^2 + dy^2)
+        two_theta = atan(dist/sdd)
+        
+        qy = qval*cos(two_theta/2)*sin(phi)
         
         return qy
 End
 
+//calculates just the z-component of the q-vector, not measured on the detector
+//input/output is the same as CalcQval()
+//ALL inputs are in detector coordinates
+//NOTE: detector locations passed in are pixel = 0.5cm real space on the Ordela detector
+//sdd is in meters
+//wavelength is in Angstroms
+//
+// not actually used, but here for completeness if anyone asks
+//
+Function CalcQZ(xaxval,yaxval,xctr,yctr,sdd,lam,pixSize)
+        Variable xaxval,yaxval,xctr,yctr,sdd,lam,pixSize
+        
+        Variable dy,qval,dx,phi,qz,dist,two_theta
+        
+        qval = CalcQval(xaxval,yaxval,xctr,yctr,sdd,lam,pixSize)
+        
+        sdd *=100               //convert to cm
+        
+        //get scattering angle to project onto flat detector => Qr = qval*cos(theta)
+        dx = (xaxval - xctr)*pixSize            //delta x in cm
+        dy = (yaxval - yctr)*pixSize            //delta y in cm
+        dist = sqrt(dx^2 + dy^2)
+        two_theta = atan(dist/sdd)
+        
+        qz = qval*sin(two_theta/2)
+        
+        return qz
+End
 
 //phi is defined from +x axis, proceeding CCW around [0,2Pi]

sans/Dev/trunk/NCNR_User_Procedures/Reduction/SANS/WriteQIS.ipf

@@ -736 +736 @@
 //      Endif
         
-        Duplicate/O data,qx_val,qy_val,z_val
+        Duplicate/O data,qx_val,qy_val,z_val,qval,qz_val
         
 //      Redimension/N=(pixelsX*pixelsY) qx_val,qy_val,z_val
@@ -743 +743 @@
         qx_val = CalcQx(p+1,q+1,rw[16],rw[17],rw[18],rw[26],rw[13]/10)          //+1 converts to detector coordinate system
         qy_val = CalcQy(p+1,q+1,rw[16],rw[17],rw[18],rw[26],rw[13]/10)
-
+        
         Redimension/N=(pixelsX*pixelsY) qx_val,qy_val,z_val
+        
+        // un-comment these if you want to write out qz_val and qval too, then use the proper save command
+//      qval = CalcQval(p+1,q+1,rw[16],rw[17],rw[18],rw[26],rw[13]/10)
+//      qz_val = CalcQz(p+1,q+1,rw[16],rw[17],rw[18],rw[26],rw[13]/10)
+//      Redimension/N=(pixelsX*pixelsY) qz_val,qval
 
         //not demo-compatible, but approx 8x faster!!   
 #if(cmpstr(stringbykey("IGORKIND",IgorInfo(0),":",";"),"pro") == 0)     
-        Save/G/M="\r\n" labelWave,qx_val,qy_val,z_val as fullpath       // /M=termStr specifies terminator
+        Save/G/M="\r\n" labelWave,qx_val,qy_val,z_val as fullpath       // /M=termStr specifies terminator      
+//      Save/G/M="\r\n" labelWave,qx_val,qy_val,qz_val,qval,z_val as fullpath   // for debugging, write out everything
 #else
         Open refNum as fullpath
@@ -757 +763 @@
 #endif
         
-        Killwaves/Z spWave,labelWave,qx_val,qy_val,z_val
+        Killwaves/Z spWave,labelWave,qx_val,qy_val,z_val,qval,qz_val
         
         Print "QxQy_Export File written: ", GetFileNameFromPathNoSemi(fullPath)