source: sans/SANSReduction/trunk/Put in User Procedures/SANS_Reduction_v5.00/WriteQIS.ipf @ 41

Last change on this file since 41 was 41, checked in by srkline, 16 years ago

change to UNIX line endings

File size: 23.9 KB
Line 
1#pragma rtGlobals=1             // Use modern global access method.
2#pragma version=5.0
3#pragma IgorVersion=4.0
4
5//************************
6// Vers 1.2 091001
7//
8//************************
9
10
11//***this function is not used - WriteWaves_W_Protocol() is used instead
12//
13//for writing out data (q-i-s) from the "type" folder
14//if fullpath is a complete HD path:filename, no dialog will be presented
15//if fullpath is just a filename, the save dialog will be presented (forced if dialog =1)
16//writes ONLY the standard header information (no protocol information)
17//
18Function WriteWaves(type,fullpath,dialog)
19        String type,fullpath
20        Variable dialog         //=1 will present dialog for name
21       
22        String destStr=""
23        destStr = "root:"+type
24       
25        Variable refNum
26        String formatStr = "%15.4g %15.4g %15.4g %15.4g %15.4g %15.4g\r\n"
27        String fname,ave="C",headerFormat = "%10.4g %8.2g %8.2g %8.2g %8.3g %8.3g %8s %5.0g\r\n"
28        Variable step=1
29       
30        //*****these waves MUST EXIST, or IGOR Pro will crash, with a type 2 error****
31        WAVE intw=$(destStr + ":integersRead")
32        WAVE rw=$(destStr + ":realsRead")
33        WAVE/T textw=$(destStr + ":textRead")
34        WAVE qvals =$(destStr + ":qval")
35        WAVE inten=$(destStr + ":aveint")
36        WAVE sig=$(destStr + ":sigave")
37        WAVE qbar = $(destStr + ":QBar")
38        WAVE sigmaq = $(destStr + ":SigmaQ")
39        WAVE fsubs = $(destStr + ":fSubS")
40
41        //check each wave
42        If(!(WaveExists(intw)))
43                Abort "intw DNExist BinaryWrite()"
44        Endif
45        If(!(WaveExists(rw)))
46                Abort "rw DNExist BinaryWrite()"
47        Endif
48        If(!(WaveExists(textw)))
49                Abort "textw DNExist BinaryWrite()"
50        Endif
51        If(!(WaveExists(qvals)))
52                Abort "qvals DNExist BinaryWrite()"
53        Endif
54        If(!(WaveExists(inten)))
55                Abort "inten DNExist BinaryWrite()"
56        Endif
57        If(!(WaveExists(sig)))
58                Abort "sig DNExist BinaryWrite()"
59        Endif
60        If(!(WaveExists(qbar)))
61                Abort "qbar DNExist BinaryWrite()"
62        Endif
63        If(!(WaveExists(sigmaq)))
64                Abort "sigmaq DNExist BinaryWrite()"
65        Endif
66        If(!(WaveExists(fsubs)))
67                Abort "fsubs DNExist BinaryWrite()"
68        Endif
69               
70        if(dialog)
71                PathInfo/S catPathName
72                fullPath = DoSaveFileDialog("Save data as")
73                If(cmpstr(fullPath,"")==0)
74                        //user cancel, don't write out a file
75                        Close/A
76                        Abort "no data file was written"
77                Endif
78                //Print "dialog fullpath = ",fullpath
79        Endif
80       
81        //actually open the file
82        Open refNum as fullpath
83       
84        fprintf refnum,"FILE: %s\t\t CREATED: %s\r\n",textw[0],textw[1]
85        fprintf refnum,"LABEL: %s\r\n",textw[6]
86        fprintf refnum,"MON CNT   LAMBDA   DET ANG   DET DIST   TRANS   THICK   AVE   STEP\r\n"
87        fprintf refnum,headerFormat,rw[0],rw[26],rw[19],rw[18],rw[4],rw[5],ave,step
88        wfprintf refnum, formatStr, qvals,inten,sig,sigmaq,qbar,fsubs
89       
90        Close refnum
91       
92        SetDataFolder root:             //(redundant)
93       
94        Return(0)
95End
96
97
98//for writing out data (q-i-s) from the "type" folder, and including reduction information
99//if fullpath is a complete HD path:filename, no dialog will be presented
100//if fullpath is just a filename, the save dialog will be presented
101//if dialog = 1, a dialog will always be presented
102//
103// root:myGlobals:Protocols:gProtoStr is the name of the currently active protocol
104//
105Function WriteWaves_W_Protocol(type,fullpath,dialog)
106        String type,fullpath
107        Variable dialog         //=1 will present dialog for name
108       
109        String destStr=""
110        destStr = "root:"+type
111       
112        Variable refNum
113        String formatStr = "%15.4g %15.4g %15.4g %15.4g %15.4g %15.4g\r\n"
114        String fname,ave="C",hdrStr1="",hdrStr2=""
115        Variable step=1
116       
117       
118       
119        //*****these waves MUST EXIST, or IGOR Pro will crash, with a type 2 error****
120        WAVE intw=$(destStr + ":integersRead")
121        WAVE rw=$(destStr + ":realsRead")
122        WAVE/T textw=$(destStr + ":textRead")
123        WAVE qvals =$(destStr + ":qval")
124        WAVE inten=$(destStr + ":aveint")
125        WAVE sig=$(destStr + ":sigave")
126        WAVE qbar = $(destStr + ":QBar")
127        WAVE sigmaq = $(destStr + ":SigmaQ")
128        WAVE fsubs = $(destStr + ":fSubS")
129
130        SVAR gProtoStr = root:myGlobals:Protocols:gProtoStr
131        Wave/T proto=$("root:myGlobals:Protocols:"+gProtoStr)
132       
133        //check each wave
134        If(!(WaveExists(intw)))
135                Abort "intw DNExist BinaryWrite_W_Protocol()"
136        Endif
137        If(!(WaveExists(rw)))
138                Abort "rw DNExist BinaryWrite_W_Protocol()"
139        Endif
140        If(!(WaveExists(textw)))
141                Abort "textw DNExist BinaryWrite_W_Protocol()"
142        Endif
143        If(!(WaveExists(qvals)))
144                Abort "qvals DNExist BinaryWrite_W_Protocol()"
145        Endif
146        If(!(WaveExists(inten)))
147                Abort "inten DNExist BinaryWrite_W_Protocol()"
148        Endif
149        If(!(WaveExists(sig)))
150                Abort "sig DNExist BinaryWrite_W_Protocol()"
151        Endif
152        If(!(WaveExists(qbar)))
153                Abort "qbar DNExist BinaryWrite_W_Protocol()"
154        Endif
155        If(!(WaveExists(sigmaq)))
156                Abort "sigmaq DNExist BinaryWrite_W_Protocol()"
157        Endif
158        If(!(WaveExists(fsubs)))
159                Abort "fsubs DNExist BinaryWrite_W_Protocol()"
160        Endif
161        If(!(WaveExists(proto)))
162                Abort "current protocol wave DNExist BinaryWrite_W_Protocol()"
163        Endif
164
165        //strings can be too long to print-- must trim to 255 chars
166        Variable ii,num=8
167        Make/O/T/N=(num) tempShortProto
168        for(ii=0;ii<num;ii+=1)
169                tempShortProto[ii] = (proto[ii])[0,240]
170        endfor
171       
172        if(dialog)
173                PathInfo/S catPathName
174                fullPath = DoSaveFileDialog("Save data as")
175                If(cmpstr(fullPath,"")==0)
176                        //user cancel, don't write out a file
177                        Close/A
178                        Abort "no data file was written"
179                Endif
180                //Print "dialog fullpath = ",fullpath
181        Endif
182       
183        hdrStr1 = num2str(rw[0])+"  "+num2str(rw[26])+"       "+num2str(rw[19])+"     "+num2str(rw[18])
184        hdrStr1 += "     "+num2str(rw[4])+"     "+num2str(rw[5]) + ave +"   "+num2str(step) + "\r\n"
185
186        hdrStr2 = num2str(rw[16])+"  "+num2str(rw[17])+"  "+num2str(rw[23])+"    "+num2str(rw[24])+"    "
187        hdrStr2 += num2str(rw[25])+"    "+num2str(rw[27])+"    "+num2str(rw[21])+"    "+textW[9] + "\r\n"
188       
189        SVAR samFiles = $("root:"+type+":fileList")
190        //actually open the file here
191        Open refNum as fullpath
192       
193        //write out the standard header information
194        fprintf refnum,"FILE: %s\t\t CREATED: %s\r\n",textw[0],textw[1]
195        fprintf refnum,"LABEL: %s\r\n",textw[6]
196        fprintf refnum,"MON CNT   LAMBDA   DET ANG   DET DIST   TRANS   THICK   AVE   STEP\r\n"
197        fprintf refnum,hdrStr1
198        fprintf refnum,"BCENT(X,Y)   A1(mm)   A2(mm)   A1A2DIST(m)   DL/L   BSTOP(mm)   DET_TYP \r\n"
199        fprintf refnum,hdrStr2
200//      fprintf refnum,headerFormat,rw[0],rw[26],rw[19],rw[18],rw[4],rw[5],ave,step
201
202        //insert protocol information here
203        //-1 list of sample files
204        //0 - bkg
205        //1 - emp
206        //2 - div
207        //3 - mask
208        //4 - abs params c2-c5
209        //5 - average params
210        fprintf refnum, "SAM: %s\r\n",samFiles
211        fprintf refnum, "BGD: %s\r\n",tempShortProto[0]
212        fprintf refnum, "EMP: %s\r\n",tempShortProto[1]
213        fprintf refnum, "DIV: %s\r\n",tempShortProto[2]
214        fprintf refnum, "MASK: %s\r\n",tempShortProto[3]
215        fprintf refnum, "ABS Parameters (3-6): %s\r\n",tempShortProto[4]
216        fprintf refnum, "Average Choices: %s\r\n",tempShortProto[5]
217       
218        //write out the data columns
219        fprintf refnum,"The 6 columns are | Q (1/A) | I(Q) (1/cm) | std. dev. I(Q) (1/cm) | sigmaQ | meanQ | ShadowFactor|\r\n"
220        wfprintf refnum, formatStr, qvals,inten,sig,sigmaq,qbar,fsubs
221       
222        Close refnum
223       
224        SetDataFolder root:             //(redundant)
225       
226        //write confirmation of write operation to history area
227        Print "Averaged File written: ", GetFileNameFromPathNoSemi(fullPath)
228        KillWaves/Z tempShortProto
229        Return(0)
230End
231
232
233//for writing out data (phi-i-s) from the "type" folder, and including reduction information
234//if fullpath is a complete HD path:filename, no dialog will be presented
235//if fullpath is just a filename, the save dialog will be presented
236//if dialog = 1, a dialog will always be presented
237//
238// root:myGlobals:Protocols:gProtoStr is the name of the currently active protocol
239//
240Function WritePhiave_W_Protocol(type,fullpath,dialog)
241        String type,fullpath
242        Variable dialog         //=1 will present dialog for name
243       
244        String destStr
245        destStr = "root:"+type
246       
247        Variable refNum
248        String formatStr = "%15.4g %15.4g %15.4g\r\n"
249        String fname,ave="C",hdrStr1,hdrStr2
250        Variable step=1
251       
252        //*****these waves MUST EXIST, or IGOR Pro will crash, with a type 2 error****
253        WAVE intw=$(destStr + ":integersRead")
254        WAVE rw=$(destStr + ":realsRead")
255        WAVE/T textw=$(destStr + ":textRead")
256        WAVE phival =$(destStr + ":phival")
257        WAVE inten=$(destStr + ":aveint")
258        WAVE sig=$(destStr + ":sigave")
259        SVAR gProtoStr = root:myGlobals:Protocols:gProtoStr
260        Wave/T proto=$("root:myGlobals:Protocols:"+gProtoStr)
261       
262        //check each wave
263        If(!(WaveExists(intw)))
264                Abort "intw DNExist BinaryWrite_W_Protocol()"
265        Endif
266        If(!(WaveExists(rw)))
267                Abort "rw DNExist BinaryWrite_W_Protocol()"
268        Endif
269        If(!(WaveExists(textw)))
270                Abort "textw DNExist BinaryWrite_W_Protocol()"
271        Endif
272        If(!(WaveExists(phival)))
273                Abort "qvals DNExist BinaryWrite_W_Protocol()"
274        Endif
275        If(!(WaveExists(inten)))
276                Abort "inten DNExist BinaryWrite_W_Protocol()"
277        Endif
278        If(!(WaveExists(sig)))
279                Abort "sig DNExist BinaryWrite_W_Protocol()"
280        Endif
281        If(!(WaveExists(proto)))
282                Abort "current protocol wave DNExist BinaryWrite_W_Protocol()"
283        Endif
284        //strings can be too long to print-- must trim to 255 chars
285        Variable ii,num=8
286        Make/O/T/N=(num) tempShortProto
287        for(ii=0;ii<num;ii+=1)
288                tempShortProto[ii] = (proto[ii])[0,240]
289        endfor
290       
291        if(dialog)
292                PathInfo/S catPathName
293                fullPath = DoSaveFileDialog("Save data as")
294                If(cmpstr(fullPath,"")==0)
295                        //user cancel, don't write out a file
296                        Close/A
297                        Abort "no data file was written"
298                Endif
299                //Print "dialog fullpath = ",fullpath
300        Endif
301       
302        hdrStr1 = num2str(rw[0])+"  "+num2str(rw[26])+"       "+num2str(rw[19])+"     "+num2str(rw[18])
303        hdrStr1 += "     "+num2str(rw[4])+"     "+num2str(rw[5]) + ave +"   "+num2str(step) + "\r\n"
304
305        hdrStr2 = num2str(rw[16])+"  "+num2str(rw[17])+"  "+num2str(rw[23])+"    "+num2str(rw[24])+"    "
306        hdrStr2 += num2str(rw[25])+"    "+num2str(rw[27])+"    "+num2str(rw[21])+"    "+textW[9] + "\r\n"
307       
308        SVAR samFiles = $("root:"+type+":fileList")
309        //actually open the file here
310        Open refNum as fullpath
311       
312        //write out the standard header information
313        fprintf refnum,"FILE: %s\t\t CREATED: %s\r\n",textw[0],textw[1]
314        fprintf refnum,"LABEL: %s\r\n",textw[6]
315        fprintf refnum,"MON CNT   LAMBDA   DET ANG   DET DIST   TRANS   THICK   AVE   STEP\r\n"
316        fprintf refnum,hdrStr1
317        fprintf refnum,"BCENT(X,Y)   A1(mm)   A2(mm)   A1A2DIST(m)   DL/L   BSTOP(mm)   DET_TYP \r\n"
318        fprintf refnum,hdrStr2
319       
320        //insert protocol information here
321        //0 - bkg
322        //1 - emp
323        //2 - div
324        //3 - mask
325        //4 - abs params c2-c5
326        //5 - average params
327        fprintf refnum, "SAM: %s\r\n",samFiles
328        fprintf refnum, "BGD: %s\r\n",tempShortProto[0]
329        fprintf refnum, "EMP: %s\r\n",tempShortProto[1]
330        fprintf refnum, "DIV: %s\r\n",tempShortProto[2]
331        fprintf refnum, "MASK: %s\r\n",tempShortProto[3]
332        fprintf refnum, "ABS Parameters (3-6): %s\r\n",tempShortProto[4]
333        fprintf refnum, "Average Choices: %s\r\n",tempShortProto[5]
334       
335        //write out the data columns
336        fprintf refnum,"The 3 columns are | Phi (deg) | I(phi) (1/cm) | std. dev. I(phi) (1/cm) |\r\n"
337        wfprintf refnum, formatStr, phival,inten,sig
338       
339        Close refnum
340       
341        SetDataFolder root:             //(redundant)
342       
343        //write confirmation of write operation to history area
344        Print "Averaged File written: ", GetFileNameFromPathNoSemi(fullPath)
345        KillWaves/Z tempShortProto
346
347        Return(0)
348End
349
350//*****************
351// saves the data after all of the desired reduction steps (average options)
352// as a 2x expanded PNG file (approx 33kb)
353//
354Function SaveAsPNG(type,fullPath,dialog)
355        String type,fullPath
356        Variable dialog
357       
358        Variable refnum
359        if(dialog)
360                PathInfo/S catPathName
361                Open/D refnum as fullpath               //won't actually open the file
362                If(cmpstr(S_filename,"")==0)
363                        //user cancel, don't write out a file
364                        Close/A
365                        Abort "no data file was written"
366                Endif
367                fullpath = S_filename
368                //Print "dialog fullpath = ",fullpath
369        Endif
370       
371        //cleanup the filename passed in from Protocol...
372        String oldStr="",newStr="",pathStr=""
373        oldStr=GetFileNameFromPathNoSemi(fullPath)      //just the filename
374        pathStr=GetPathStrFromfullName(fullPath)        //just the path
375       
376        newStr = CleanupName(oldStr, 0 )                                //filename with _EXT rather than .EXT
377        fullPath=pathStr+newStr+".png"                          //tack on the png extension
378       
379        print "type=",type
380        //graph the current data and save a little graph
381        Wave data =  $("root:"+type+":data")
382        Wave q_x_axis = $"root:myGlobals:q_x_axis"
383        Wave q_y_axis = $"root:myGlobals:q_y_axis"
384        Wave NIHColors = $"root:myGlobals:NIHColors"
385       
386        NewImage/F data
387        DoWindow/C temp_PNG
388        ModifyImage data cindex= NIHColors
389        AppendToGraph/R q_y_axis
390        ModifyGraph tkLblRot(right)=90,lowTrip(right)=0.001
391        AppendToGraph/T q_x_axis
392        ModifyGraph lowTrip(top)=0.001,standoff=0,mode=2
393        ModifyGraph fSize(right)=9,fSize(top)=9,btLen=3
394       
395//      ModifyGraph nticks=0
396       
397//      WaveStats/Q data
398//      ScaleColorsToData(V_min, V_max, NIHColors)
399
400// ***comment out for DEMO_MODIFIED version
401        SavePict/Z/E=-5/B=144 as fullPath                       //PNG at 2x screen resolution
402//***
403
404        Print "Saved graphic as ",newStr+".png"
405        DoWindow/K temp_PNG
406End
407
408//****************
409//Testing only , not called
410Proc Fast_ASCII_2D_Export(type,term)
411        String type,term
412        Prompt type,"2-D data type for Export",popup,"SAM;EMP;BGD;DIV;COR;CAL;RAW;ABS;MSK;"
413        Prompt term,"line termination",popup,"CR;LF;CRLF;"
414       
415        //terminator is currently ignored
416        Fast2dExport(type,"",1)
417       
418End
419
420//the default termination for the platform is used...
421//if RAW export, sets "dummy" protocol to "RAW data export"
422Function Fast2dExport(type,fullpath,dialog)
423        String type,fullpath
424        Variable dialog         //=1 will present dialog for name
425       
426        String destStr="",ave="C",typeStr=""
427        Variable step=1,refnum
428        destStr = "root:"+type
429       
430        //must select the linear_data to export
431        // can't export log data if there are -ve intensities from a subtraction
432        NVAR isLog = $(destStr+":gIsLogScale")
433        if(isLog==1)
434                typeStr = ":linear_data"
435        else
436                typeStr = ":data"
437        endif
438       
439        NVAR pixelsX = root:myGlobals:gNPixelsX
440        NVAR pixelsY = root:myGlobals:gNPixelsY
441       
442        Wave data=$(destStr+typeStr)
443        WAVE intw=$(destStr + ":integersRead")
444        WAVE rw=$(destStr + ":realsRead")
445        WAVE/T textw=$(destStr + ":textRead")
446
447        SVAR gProtoStr = root:myGlobals:Protocols:gProtoStr
448        String rawTag=""
449        if(cmpstr(type,"RAW")==0)
450                Make/O/T/N=8 proto={"none","none","none","none","none","none","none","none"}
451                RawTag = "RAW Data File: "     
452        else
453                Wave/T proto=$("root:myGlobals:Protocols:"+gProtoStr)
454        endif
455        SVAR samFiles = $("root:"+type+":fileList")
456        //check each wave - MUST exist, or will cause a crash
457        If(!(WaveExists(data)))
458                Abort "data DNExist AsciiExport()"
459        Endif
460        If(!(WaveExists(intw)))
461                Abort "intw DNExist AsciiExport()"
462        Endif
463        If(!(WaveExists(rw)))
464                Abort "rw DNExist AsciiExport()"
465        Endif
466        If(!(WaveExists(textw)))
467                Abort "textw DNExist AsciiExport()"
468        Endif
469        If(!(WaveExists(proto)))
470                Abort "current protocol wave DNExist AsciiExport()"
471        Endif
472       
473        if(dialog)
474                PathInfo/S catPathName
475                fullPath = DoSaveFileDialog("Save data as")
476                If(cmpstr(fullPath,"")==0)
477                        //user cancel, don't write out a file
478                        Close/A
479                        Abort "no data file was written"
480                Endif
481                //Print "dialog fullpath = ",fullpath
482        Endif
483       
484/////////
485        Variable numTextLines=18
486        Make/O/T/N=(numTextLines) labelWave
487        labelWave[0] = "FILE: "+textw[0]+"   CREATED: "+textw[1]
488        labelWave[1] = "LABEL: "+textw[6]
489        labelWave[2] = "MON CNT   LAMBDA(A)   DET_OFF(cm)   DET_DIST(m)   TRANS   THICK(cm)"
490        labelWave[3] = num2str(rw[0])+"  "+num2str(rw[26])+"       "+num2str(rw[19])+"     "+num2str(rw[18])
491        labelWave[3] += "     "+num2str(rw[4])+"     "+num2str(rw[5])
492        labelWave[4] = "BCENT(X,Y)   A1(mm)   A2(mm)   A1A2DIST(m)   DL/L   BSTOP(mm)   DET_TYP  "
493        labelWave[5] = num2str(rw[16])+"  "+num2str(rw[17])+"  "+num2str(rw[23])+"  "+num2str(rw[24])+"  "
494        labelWave[5] += num2str(rw[25])+"  "+num2str(rw[27])+"  "+num2str(rw[21])+"  "+textW[9]
495        labelWave[6] =  "SAM: "+rawTag+samFiles
496        labelWave[7] =  "BGD: "+proto[0]
497        labelWave[8] =  "EMP: "+proto[1]
498        labelWave[9] =  "DIV: "+proto[2]
499        labelWave[10] =  "MASK: "+proto[3]
500        labelWave[11] =  "ABS Parameters (3-6): "+proto[4]
501        labelWave[12] = "Average Choices: "+proto[5]
502        labelWave[13] = ""
503        labelWave[14] = "*** Data written from "+type+" folder and may not be a fully corrected data file ***"
504        labelWave[15] = "The detector image is a standard X-Y coordinate system"
505        labelWave[16] = "Data is written by row, starting with Y=1 and X=(1->128)"
506        labelWave[17] = "ASCII data created " +date()+" "+time()
507        //strings can be too long to print-- must trim to 255 chars
508        Variable ii
509        for(ii=0;ii<numTextLines;ii+=1)
510                labelWave[ii] = (labelWave[ii])[0,240]
511        endfor
512//      If(cmpstr(term,"CR")==0)
513//              termStr = "\r"
514//      Endif
515//      If(cmpstr(term,"LF")==0)
516//              termStr = "\n"
517//      Endif
518//      If(cmpstr(term,"CRLF")==0)
519//              termStr = "\r\n"
520//      Endif
521       
522        Duplicate/O data,spWave         
523        Redimension/S/N=(pixelsX*pixelsY) spWave                //single precision (/S)
524       
525        //not demo- compatible, but approx 6x faster!!
526//      Save/G/M="\r\n" labelWave,spWave as fullPath            // /M=termStr specifies terminator
527        ///
528       
529        Open refNum as fullpath
530        wfprintf refNum,"%s\r\n",labelWave
531        fprintf refnum,"\r\n"
532        wfprintf refNum,"%g\r\n",spWave
533        Close refNum
534       
535        Killwaves/Z spWave,labelWave            //don't delete proto!
536       
537        Print "2D ASCII File written: ", GetFileNameFromPathNoSemi(fullPath)
538        return(0)
539End
540
541//**********************
542// Resolution calculation - used by the averaging routines
543// to calculate the resolution function at each q-value
544// - the return value is not used
545//
546// equivalent to John's routine on the VAX Q_SIGMA_AVE.FOR
547// Incorporates eqn. 3-15 from J. Appl. Cryst. (1995) v. 28 p105-114
548//
549Function/S getResolution(inQ,lambda,lambdaWidth,DDet,apOff,S1,S2,L1,L2,BS,del_r,SigmaQ,QBar,fSubS)
550        Variable inQ, lambda, lambdaWidth, DDet, apOff, S1, S2, L1, L2, BS, del_r
551        Variable &fSubS, &QBar, &SigmaQ         //these are the output quantities at the input Q value
552       
553        //lots of calculation variables
554        Variable a2, q_small, lp, v_lambda, v_b, v_d, vz, yg, v_g
555        Variable r0, delta, inc_gamma, fr, fv, rmd, v_r1, rm, v_r
556
557        //Constants
558        //Variable del_r = .1
559        Variable vz_1 = 3.956e5         //velocity [cm/s] of 1 A neutron
560        Variable g = 981.0                              //gravity acceleration [cm/s^2]
561
562        String results
563        results ="Failure"
564
565        //rename for working variables,  these must be gotten from global
566        //variables
567
568//      Variable wLam, wLW, wL1, wL2, wS1, wS2
569//      Variable wBS, wDDet, wApOff
570//      wLam = lambda
571//      wLW = lambdaWidth
572//      wDDet = DDet
573//      wApOff = apOff
574        S1 *= 0.5*0.1                   //convert to radius and [cm]
575        S2 *= 0.5*0.1
576
577        L1 *= 100.0                     // [cm]
578        L1 -= apOff                             //correct the distance
579
580        L2 *= 100.0
581        L2 += apOff
582
583        BS *= 0.5*0.1                   //convert to radius and [cm]
584        del_r *= 0.1                            //width of annulus, convert mm to [cm]
585       
586        //Start resolution calculation
587        a2 = S1*L2/L1 + S2*(L1+L2)/L1
588        q_small = 2.0*Pi*(BS-a2)*(1.0-lambdaWidth)/(lambda*L2)
589        lp = 1.0/( 1.0/L1 + 1.0/L2)
590
591        v_lambda = lambdaWidth^2/6.0
592        v_b = 0.25*(S1*L2/L1)^2 +0.25*(S2*L2/lp)^2
593        v_d = (DDet/2.3548)^2 + del_r^2/12.0
594        vz = vz_1 / lambda
595        yg = 0.5*g*L2*(L1+L2)/vz^2
596        v_g = 2.0*yg^2*v_lambda
597
598        r0 = L2*tan(2.0*asin(lambda*inQ/(4.0*Pi) ))
599        delta = 0.5*(BS - r0)^2/v_d
600
601        if (r0 < BS)
602                inc_gamma=exp(gammln(1.5))*(1-gammp(1.5,delta))
603        else
604                inc_gamma=exp(gammln(1.5))*(1+gammp(1.5,delta))
605        endif
606
607        fSubS = 0.5*(1.0+erf( (r0-BS)/sqrt(2.0*v_d) ) )
608        if (fSubS <= 0.0)
609                fSubS = 1.e-10
610        endif
611        fr = 1.0 + sqrt(v_d)*exp(-1.0*delta) /(r0*fSubS*sqrt(2.0*Pi))
612        fv = inc_gamma/(fSubS*sqrt(Pi)) - r0^2*(fr-1.0)^2/v_d
613
614        rmd = fr*r0
615        v_r1 = v_b + fv*v_d +v_g
616
617        rm = rmd + 0.5*v_r1/rmd
618        v_r = v_r1 - 0.5*(v_r1/rmd)^2
619        if (v_r < 0.0)
620                v_r = 0.0
621        endif
622        QBar = (4.0*Pi/lambda)*sin(0.5*atan(rm/L2))
623        SigmaQ = QBar*sqrt(v_r/rmd^2 +v_lambda)
624
625        results = "success"
626        Return results
627End
628
629//ASCII export of data as 3-columns qx-qy-Intensity
630//limited header information?
631//
632// - creates the qx and qy data here, based on the data and header information
633//
634// Need to ensure that the data being exported is the linear copy of the dataset - check the global
635//
636Function QxQy_Export(type,fullpath,dialog)
637        String type,fullpath
638        Variable dialog         //=1 will present dialog for name
639       
640        String destStr="",typeStr=""
641        Variable step=1,refnum
642        destStr = "root:"+type
643       
644        //must select the linear_data to export
645        NVAR isLog = $(destStr+":gIsLogScale")
646        if(isLog==1)
647                typeStr = ":linear_data"
648        else
649                typeStr = ":data"
650        endif
651       
652        NVAR pixelsX = root:myGlobals:gNPixelsX
653        NVAR pixelsY = root:myGlobals:gNPixelsY
654       
655        Wave data=$(destStr+typeStr)
656        WAVE intw=$(destStr + ":integersRead")
657        WAVE rw=$(destStr + ":realsRead")
658        WAVE/T textw=$(destStr + ":textRead")
659
660        SVAR gProtoStr = root:myGlobals:Protocols:gProtoStr
661        String rawTag=""
662        if(cmpstr(type,"RAW")==0)
663                Make/O/T/N=8 proto={"none","none","none","none","none","none","none","none"}
664                RawTag = "RAW Data File: "     
665        else
666                Wave/T proto=$("root:myGlobals:Protocols:"+gProtoStr)
667        endif
668        SVAR samFiles = $("root:"+type+":fileList")
669        //check each wave - MUST exist, or will cause a crash
670        If(!(WaveExists(data)))
671                Abort "data DNExist QxQy_Export()"
672        Endif
673        If(!(WaveExists(intw)))
674                Abort "intw DNExist QxQy_Export()"
675        Endif
676        If(!(WaveExists(rw)))
677                Abort "rw DNExist QxQy_Export()"
678        Endif
679        If(!(WaveExists(textw)))
680                Abort "textw DNExist QxQy_Export()"
681        Endif
682        If(!(WaveExists(proto)))
683                Abort "current protocol wave DNExist QxQy_Export()"
684        Endif
685       
686        if(dialog)
687                PathInfo/S catPathName
688                fullPath = DoSaveFileDialog("Save data as")
689                If(cmpstr(fullPath,"")==0)
690                        //user cancel, don't write out a file
691                        Close/A
692                        Abort "no data file was written"
693                Endif
694                //Print "dialog fullpath = ",fullpath
695        Endif
696       
697/////////
698        Variable numTextLines=18
699        Make/O/T/N=(numTextLines) labelWave
700        labelWave[0] = "FILE: "+textw[0]+"   CREATED: "+textw[1]
701        labelWave[1] = "LABEL: "+textw[6]
702        labelWave[2] = "MON CNT   LAMBDA (A)  DET_OFF(cm)   DET_DIST(m)   TRANS   THICK(cm)"
703        labelWave[3] = num2str(rw[0])+"  "+num2str(rw[26])+"       "+num2str(rw[19])+"     "+num2str(rw[18])
704        labelWave[3] += "     "+num2str(rw[4])+"     "+num2str(rw[5])
705        labelWave[4] = "BCENT(X,Y)   A1(mm)   A2(mm)   A1A2DIST(m)   DL/L   BSTOP(mm)   DET_TYP  "
706        labelWave[5] = num2str(rw[16])+"  "+num2str(rw[17])+"  "+num2str(rw[23])+"    "+num2str(rw[24])+"    "
707        labelWave[5] += num2str(rw[25])+"    "+num2str(rw[27])+"    "+num2str(rw[21])+"    "+textW[9]
708        labelWave[6] =  "SAM: "+rawTag+samFiles
709        labelWave[7] =  "BGD: "+proto[0]
710        labelWave[8] =  "EMP: "+proto[1]
711        labelWave[9] =  "DIV: "+proto[2]
712        labelWave[10] =  "MASK: "+proto[3]
713        labelWave[11] =  "ABS Parameters (3-6): "+proto[4]
714        labelWave[12] = "Average Choices: "+proto[5]
715        labelWave[13] = ""
716        labelWave[14] = "*** Data written from "+type+" folder and may not be a fully corrected data file ***"
717        labelWave[15] = "Data columns are Qx - Qy - I(Qx,Qy)"
718        labelWave[16] = ""
719        labelWave[17] = "ASCII data created " +date()+" "+time()
720        //strings can be too long to print-- must trim to 255 chars
721        Variable ii,jj
722        for(ii=0;ii<numTextLines;ii+=1)
723                labelWave[ii] = (labelWave[ii])[0,240]
724        endfor
725//      If(cmpstr(term,"CR")==0)
726//              termStr = "\r"
727//      Endif
728//      If(cmpstr(term,"LF")==0)
729//              termStr = "\n"
730//      Endif
731//      If(cmpstr(term,"CRLF")==0)
732//              termStr = "\r\n"
733//      Endif
734       
735        Duplicate/O data,qx_val,qy_val,z_val
736        Redimension/N=(pixelsX*pixelsY) qx_val,qy_val,z_val
737        MyMat2XYZ(data,qx_val,qy_val,z_val)             //x and y are [p][q] indexes, not q-vals yet
738       
739        qx_val = CalcQx(qx_val+1,rw[16],rw[18],rw[26])          //+1 converts to detector coordinate system
740        qy_val = CalcQy(qy_val+1,rw[17],rw[18],rw[26])
741
742        //not demo-compatible, but approx 8x faster!!           
743//      Save/G/M="\r\n" labelWave,qx_val,qy_val,z_val as fullpath       // /M=termStr specifies terminator
744        //
745       
746        Open refNum as fullpath
747        wfprintf refNum,"%s\r\n",labelWave
748        fprintf refnum,"\r\n"
749        wfprintf refNum,"%8g\t%8g\t%8g\r\n",qx_val,qy_val,z_val
750        Close refNum
751       
752        Killwaves/Z spWave,labelWave,qx_val,qy_val,z_val
753       
754        Print "QxQy_Export File written: ", GetFileNameFromPathNoSemi(fullPath)
755        return(0)
756
757End
758
759
760Function MyMat2XYZ(mat,xw,yw,zw)
761        WAVE mat,xw,yw,zw
762
763        NVAR pixelsX = root:myGlobals:gNPixelsX
764        NVAR pixelsY = root:myGlobals:gNPixelsY
765       
766        xw= mod(p,pixelsX)              // X varies quickly
767        yw= floor(p/pixelsY)    // Y varies slowly
768        zw= mat(xw[p])(yw[p])
769
770End
771
772//converts xyz triple to a matrix
773//MAJOR assumption is that the x and y-spacings are LINEAR
774// (ok for small-angle approximation)
775//
776// currently unused
777//
778Function LinXYZToMatrix(xw,yw,zw,matStr)
779        WAVE xw,yw,zw
780        String matStr
781       
782        NVAR pixelsX = root:myGlobals:gNPixelsX
783        NVAR pixelsY = root:myGlobals:gNPixelsY
784        //mat is "zw" redimensioned to a matrix
785        Make/O/N=(pixelsX*pixelsY) $matStr
786        WAVE mat=$matStr
787        mat=zw
788        Redimension/N=(pixelsX,pixelsY) mat
789        WaveStats/Q xw
790        SetScale/I x, V_min, V_max, "",mat
791        WaveStats/Q yw
792        SetScale/I y, V_min, V_max, "",mat
793       
794        Display;Appendimage mat
795        ModifyGraph lowTrip=0.0001
796        ModifyGraph width={Plan,1,bottom,left},height={Plan,1,left,bottom}
797        ModifyImage $matStr ctab={*,*,YellowHot,0}
798       
799        return(0)
800End
Note: See TracBrowser for help on using the repository browser.