source: sans/Dev/trunk/NCNR_User_Procedures/Reduction/USANS/BT5_Loader.ipf @ 734

Last change on this file since 734 was 570, checked in by srkline, 13 years ago

Change (1):
In preparation for release, updated pragma IgorVersion?=6.1 in all procedures

Change (2):
As a side benefit of requiring 6.1, we can use the MultiThread? keyword to thread any model function we like. The speed benefit is only noticeable on functions that require at least one integration and at least 100 points (resolution smearing is NOT threaded, too many threadSafe issues, too little benefit). I have chosen to use the MultiThread? only on the XOP assignment. In the Igor code there are too many functions that are not explicitly declared threadsafe, making for a mess.

  • Property eol-style set to native
  • Property svn:executable set to *
File size: 7.7 KB
Line 
1#pragma rtGlobals=1             // Use modern global access method.
2#pragma Version=2.20
3#pragma IgorVersion=6.1
4
5/////////////////////////
6// 101001 Vers. 1
7//
8// functions to load and parse the ASCII ICP files generated at BT5 USANS
9// Nearly all of the important information about the data is dragged around
10// with the DetCts wave as a wave note
11//
12// - thes wave note is a string of KEY:value items
13//
14//      str = "FILE:"+filen+";"
15//      str += "TIMEPT:"+num2str(counttime)+";"
16//      str += "PEAKANG:;"              //no value yet
17//      str += "STATUS:;"               //no value yet
18//      str += "LABEL:"+fileLabel+";"
19//      str += "PEAKVAL:;"              //no value yet
20//      str += "TWIDE:0;"               //no value yet
21//
22/////////////////////////
23
24
25//given the filename, loads a single data file into a "TYPE" folder, into a
26//series of 1D waves. counting time is associated with the DetCts wave as a note
27//"blank" note entries are set up here as well, so later they can be polled/updated
28Function LoadBT5File(fname,type)
29        String fname,type
30
31        SVAR USANSFolder = root:Packages:NIST:USANS:Globals:gUSANSFolder
32       
33        Variable num=200,err=0,refnum
34        Make/O/D/N=(num) $(USANSFolder+":"+type+":Angle")
35        Make/O/D/N=(num) $(USANSFolder+":"+type+":DetCts")
36        Make/O/D/N=(num) $(USANSFolder+":"+type+":ErrDetCts")
37        Make/O/D/N=(num) $(USANSFolder+":"+type+":MonCts")
38        Make/O/D/N=(num) $(USANSFolder+":"+type+":TransCts")
39        Wave Angle = $(USANSFolder+":"+type+":Angle")
40        Wave DetCts = $(USANSFolder+":"+type+":DetCts")
41        Wave ErrDetCts = $(USANSFolder+":"+type+":ErrDetCts")
42        Wave MonCts = $(USANSFolder+":"+type+":MonCts")
43        Wave TransCts = $(USANSFolder+":"+type+":TransCts")
44       
45        Open/R refNum as fname          //if fname is "", a dialog will be presented
46        if(refnum==0)
47                return(1)               //user cancelled
48        endif
49        //read in the ASCII data line-by-line
50        Variable numLinesLoaded = 0,firstchar,countTime
51        Variable v1,v2,v3,v4,v5,v6,v7,v8,v9,v10,v11,v12,v13,v14,v15,v16,ii,valuesRead
52        String buffer ="",s1,s2,s3,s4,s5,s6,s7,s8,s9,s10
53        String filen="",fileLabel=""
54       
55        //parse the first line
56        FReadLine refnum,buffer
57        sscanf buffer, "%s%s%s%s%s%s%g%g%s%g%s",s1,s2,s3,s4,s5,s6,v1,v2,s7,v3,s8
58        ii=strlen(s1)
59        filen=s1[1,(ii-2)]                      //remove the ' from beginning and end
60        //Print "filen = ",filen,strlen(filen)
61        //v2 is the monitor prefactor. multiply monitor time by prefactor. AJJ 5 March 07
62        countTime=v1*v2
63        //Print "time (sec) = ",countTime
64       
65        //skip line 2
66        FReadLine refnum,buffer
67        //the next line is the sample label, use it all, minus the terminator
68        FReadLine refnum,filelabel
69       
70        //skip the next 10 lines
71        For(ii=0;ii<10;ii+=1)
72                FReadLine refnum,buffer
73        EndFor
74       
75        //read the data until EOF - assuming always a pair or lines
76        do
77                FReadLine refNum, buffer
78                if (strlen(buffer) == 0)
79                        break                                                   // End of file
80                endif
81                firstChar = char2num(buffer[0])
82                if ( (firstChar==10) || (firstChar==13) )
83                        break                                                   // Hit blank line. End of data in the file.
84                endif
85                sscanf buffer,"%g%g%g%g%g",v1,v2,v3,v4,v5               // 5 values here now
86                angle[numlinesloaded] = v1              //[0] is the ANGLE
87               
88                FReadLine refNum,buffer //assume a 2nd line is there, w/16 values
89                sscanf buffer,"%g,%g,%g,%g,%g,%g,%g,%g,%g,%g,%g,%g,%g,%g,%g,%g",v1,v2,v3,v4,v5,v6,v7,v8,v9,v10,v11,v12,v13,v14,v15,v16
90                //valuesRead = V_flag
91                //print valuesRead
92                MonCts[numlinesloaded] = v1             //monitor
93                DetCts[numlinesloaded] = v2 + v3 + v5 + v6 + v7         //5 detectors
94                TransCts[numlinesloaded] = v4           //trans detector
95                ErrDetCts[numlinesloaded] = sqrt(detCts[numlinesloaded])
96                //values 8-16 are always zero
97                numlinesloaded += 1             //2 lines in file read, one "real" line of data
98        while(1)
99               
100        Close refNum            // Close the file.
101        //Print "numlines = ",numlinesloaded
102        //trim the waves to the correct number of points
103        Redimension/N=(numlinesloaded) Angle,MonCts,DetCts,TransCts,ErrDetCts
104       
105        //remove LF from end of filelabel
106        filelabel=fileLabel[0,(strlen(fileLabel)-2)]
107       
108        //set the wave note for the DetCts
109        String str=""
110        str = "FILE:"+filen+";"
111        str += "TIMEPT:"+num2str(counttime)+";"
112        str += "PEAKANG:;"              //no value yet
113        str += "STATUS:;"               //no value yet
114        str += "LABEL:"+fileLabel+";"
115        str += "PEAKVAL:;"              //no value yet
116        str += "TWIDE:0;"               //no value yet
117        Note DetCts,str
118       
119        String/G fileLbl = filelabel
120        return err                      // Zero signifies no error.     
121End
122
123
124//will search for peak in root:type:DetCts,  find and returns the angle
125// of the peak poistion (= q=0)
126//returns -9999 if error, appends value to wavenote
127// !can't return -1,0,1, since these may be the peak angle!
128//
129Function FindZeroAngle(type)
130        String type
131       
132        SVAR USANSFolder = root:Packages:NIST:USANS:Globals:gUSANSFolder       
133       
134        Variable pkNotFound,pkPt,pkAngle,pkVal,temp
135        Wave angle = $(USANSFolder+":"+type+":Angle")
136        Wave detCts = $(USANSFolder+":"+type+":DetCts")
137
138
139        WaveStats/Q detcts
140        temp=V_Maxloc           //in points
141        FindPeak/P/Q/M=10/R=[(temp-10),(temp+10)] DetCts                //+/- 10 pts from maximum
142       
143        //FindPeak/P/Q/M=(1000) detCts          // /M=min ht -- peak must be at least 100x higher than 1st pt
144        pkNotFound=V_Flag               //V_Flag==1 if no pk found
145        pkPt = V_PeakLoc
146        pkVal = V_PeakVal               //for calc of T_rock
147        if(pkNotFound)
148                //DoAlert 0, "Peak not found"
149                Return (-9999)          //fatal error
150        Endif
151        pkAngle = Angle[pkPt]
152        //Print "Peak Angle = ",pkAngle
153        //update the note
154        String str=""
155        str=note(DetCts)
156        str = ReplaceNumberByKey("PEAKANG",str,pkangle,":",";")
157        str = ReplaceNumberByKey("PEAKVAL",str,pkVal,":",";")
158        Note/K DetCts
159        Note detCts,str
160       
161        return(pkAngle)
162End
163
164// given the peakAngle (q=0) and the "type" folder of scattering angle,
165// convert the angle to Q-values [degrees to (1/A)]
166// makes a new Qvals wave, duplicating the Angle wave, which is assumed to exist
167// Uses a conversion constant supplied by John Barker, and is hard-wired in
168//
169Function ConvertAngle2Qvals(type,pkAngle)
170        String type
171        Variable pkAngle
172
173        SVAR USANSFolder = root:Packages:NIST:USANS:Globals:gUSANSFolder
174       
175        Wave angle = $(USANSFolder+":"+type+":Angle")
176        Variable num=numpnts(angle)
177        Variable deg2QConv=5.55e-5              //JGB -- 2/24/01
178       
179        Make/O/N=(num) $(USANSFolder+":"+type+":Qvals")
180        Wave qvals = $(USANSFolder+":"+type+":Qvals")   
181        Qvals = deg2QConv*(angle[p] - pkAngle)
182       
183        return(0)       //no error
184End
185
186//updates the wavenote with the average Trans det cts for calculation of T_Wide
187//  finds the average number of counts on the transmission detector at angles
188// greater than 2 deg. (per John, practical experience where the trans data levels off)
189//
190// error value of 1 is returned and wavenote not updated if level is not found
191//
192// now normalizes to the monitor counts
193//
194Function FindTWideCts(type)
195        String type
196       
197        SVAR USANSFolder = root:Packages:NIST:USANS:Globals:gUSANSFolder
198       
199        Variable levNotFound,levPt,Cts,num,ii,sumMonCts
200       
201        Wave angle = $(USANSFolder+":"+type+":Angle")
202        Wave detCts = $(USANSFolder+":"+type+":DetCts")
203        Wave TransCts = $(USANSFolder+":"+type+":TransCts")
204        Wave monCts = $(USANSFolder+":"+type+":MonCts")
205       
206        num=numpnts(TransCts)
207
208
209        FindLevel/Q/P angle,2           //use angles greater than 2 deg
210        levNotFound=V_Flag              //V_Flag==1 if no pk found
211       
212        if(levNotFound)
213                //post a warning, and use just the last 4 points...
214                DoAlert 0,"You don't have "+type+" data past 2 degrees - so Twide may not be reliable"
215                levPt = num-4
216        else
217                levPt = trunc(V_LevelX)         // in points, force to integer
218        endif
219       
220        //average the trans cts from levPt to the end of the dataset
221        Cts=0
222        sumMonCts=0
223       
224        for(ii=levPt;ii<num;ii+=1)
225                Cts += transCts[ii]
226                sumMonCts += monCts[ii]
227        endfor
228        // get the average over that number of data points
229        sumMonCts /= (num-levPt-1)
230        Cts /= (num-levPt-1)
231       
232//      Print "cts = ",cts
233//      Print "sumMoncts = ",sumMoncts
234       
235        //normalize to the average monitor counts
236        Cts /= sumMonCts
237       
238//      Print "normalized counts = ",cts
239
240
241        //update the note
242        Wave DetCts = $(USANSFolder+":"+type+":DetCts")
243        String str,strVal
244        str=note(DetCts)
245        str = ReplaceNumberByKey("TWIDE",str,Cts,":",";")
246        Note/K DetCts
247        Note detCts,str
248
249        return(0)
250End
Note: See TracBrowser for help on using the repository browser.