Changeset 1120


Ignore:
Timestamp:
Feb 1, 2019 3:47:27 PM (4 years ago)
Author:
srkline
Message:

bug fix in the QxQy? export function

Location:
sans/Dev/trunk/NCNR_User_Procedures/Reduction/VSANS
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • sans/Dev/trunk/NCNR_User_Procedures/Reduction/VSANS/V_Instrument_Resolution.ipf

    r1119 r1120  
    709709 
    710710 
    711 //********************** 
    712 // Resolution calculation - used by the averaging routines 
    713 // to calculate the resolution function at each q-value 
    714 // - the return value is not used 
    715 // 
    716 // equivalent to John's routine on the VAX Q_SIGMA_AVE.FOR 
    717 // Incorporates eqn. 3-15 from J. Appl. Cryst. (1995) v. 28 p105-114 
    718 // 
    719 // - 21 MAR 07 uses projected BS diameter on the detector 
    720 // - APR 07 still need to add resolution with lenses. currently there is no flag in the  
    721 //          raw data header to indicate the presence of lenses. 
    722 // 
    723 // - Aug 07 - added input to switch calculation based on lenses (==1 if in) 
    724 // 
    725 // - SANS -- called by CircSectAvg.ipf and RectAnnulAvg.ipf 
    726 // 
    727 // - VSANS -- called in VC_fDoBinning_QxQy2D(folderStr, binningType) 
    728 // 
    729 // DDet is the detector pixel resolution 
    730 // apOff is the offset between the sample aperture and the sample position 
    731 // 
    732 // 
    733 // INPUT: 
    734 // inQ = q-value [1/A] 
    735 // lambda = wavelength [A] 
    736 // lambdaWidth = [dimensionless] 
    737 // DDet = detector pixel resolution [cm]        **assumes square pixel 
    738 // apOff = sample aperture to sample distance [cm] 
    739 // S1 = source aperture diameter [mm] 
    740 // S2 = sample aperture diameter [mm] 
    741 // L1 = source to sample distance [m]  
    742 // L2 = sample to detector distance [m] 
    743 // BS = beam stop diameter [mm] 
    744 // del_r = step size [mm] = binWidth*(mm/pixel)  
    745 // usingLenses = flag for lenses = 0 if no lenses, non-zero if lenses are in-beam 
    746 // 
    747 // OUPUT: 
    748 // SigmaQ 
    749 // QBar 
    750 // fSubS 
    751 // 
    752 // 
    753 Function V_getResolution_old(inQ,lambda,lambdaWidth,DDet,apOff,S1,S2,L1,L2,BS,del_r,usingLenses,SigmaQ,QBar,fSubS) 
    754         Variable inQ, lambda, lambdaWidth, DDet, apOff, S1, S2, L1, L2, BS, del_r,usingLenses 
    755         Variable &SigmaQ, &QBar, &fSubS         //these are the output quantities at the input Q value 
    756          
    757         //lots of calculation variables 
    758         Variable a2, q_small, lp, v_lambda, v_b, v_d, vz, yg, v_g 
    759         Variable r0, delta, inc_gamma, fr, fv, rmd, v_r1, rm, v_r 
    760  
    761         //Constants 
    762         Variable vz_1 = 3.956e5         //velocity [cm/s] of 1 A neutron 
    763         Variable g = 981.0                              //gravity acceleration [cm/s^2] 
    764  
    765  
    766         S1 *= 0.5*0.1                   //convert to radius and [cm] 
    767         S2 *= 0.5*0.1 
    768  
    769         L1 *= 100.0                     // [cm] 
    770         L1 -= apOff                             //correct the distance 
    771  
    772         L2 *= 100.0 
    773         L2 += apOff 
    774         del_r *= 0.1                            //width of annulus, convert mm to [cm] 
    775          
    776         BS *= 0.5*0.1                   //nominal BS diameter passed in, convert to radius and [cm] 
    777         // 21 MAR 07 SRK - use the projected BS diameter, based on a point sample aperture 
    778         Variable LB 
    779         LB = 20.1 + 1.61*BS                     //distance in cm from beamstop to anode plane (empirical) 
    780         BS = bs + bs*lb/(l2-lb)         //adjusted diameter of shadow from parallax 
    781          
    782         //Start resolution calculation 
    783         a2 = S1*L2/L1 + S2*(L1+L2)/L1 
    784         q_small = 2.0*Pi*(BS-a2)*(1.0-lambdaWidth)/(lambda*L2) 
    785         lp = 1.0/( 1.0/L1 + 1.0/L2) 
    786  
    787         v_lambda = lambdaWidth^2/6.0 
    788          
    789 //      if(usingLenses==1)                      //SRK 2007 
    790         if(usingLenses != 0)                    //SRK 2008 allows for the possibility of different numbers of lenses in header 
    791                 v_b = 0.25*(S1*L2/L1)^2 +0.25*(2/3)*(lambdaWidth/lambda)^2*(S2*L2/lp)^2         //correction to 2nd term 
    792         else 
    793                 v_b = 0.25*(S1*L2/L1)^2 +0.25*(S2*L2/lp)^2              //original form 
    794         endif 
    795          
    796         v_d = (DDet/2.3548)^2 + del_r^2/12.0                    //the 2.3548 is a conversion from FWHM->Gauss, see http://mathworld.wolfram.com/GaussianFunction.html 
    797         vz = vz_1 / lambda 
    798         yg = 0.5*g*L2*(L1+L2)/vz^2 
    799         v_g = 2.0*(2.0*yg^2*v_lambda)                                   //factor of 2 correction, B. Hammouda, 2007 
    800  
    801         r0 = L2*tan(2.0*asin(lambda*inQ/(4.0*Pi) )) 
    802         delta = 0.5*(BS - r0)^2/v_d 
    803  
    804         if (r0 < BS)  
    805                 inc_gamma=exp(gammln(1.5))*(1-gammp(1.5,delta)) 
    806         else 
    807                 inc_gamma=exp(gammln(1.5))*(1+gammp(1.5,delta)) 
    808         endif 
    809  
    810         fSubS = 0.5*(1.0+erf( (r0-BS)/sqrt(2.0*v_d) ) ) 
    811         if (fSubS <= 0.0)  
    812                 fSubS = 1.e-10 
    813         endif 
    814         fr = 1.0 + sqrt(v_d)*exp(-1.0*delta) /(r0*fSubS*sqrt(2.0*Pi)) 
    815         fv = inc_gamma/(fSubS*sqrt(Pi)) - r0^2*(fr-1.0)^2/v_d 
    816  
    817         rmd = fr*r0 
    818         v_r1 = v_b + fv*v_d +v_g 
    819  
    820         rm = rmd + 0.5*v_r1/rmd 
    821         v_r = v_r1 - 0.5*(v_r1/rmd)^2 
    822         if (v_r < 0.0)  
    823                 v_r = 0.0 
    824         endif 
    825         QBar = (4.0*Pi/lambda)*sin(0.5*atan(rm/L2)) 
    826         SigmaQ = QBar*sqrt(v_r/rmd^2 +v_lambda) 
    827  
    828  
    829 // more readable method for calculating the variance in Q 
    830 // EXCEPT - this is calculated for Qo, NOT qBar 
    831 // (otherwise, they are nearly equivalent, except for close to the beam stop) 
    832 //      Variable kap,a_val,a_val_l2,m_h 
    833 //      g = 981.0                               //gravity acceleration [cm/s^2] 
    834 //      m_h     = 252.8                 // m/h [=] s/cm^2 
    835 //       
    836 //      kap = 2*pi/lambda 
    837 //      a_val = L2*(L1+L2)*g/2*(m_h)^2 
    838 //      a_val_L2 = a_val/L2*1e-16               //convert 1/cm^2 to 1/A^2 
    839 // 
    840 //      sigmaQ = 0 
    841 //      sigmaQ = 3*(S1/L1)^2 
    842 //       
    843 //      if(usingLenses != 0) 
    844 //              sigmaQ += 2*(S2/lp)^2*(lambdaWidth)^2   //2nd term w/ lenses 
    845 //      else     
    846 //              sigmaQ += 2*(S2/lp)^2                                           //2nd term w/ no lenses 
    847 //      endif 
    848 //       
    849 //      sigmaQ += (del_r/L2)^2 
    850 //      sigmaQ += 2*(r0/L2)^2*(lambdaWidth)^2 
    851 //      sigmaQ += 4*(a_val_l2)^2*lambda^4*(lambdaWidth)^2 
    852 //       
    853 //      sigmaQ *= kap^2/12 
    854 //      sigmaQ = sqrt(sigmaQ) 
    855  
    856  
    857         Return (0) 
    858 End 
    859  
    860  
     711 
     712 
     713 
  • sans/Dev/trunk/NCNR_User_Procedures/Reduction/VSANS/V_Write_VSANS_QIS.ipf

    r1119 r1120  
    518518        //      qval = CalcQval(p+1,q+1,rw[16],rw[17],rw[18],rw[26],rw[13]/10) 
    519519                Duplicate/O qTot,phi,r_dist 
    520                 Variable pixSizeX,pixSizeY 
     520                Variable pixSizeX,pixSizeY,xctr,yctr 
    521521                pixSizeX = V_getDet_x_pixel_size(type,detStr) 
    522522                pixSizeY = V_getDet_y_pixel_size(type,detStr) 
    523523 
    524                 phi = V_FindPhi( pixSizeX*((p+1)-bCentX) , pixSizeY*((q+1)-bCentY)+(2)*yg_d)            //(dx,dy+yg_d) 
    525                 r_dist = sqrt(  (pixSizeX*((p+1)-bCentX))^2 +  (pixSizeY*((q+1)-bCentY)+(2)*yg_d)^2 )           //radial distance from ctr to pt 
     524                xctr = V_getDet_beam_center_x_pix(type,detStr) 
     525                yctr = V_getDet_beam_center_y_pix(type,detStr) 
     526                phi = V_FindPhi( pixSizeX*((p+1)-xctr) , pixSizeY*((q+1)-yctr)+(2)*yg_d)                //(dx,dy+yg_d) 
     527                r_dist = sqrt(  (pixSizeX*((p+1)-xctr))^2 +  (pixSizeY*((q+1)-yctr)+(2)*yg_d)^2 )               //radial distance from ctr to pt 
    526528         
    527529                //make everything in 1D now 
     
    591593                KillWaves/Z qx_val_s,qy_val_s,z_val_s,qz_val_s,SigmaQx_s,SigmaQy_s,fSubS_s,sw,sw_s 
    592594                 
    593                 Killwaves/Z qx_val,qy_val,z_val,qval,qz_val,sigmaQx,SigmaQy,fSubS 
     595                Killwaves/Z qval,sigmaQx,SigmaQy,fSubS,phi,r_dist 
    594596                 
    595597                Print "QxQy_Export File written: ", V_GetFileNameFromPathNoSemi(detSavePath) 
Note: See TracChangeset for help on using the changeset viewer.