Ignore:
Timestamp:
Jun 21, 2011 1:07:46 PM (12 years ago)
Author:
srkline
Message:

Fixed comparison bug in NSORT / MRED / CatTable?

Fixed bug in writing of RAW VAX files (only used in simulation) that incorrectly converted 32 bit integer to 16 bit integer before compression, leading to negative counts for n ~> 37000.

Added source aperture and gravity to the MC simulation for a realistic beam profile.

Location:
sans/Dev/trunk/NCNR_User_Procedures/Reduction/SANS
Files:
7 edited

Legend:

Unmodified
Added
Removed
  • sans/Dev/trunk/NCNR_User_Procedures/Reduction/SANS/CatVSTable.ipf

    r707 r811  
    151151                        //write to notebook that file was not found 
    152152                        //if string is not a number, report the error 
    153                         if(str2num(partialName) == NaN) 
     153                        if(numtype(str2num(partialName)) == 2) 
    154154                                str = "this file was not found: "+partialName+"\r\r" 
    155155                                //Notebook CatWin,font="Times",fsize=12,text=str 
     
    508508                        //write to notebook that file was not found 
    509509                        //if string is not a number, report the error 
    510                         if(str2num(partialName) == NaN) 
     510                        if(numtype(str2num(partialName)) == 2) 
    511511                                str = "this file was not found: "+partialName+"\r\r" 
    512512                                Notebook CatWin,font="Times",fsize=12,text=str 
     
    677677                        //write to notebook that file was not found 
    678678                        //if string is not a number, report the error 
    679                         if(str2num(partialName) == NaN) 
     679                        if(numtype(str2num(partialName)) == 2) 
    680680                                str = "this file was not found: "+partialName+"\r\r" 
    681681                                Notebook CatWin,font="Times",fsize=12,text=str 
  • sans/Dev/trunk/NCNR_User_Procedures/Reduction/SANS/MultScatter_MonteCarlo_2D.ipf

    r800 r811  
    2020// whatever XOP crashing was happining during threading is really unlikely to be from the RNG 
    2121// 
    22 // -- June 2010 - calls from different threads to the same RNG really seems to cause a crash. Probably as soon 
     22// -- June 2010 - calls from different threads to the same RNG will absolutely cause a crash. Probably as soon 
    2323//                              as the different threads try to call at the same time. Found this out by accident doing the  
    2424//                              wavelength spread. Each thread called ran3 at that point, and the crash came quickly. Went 
     
    8585// --- TO ADD --- 
    8686// X- wavelength distribution = another RNG to select the wavelength 
    87 // ---- done Jun 2010, approximating the wavelength distribution as a Gaussian, based on the triangular 
     87// DONE Jun 2010, approximating the wavelength distribution as a Gaussian, based on the triangular 
    8888//                      FWHM. Wavelength distribution added to XOP too, and now very accurately matches the shape of the 1D 
    8989//                      simulation. 
     
    9696//               data can be plotted, but only by hand right now. EC and blocked beam are combined. 
    9797// 
    98 // -- divergence / size of the incoming beam. Currently everything is parallel, and anything that is transmitted 
     98 
     99// X- divergence / size of the incoming beam. Currently everything is parallel, and anything that is transmitted 
    99100//              simply ends up in (xCtr,yCtr), and the "real" profile of the beam is not captured. 
    100  
     101// DONE, 21 JUN 11 SRK 
     102// point is picked in source aperture, then sample aperture. This sets the initial direction of the neutron. 
     103// Gravity is included, lowering every neutron by yg_d. This affects qy only.Simulated beam center makes sense 
     104// now both in size and xy positon. Gravity effects can be seen in the beam spot itself (fall + spread) and in the scattering 
     105// pattern at extreme cases (sharp peaks, wide spread, long SDD, long wavelength, small source aperture) 
    101106 
    102107 
     
    138143        nthreads = 1 
    139144#endif 
    140          
    141 //      nthreads = 1 
     145 
     146         
     147//      nthreads = 2 
     148 
     149 
    142150        NVAR mt=root:myGlobals:gThreadGroupID 
    143151        mt = ThreadGroupCreate(nthreads) 
     
    382390        Variable NMultipleCoherent,NCoherentEvents 
    383391        Variable deltaLam,v1,v2,currWavelength,rsq,fac          //for simulating wavelength distribution 
     392        Variable SSD, sourAp, souXX, souYY, magn                //source-to-sample, and source Ap radius for initlal trajectory 
     393 
     394        Variable vz_1 = 3.956e5         //velocity [cm/s] of 1 A neutron 
     395        Variable g = 981.0                              //gravity acceleration [cm/s^2] 
     396        Variable yg_d           //fall due to gravity 
     397         
    384398         
    385399        // don't set to other than one here. Detector efficiency is handled outside, only passing the number of  
     
    399413        sig_sas = inputWave[10] 
    400414        deltaLam = inputWave[11] 
     415        SSD = inputWave[12]                     // in cm, like SDD 
     416        sourAp = inputWave[13]          // radius, in cm, like r1 and r2 
     417         
     418//      print SSD, sourAp 
    401419         
    402420//      SetRandomSeed 0.1               //to get a reproduceable sequence 
     
    462480// NOW, start the loop, throwing neutrons at the sample. 
    463481        do 
    464                 Vx = 0.0                        // Initialize direction vector. 
    465                 Vy = 0.0 
    466                 Vz = 1.0 
     482//              Vx = 0.0                        // Initialize direction vector. 
     483//              Vy = 0.0 
     484//              Vz = 1.0 
    467485                 
    468486                Theta = 0.0             //      Initialize scattering angle. 
     
    474492                incoherentEvent = 0 
    475493                coherentEvent = 0 
    476                  
     494         
     495                // pick point in source aperture area 
     496                do                                      //      Makes sure position is within circle. 
     497                        ran = abs(enoise(1))            //[0,1] 
     498                        souxx = 2.0*sourAp*(Ran-0.5)            //X beam position of neutron entering sample. 
     499                        ran = abs(enoise(1))            //[0,1] 
     500                        souyy = 2.0*sourAp*(Ran-0.5)            //Y beam position ... 
     501                        RR = SQRT(souxx*souxx+souyy*souyy)              //Radial position of neutron in incident beam. 
     502                while(rr>sourAp) 
     503                                                 
     504                // pick point in sample aperture 
    477505                do                                      //      Makes sure position is within circle. 
    478506                        ran = abs(enoise(1))            //[0,1] 
     
    497525                currWavelength = (v2*fac)*deltaLam*wavelength/2.35 + wavelength 
    498526                 
     527//              if(n1 == 1) 
     528//                      Print "lambda ", currWavelength 
     529//              endif    
     530//               
     531                magn = sqrt((souXX - xx)^2 + (souYY - yy)^2 + ssd^2) 
     532                Vx = (souXX - xx)/magn          // Initialize direction vector. 
     533                Vy = (souYY - yy)/magn 
     534                Vz = (ssd - 0)/magn 
     535                 
     536//              Vx = 0.0                        // Initialize direction vector. 
     537//              Vy = 0.0 
     538//              Vz = 1.0 
     539                 
     540                 
     541// 
     542//              if(n1 == 1) 
     543//                      Print "vx, vy, vz, mag",vx,vy,vz,sqrt(vx^2+vy^2+vz^2) 
     544//              endif            
     545                 
    499546                 
    500547                do    //Scattering Loop, will exit when "done" == 1 
     
    505552                        err = NewDirection(vx,vy,vz,Theta,Phi)          //vx,vy,vz is updated, theta, phi unchanged by function 
    506553 
     554//                      if(vx != 0) 
     555//                              Print "vx, vy, vz, mag" = vx,vy,vz,sqrt(vx^2+vy^2+vz^2) 
     556//                      endif    
     557                         
     558                         
    507559                        //X,Y,Z-POSITION OF SCATTERING EVENT. 
    508560                        xx += ll*vx 
     
    578630                                Endif 
    579631                                 
     632                                // calculate fall due to gravity (in cm) (note that it is negative) 
     633                                YG_d = -0.5*g*SDD*(SSD+SDD)*(currWavelength/vz_1)^2 
     634                                 
     635//                              yg_d=0 
     636                                 
     637                                if(n1 == 1) 
     638//                                      Print "gravity fall (cm) = ",Yg_d 
     639                                        Print "gravity fall at mean lam (cm) = ",-0.5*g*SDD*(SSD+SDD)*(wavelength/vz_1)^2 
     640                                endif    
     641                 
    580642                                if(index != 0)          //the neutron interacted at least once, figure out where it ends up 
    581643 
     
    590652                                         
    591653                                        // is it on the detector?        
    592                                         FindPixel(testQ,testPhi,currWavelength,sdd,pixSize,xCtr,yCtr,xPixel,yPixel) 
     654                                        FindPixel(testQ,testPhi,currWavelength,yg_d,sdd,pixSize,xCtr,yCtr,xPixel,yPixel) 
    593655                                         
    594656                                        if(xPixel != -1 && yPixel != -1) 
    595657                                                //if(index==1)  // only the single scattering events 
     658                                                        if( xPixel > 127 || yPixel > 127) 
     659                                                                print "error XY=",xPixel,yPixel 
     660                                                        endif 
    596661                                                        MC_linear_data[xPixel][yPixel] += 1             //this is the total scattering, including multiple scattering 
    597662                                                //endif 
     
    641706                                 
    642707                                        // then it must be a transmitted neutron 
    643                                         // don't need to calculate, just increment the proper counters 
     708                                        // Now, calculate where it lands 
    644709                                         
    645                                         MC_linear_data[xCtr+xx/pixsize][yCtr+yy/pixsize] += 1 
    646                                         isOn += 1 
    647                                         nt[0] += 1 
     710                                        Theta_z = acos(Vz)              // Angle WITH respect to z axis. 
     711                                        testQ = 2*pi*sin(theta_z)/currWavelength 
     712                                         
     713                                        testPhi = MC_FindPhi(Vx,Vy)             //use the exiting phi value as defined by Vx and Vy 
     714                                         
     715                                        // is it on the detector?        
     716                                        FindPixel(testQ,testPhi,currWavelength,yg_d,sdd,pixSize,xCtr,yCtr,xPixel,yPixel) 
     717                                         
     718                                        if(xPixel != -1 && yPixel != -1) 
     719                                                //if(index==1)  // only the single scattering events 
     720                                                        if( xPixel > 127 || yPixel > 127) 
     721                                                                print "error XY=",xPixel,yPixel 
     722                                                        endif 
     723                                                        MC_linear_data[xPixel][yPixel] += 1             //this is the total scattering, including multiple scattering 
     724                                                //endif 
     725                                                        isOn += 1               // neutron that lands on detector 
     726                                        endif 
     727                                         
     728//                                      MC_linear_data[xCtr+xx/pixsize][yCtr+yy/pixsize] += 1 
     729//                                      isOn += 1 
     730                                        If(theta_z < theta_max) 
     731                                                //Choose index for scattering angle array. 
     732                                                //IND = NINT(THETA_z/DTH + 0.4999999) 
     733                                                ind = round(THETA_z/DTH + 0.4999999)            //round is eqivalent to nint() 
     734                                                NT[ind] += 1                    //Increment bin for angle. 
     735                                        endif 
     736                                        //nt[0] += 1            // may not be at zero angle anmore 
    648737                                         
    649738                                endif           //if interacted 
     
    652741        while(n1 < imon) 
    653742 
    654 //      Print "Monte Carlo Done" 
    655         results[0] = n1 
    656         results[1] = n2 
    657         results[2] = isOn 
    658         results[3] = NScatterEvents             //sum of # of times that neutrons scattered (coh+incoh) 
    659         results[4] = NSingleCoherent            //# of events that are single, coherent 
    660         results[5] = NMultipleCoherent  //# of scattered neutrons that are coherently scattered more than once 
    661         results[6] = NMultipleScatter           //# of scattered neutrons that are scattered more than once (coh and/or incoh) 
    662         results[7] = NCoherentEvents            //# of scattered neutrons that are scattered coherently one or more times 
    663          
     743        Print "Non-XOP Monte Carlo Done" 
     744//      results[0] = n1 
     745//      results[1] = n2 
     746//      results[2] = isOn 
     747//      results[3] = NScatterEvents             //sum of # of times that neutrons scattered (coh+incoh) 
     748//      results[4] = NSingleCoherent            //# of events that are single, coherent 
     749//      results[5] = NMultipleCoherent  //# of scattered neutrons that are coherently scattered more than once 
     750//      results[6] = NMultipleScatter           //# of scattered neutrons that are scattered more than once (coh and/or incoh) 
     751//      results[7] = NCoherentEvents            //# of scattered neutrons that are scattered coherently one or more times 
     752 
     753 
     754        Variable xc,yc 
     755        xc=inputWave[3] 
     756        yc=inputWave[4] 
     757        results[0] = inputWave[9]+inputWave[10]         //total XS 
     758        results[1] = inputWave[10]                                              //SAS XS 
     759        results[2] = n2                                         //number that interact n2 
     760        results[3] = isOn       - MC_linear_data[xc][yc]                                //# reaching detector minus Q(0) 
     761        results[4] = NScatterEvents/n2                          //avg# times scattered 
     762        results[5] = NSingleCoherent/NCoherentEvents                                            //single coherent fraction 
     763        results[6] = NMultipleCoherent/NCoherentEvents                          //multiple coherent fraction 
     764        results[7] = NMultipleScatter/n2                                //multiple scatter fraction 
     765        results[8] = (n1-n2)/n1                 //transmitted fraction 
     766         
     767                 
    664768//      Print "# absorbed = ",n3 
    665769 
     
    790894End 
    791895 
    792 ThreadSafe Function FindPixel(testQ,testPhi,lam,sdd,pixSize,xCtr,yCtr,xPixel,yPixel) 
    793         Variable testQ,testPhi,lam,sdd,pixSize,xCtr,yCtr,&xPixel,&yPixel 
     896 
     897 
     898// xCtr and yCtr here are the "optical" center of the detector ~(65,65) and the full fall due to  
     899// gravity is calculated from this horizontal axis 
     900// 
     901ThreadSafe Function FindPixel(testQ,testPhi,lam,yg_d,sdd,pixSize,xCtr,yCtr,xPixel,yPixel) 
     902        Variable testQ,testPhi,lam,yg_d,sdd,pixSize,xCtr,yCtr,&xPixel,&yPixel 
    794903 
    795904        Variable theta,dy,dx,qx,qy 
     
    798907        qy = testQ*sin(testPhi) 
    799908 
     909// correct qy for gravity 
     910// qy = 4*pi/lam * (theta/2) 
     911        qy += 4*pi/lam*(yg_d/sdd/2) 
     912         
     913         
    800914        //convert qx,qy to pixel locations relative to # of pixels x, y from center 
    801915        theta = 2*asin(qy*lam/4/pi) 
     
    811925         
    812926        //if on detector, return xPix and yPix values, otherwise -1 
    813         if(yPixel > pixelsY || yPixel < 0) 
     927        if(yPixel >= pixelsY || yPixel < 0) 
    814928                yPixel = -1 
    815929        endif 
    816         if(xPixel > pixelsX || xPixel < 0) 
     930        if(xPixel >= pixelsX || xPixel < 0) 
    817931                xPixel = -1 
    818932        endif 
     
    13321446        inputWave[10] = sig_sas 
    13331447        inputWave[11] = deltaLam 
    1334 //      inputWave[] 12-14 are currently unused 
     1448         
     1449        inputWave[12] = sourceToSampleDist()            // function returns dist in cm 
     1450        inputWave[13] = sourceApertureDiam()/2          //      function returns diam in cm, this is radius 
     1451//      inputWave[] 14 are currently unused 
    13351452 
    13361453        linear_data = 0         //initialize 
     
    14021519         
    14031520        Variable rad=beamstopDiam()/2           //beamstop radius in cm 
     1521        // the function detectorOffset() does not take gravity into account, and resets the beam center to (64.5,64.5) 
     1522        // after the MC simulation is done (displayConfigurationText is always called) 
     1523// 
     1524// so this (locally) changes the beam center, just for the correct placement of the beam stop 
     1525// the data files will have the wrong 
     1526        Variable vz_1 = 3.956e5         //velocity [cm/s] of 1 A neutron 
     1527        Variable g = 981.0                              //gravity acceleration [cm/s^2] 
     1528        Variable yg_d,ssd               //fall due to gravity 
     1529 
     1530        ssd = sourceToSampleDist() 
     1531        YG_d = -0.5*g*SDD*(SSD+SDD)*(wavelength/vz_1)^2         // fall in cm (negative value) 
     1532         
     1533        xCtr = 64.5      + round(2*rw[19])              // I'm always off by one for some reason, so start at 65.5? 
     1534        yCtr = 65 + yg_d/0.5                                    // this will lower the beam center 
     1535//      rw[16] = xCtr 
     1536//      rw[17] = yCtr 
     1537         
     1538        Print "Gravity q* = ",-2*pi/wavelength*2*yg_d/sdd 
     1539 
     1540                 
    14041541        if(MC_BS_in) 
    14051542                rad /= 0.5                              //convert cm to pixels 
     
    14071544                Duplicate/O linear_data,root:Packages:NIST:SAS:tmp_mask//,root:Packages:NIST:SAS:MC_linear_data 
    14081545                WAVE tmp_mask = root:Packages:NIST:SAS:tmp_mask 
     1546                 
    14091547                tmp_mask = (sqrt((p-xCtr)^2+(q-yCtr)^2) < rad) ? 0 : 1          //behind beamstop = 0, away = 1 
    14101548                 
     
    14581596        // re-average the 2D data 
    14591597        S_CircularAverageTo1D("SAS") 
    1460          
     1598                 
    14611599        // put the new result into the simulation folder 
    14621600        Fake1DDataFolder(qval,aveint,sigave,sigmaQ,qbar,fSubs,"Simulation")      
     
    14741612                endif 
    14751613        endif 
    1476  
    14771614 
    14781615        return(0) 
  • sans/Dev/trunk/NCNR_User_Procedures/Reduction/SANS/MultipleReduce.ipf

    r764 r811  
    607607                        //write to notebook that file was not found 
    608608                        //if string is not a number, report the error 
    609                         if(str2num(partialName) == NaN) 
     609                        if(numtype(str2num(partialName)) == 2) 
    610610                                str = "this file was not found: "+partialName+"\r\r" 
    611611                                //Notebook CatWin,font="Times",fsize=12,text=str 
     
    664664                loc = strsearch(labels[ii], findThisStr, 0 ,2)          //2==case insensitive, but Igor 5 specific 
    665665                if(loc != -1) 
    666                         Print "Remove w[ii] = ",num,"  ",labels[ii] 
     666                        Print "Remove w[ii] = ",runnum[ii],"  ",labels[ii] 
    667667                        DeletePoints ii, 1, filenames,suffix,labels,sdd,runnum,isTrans 
    668668                endif 
  • sans/Dev/trunk/NCNR_User_Procedures/Reduction/SANS/NCNR_DataReadWrite.ipf

    r806 r811  
    20042004        WAVE w=tmp_data 
    20052005 
     2006        // check for data values that are too large. the maximum VAX compressed data value is 2767000 
     2007        // 
     2008        WaveStats/Q w 
     2009        if(V_max > 2767000) 
     2010                Abort "Some individual pixel values are > 2767000 and the data can't be saved in VAX format" 
     2011        Endif 
    20062012         
    20072013        //check each wave 
     
    20442050        tmpFile=0 
    20452051         
    2046         Make/O/W/N=16401 dataWRecMarkers 
     2052//      Make/O/W/N=16401 dataWRecMarkers                        // don't convert to 16 bit here, rather write to file as 16 bit later 
     2053        Make/O/I/N=16401 dataWRecMarkers 
    20472054        AddRecordMarkers(w,dataWRecMarkers) 
    2048          
     2055                 
    20492056        // need to re-compress?? maybe never a problem, but should be done for the odd case 
    20502057        dataWRecMarkers = CompressI4toI2(dataWRecMarkers)               //unless a pixel value is > 32767, the same values are returned 
  • sans/Dev/trunk/NCNR_User_Procedures/Reduction/SANS/NCNR_Utils.ipf

    r801 r811  
    320320 
    321321///////////////////////////////////////////////// 
    322 ///      
    323         ////// this is all experimental, inclusion of gravity effect into the parallel component 
    324         // 
    325         Variable yg_d,acc,sdd,ssd,lambda0,DL_L,sig_l 
    326         Variable var_qlx,var_qly,var_ql,qx,qy,sig_perp,sig_para 
     322/////    
     323//      ////// this is all experimental, inclusion of gravity effect into the parallel component 
     324////    // 
     325//      Variable yg_d,acc,sdd,ssd,lambda0,DL_L,sig_l 
     326//      Variable var_qlx,var_qly,var_ql,qx,qy,sig_perp,sig_para 
    327327//      G = 981.  //!   ACCELERATION OF GRAVITY, CM/SEC^2 
    328         acc = vz_1              //      3.956E5 //!     CONVERT WAVELENGTH TO VELOCITY CM/SEC 
    329         SDD = L2                //1317 
    330         SSD = L1                //1627          //cm 
    331         lambda0 = lambda                //              15 
    332         DL_L = lambdaWidth              //0.236 
    333         SIG_L = DL_L/sqrt(6) 
    334         YG_d = -0.5*G*SDD*(SSD+SDD)*(LAMBDA0/acc)^2 
    335 //      Print "DISTANCE BEAM FALLS DUE TO GRAVITY (CM) = ",YG 
    336          
    337         sig_perp = kap*kap/12 * (3*(S1/L1)^2 + 3*(S2/LP)^2 + (proj_DDet/L2)^2) 
    338         sig_perp = sqrt(sig_perp) 
    339          
    340          
    341         FindQxQy(inQ,phi,qx,qy) 
    342          
    343         VAR_QLY = SIG_L^2 * (QY+4*PI*YG_d/(SDD*LAMBDA0))^2 
    344         VAR_QLX = (SIG_L*QX)^2 
    345         VAR_QL = VAR_QLY + VAR_QLX  //! WAVELENGTH CONTRIBUTION TO VARIANCE 
    346         sig_para = (sig_perp^2 + VAR_QL)^0.5 
    347          
    348 //return values PBR      
    349         SigmaQX = sig_para 
    350         SigmaQy = sig_perp 
    351          
     328//      acc = vz_1              //      3.956E5 //!     CONVERT WAVELENGTH TO VELOCITY CM/SEC 
     329//      SDD = L2                //1317 
     330//      SSD = L1                //1627          //cm 
     331//      lambda0 = lambda                //              15 
     332//      DL_L = lambdaWidth              //0.236 
     333//      SIG_L = DL_L/sqrt(6) 
     334//      YG_d = -0.5*G*SDD*(SSD+SDD)*(LAMBDA0/acc)^2 
     335/////// Print "DISTANCE BEAM FALLS DUE TO GRAVITY (CM) = ",YG 
    352336//       
     337//      sig_perp = kap*kap/12 * (3*(S1/L1)^2 + 3*(S2/LP)^2 + (proj_DDet/L2)^2) 
     338//      sig_perp = sqrt(sig_perp) 
     339//       
     340//       
     341//      FindQxQy(inQ,phi,qx,qy) 
     342//       
     343//      VAR_QLY = SIG_L^2 * (QY+4*PI*YG_d/(SDD*LAMBDA0))^2 
     344//      VAR_QLX = (SIG_L*QX)^2 
     345//      VAR_QL = VAR_QLY + VAR_QLX  //! WAVELENGTH CONTRIBUTION TO VARIANCE 
     346//      sig_para = (sig_perp^2 + VAR_QL)^0.5 
     347//       
     348/////// return values PBR        
     349//      SigmaQX = sig_para 
     350//      SigmaQy = sig_perp 
     351//       
     352//////   
    353353         
    354354        results = "success" 
  • sans/Dev/trunk/NCNR_User_Procedures/Reduction/SANS/NSORT.ipf

    r741 r811  
    15351535                        //write to notebook that file was not found 
    15361536                        //if string is not a number, report the error 
    1537                         if(str2num(partialName) == NaN) 
     1537                        if(numtype(str2num(partialName)) == 2) 
    15381538                                str = "this file was not found: "+partialName+"\r\r" 
    15391539                                //Notebook CatWin,font="Times",fsize=12,text=str 
  • sans/Dev/trunk/NCNR_User_Procedures/Reduction/SANS/SASCALC.ipf

    r800 r811  
    294294        // real wave 
    295295        rw = 0 
    296         rw[16] = 64             // beamcenter X (pixels) 
     296        rw[16] = 65             // beamcenter X (pixels) 
    297297        rw[17] = 64             // beamcenter Y 
    298298         
     
    754754                NVAR do2D = root:Packages:NIST:SAS:gDoMonteCarlo 
    755755                 
    756                 if(CB_Struct.eventMod == 2)             //if the shift key is down - go to 2D mode 
     756                if((CB_Struct.eventMod & 2^1) != 0 )            //if the shift key is down (bit 1)- go to 2D mode 
    757757                        do2D = 1 
    758758                endif 
     
    946946                if(doMonteCarlo == 1) 
    947947                        //2D simulation (in MultiScatter_MonteCarlo_2D.ipf) 
     948                         
     949                        // may want to take this back out--- 20 JUN 2011 SRK 
     950                        // this call to detectorOffset() sets the (xCtr,yCtr) in RealsRead to default values of (64.5,64.5) 
     951                        // -- if user changes the values in RealsRead to be able to do a proper average, the next sim starts wtih the 
     952                        // wrong values... 
     953                         
     954                        detectorOffset() 
    948955                         
    949956                        Simulate_2D_MC(funcStr,aveint,qval,sigave,sigmaq,qbar,fsubs) 
     
    17991806        rw[16] = 64  + round(2*rw[19]) + 0.5            //approximate beam X is 64 w/no offset, 114 w/25 cm offset  
    18001807        rw[17] = 64     + 0.5 //typical value 
    1801          
     1808 
    18021809        return(val) 
    18031810end 
Note: See TracChangeset for help on using the changeset viewer.