source: sans/USANSReduction/trunk/Put in User Procedures Folder/USANS Procedures v2.20/LakeDesmearing_JB.ipf @ 35

Last change on this file since 35 was 35, checked in by ajj, 16 years ago

Set pragma = 2.20 and required version >= 2.20 for includes in USANS_Includes.ipf

  • Property eol-style set to native
File size: 49.2 KB
Line 
1#pragma rtGlobals=1             // Use modern global access method.
2#pragma Version=2.20
3
4//////////////////////////////////////////////
5// Igor conversion:     03 FEB 06 SRK
6//                      Revision:       03 MAR 06 SRK
7//
8//      Program uses Lake's iterative technique, J. A. Lake, Acta Cryst. 23 (1967) 191-4.
9//      to DESMEAR Infinite slit smeared USANS DATA.
10//
11//      TO SEE DESCRIPTION, CHECK COMPUTER SOFTWARE LOGBOOK, P13,41-46, 50
12//      J. BARKER, August, 2001
13//      Switches from fast to slow convergence near target chi JGB 2/2003
14//
15// steps are:
16// load
17// mask
18// extrapolate (find the power law of the desmeared set = smeared -1)
19// smooth (optional)
20// desmear based on dQv and chi^2 target
21// save result
22//
23/////////
24// Waves produced at each step: (Q/I/S with the following extensions)
25//
26// Load:                Uses -> nothing (Kills old waves)
27//                              Produces -> "_exp" and "_exp_orig"
28//
29// Mask:                Uses -> "_exp"
30//                              Produces -> "_msk"
31//
32// Extrapolate: Uses -> nothing
33//                                      Produces -> nothing
34//
35// Smooth:      Uses -> "_smth" OR "_msk" OR "_exp" (in that order)
36//                              Produces -> "_smth"
37//
38// Desmear:     Uses ->  "_smth" OR "_msk" OR "_exp" (in that order)
39//                              Produces -> "_dsm"
40//
41////////
42
43///***** TO FIX *******
44// - ? power law + background extrapolation (only useful for X-ray data...)
45//
46// see commented code lines for Igor 4 or Igor 5
47// - there are only a few options and calls that are not Igor 4 compatible.
48// Igor 4 routines are currently used.
49
50
51////////////////////////////////////////////////////////////////////////
52
53// main entry routine
54Proc Desmear()
55        //check for the correct data folder, initialize if necessary
56        //
57        if(DataFolderExists("root:DSM") == 0)
58                Execute "Init_Desmearing()"
59        endif
60       
61        //always initialize these
62        String/G root:DSM:gStr1 = ""
63        String/G root:DSM:gStr2 = ""
64        String/G root:DSM:gIterStr = ""
65       
66        DoWindow/F Desmear_Graph
67        if(V_flag==0)
68                Execute "Desmear_Graph()"
69        endif
70End
71
72Proc Init_Desmearing()
73        //set up the folder(s) needed
74        NewDataFolder/O root:DSM
75        NewDataFolder/O root:myGlobals          //in case it wasn't created elsewhere
76       
77        String/G root:DSM:gCurFile=""
78       
79        Variable/G root:DSM:gMaxFastIter = 100                  //max number of iter in Fast convergence
80        Variable/G root:DSM:gMaxSlowIter = 10000
81       
82        Variable/G root:DSM:gNptsExtrap = 10            //points for high q extrapolation
83        Variable/G root:DSM:gChi2Target = 1             //chi^2 target
84        Variable/G root:DSM:gPowerM = -4
85        Variable/G root:DSM:gDqv = 0.117                        //2005 measured slit height - see John
86        Variable/G root:DSM:gNq = 1
87        Variable/G root:DSM:gS = 0              // global varaible for Midpnt()
88        Variable/G root:DSM:gSmoothFac=0.03
89       
90        Variable/G root:DSM:gChi2Final = 0              //chi^2 final
91        Variable/G root:DSM:gIterations = 0             //total number of iterations
92       
93        String/G root:DSM:gStr1 = ""                            //information strings
94        String/G root:DSM:gStr2 = ""
95        String/G root:DSM:gIterStr = ""
96        Variable/G root:DSM:gChi2Smooth = 0
97       
98        Variable/G root:DSM:gFreshMask=1
99End
100
101//////////// Lake desmearing routines
102
103//      Smears guess Y --> YS using weighting array
104Function DSM_Smear(Y_wave,Ys,NQ,FW)
105        Wave Y_wave,Ys
106        Variable nq
107        Wave FW
108       
109        Variable ii,jj,summ
110        for(ii=0;ii<nq;ii+=1)
111                summ=0
112                for(jj=0;jj<nq;jj+=1)
113                        summ = summ + Y_wave[jj]*FW[ii][jj]
114                endfor
115                Ys[ii] = summ
116        endfor
117       
118        Return (0)
119End
120
121//      CALCULATES CHI^2
122FUNCTION CHI2(ys_new,ys_exp,w,NQ)
123        Wave ys_new,ys_exp,w
124        Variable nq
125
126        Variable CHI2,summ,ii
127
128        SUMm = 0.0
129        for(ii=0;ii<nq;ii+=1)
130//              if(numtype(YS_EXP[ii]) == 0 && numtype(YS_NEW[ii]) == 0)
131                SUMm=SUMm+W[ii]*(YS_EXP[ii]-YS_NEW[ii])^2
132//              endif
133        endfor
134       
135        CHI2=SUMm/(NQ-1)
136       
137        RETURN CHI2
138END
139
140//      Routine calculates the weights needed to convert a table
141//      representation of a scattering function I(q) into the infinite
142//      slit smeared table represented as Is(q)
143//      Is(qi) = Sum(j=i-n) Wij*I(qj)
144//      Program assumes data is extrapolated from qn to dqv using powerlaw
145//      I(q) = Aq^m
146//      Required input:
147//      N       : Number of data points {in Global common block }
148//      q(N)    : Array of q values     {in Global common block }
149//      dqv     : Limit of slit length  {in Global common block }
150//      m       : powerlaw for extrapolation.
151//
152//      Routine output:
153//      W(N,N)  : Weights array
154//
155//      The weights obtained from this routine can be used to calculate
156//      Is(q) for any I(q) function.
157//
158//      Calculation based upon linear interpolation between points.
159//      See p.44, Computer Software Log book.
160//      11/00 JGB
161Function Weights_L(m,FW,Q_exp)          //changed input w to FW (FW is two dimensional)
162        Variable m
163        Wave FW,Q_exp
164       
165        NVAR dqv = root:DSM:gDqv
166        NVAR NN = root:DSM:gNq
167
168//      Calculate Remainder fractions and put into separate array.
169        Variable lower,ii,ss,jj
170        Duplicate/O Q_exp root:DSM:R_wave
171        wave r_wave = root:DSM:R_wave
172       
173//      Make/O/D/N=75 root:DSM:SS_save  //debug
174//      wave SS_save = root:DSM:SS_save
175       
176        Print "calculating remainders by integration..."
177        for(ii=0;ii<NN-1;ii+=1)
178                lower = sqrt(Q_exp[NN-1]^2-Q_exp[ii]^2)
179                ss = Qromo(lower,dqv,Q_exp[ii],m)               //new parameter list for Qromo call
180//              SS_save[ii] = ss
181                R_wave[ii] = (Q_exp[NN-1]^(-m)/dqv)*SS
182//              Printf "I = %d R_wave[ii] =%g\r",ii,R_wave[ii]
183        endfor
184       
185        lower = 0.0
186        ss = Qromo(lower,dqv,Q_exp[NN-1],m)             //new parameter list for Qromo call
187        R_wave[NN-1] = (Q_exp[NN-1]^(-m)/dqv)*SS
188//      Printf "I = %d R_wave[ii] =%g\r",NN-1,R_wave[NN-1]
189//      SS_save[ii] = ss
190
191//      Make/O/D/N=(75,75) root:DSM:IG_save             //debug
192//      wave IG_save=root:DSM:IG_save
193       
194//      Zero weight matrix... then fill it
195        FW = 0
196        Print "calculating full weights...."
197//      Calculate weights
198        for(ii=0;ii<NN;ii+=1)
199                for(jj=ii+1;jj<NN-1;jj+=1)
200                        FW[ii][jj] = DU(ii,jj)*(1.0+Q_exp[jj]/(Q_exp[jj+1]-Q_exp[jj]))
201                        FW[ii][jj] -= (1.0/(Q_exp[jj+1]-Q_exp[jj]))*IG(ii,jj)
202                        FW[ii][jj] -= DU(ii,jj-1)*Q_exp[jj-1]/(Q_exp[jj]-Q_exp[jj-1])
203                        FW[ii][jj] += (1.0/(Q_exp[jj]-Q_exp[jj-1]))*IG(ii,jj-1)
204                        FW[ii][jj] *= (1.0/dqv)                 
205//              Printf "FW[%d][%d] = %g\r",ii,jj,FW[ii][jj]
206                endfor
207        endfor
208//
209//      special case: I=J,I=N
210        for(ii=0;ii<NN-1;ii+=1)
211                FW[ii][ii] = DU(ii,ii)*(1.0+Q_exp[ii]/(Q_exp[ii+1]-Q_exp[ii]))
212                FW[ii][ii] -= (1.0/(Q_exp[ii+1]-Q_exp[ii]))*IG(ii,ii)
213                FW[ii][ii] *= (1.0/dqv)
214//              Printf "FW[%d][%d] = %g\r",ii,jj,FW[ii][jj]
215//     following line corrected for N -> NN-1 since Igor indexes from 0! (Q_exp[NN] DNE!)
216                FW[ii][NN-1] = -DU(ii,NN-2)*Q_exp[NN-2]/(Q_exp[NN-1]-Q_exp[NN-2])
217                FW[ii][NN-1] += (1.0/(Q_exp[NN-1]-Q_exp[NN-2]))*IG(ii,NN-2)
218                FW[ii][NN-1] *= (1.0/dqv)
219                FW[ii][NN-1] += R_wave[ii]
220//              Printf "FW[%d][%d] = %g\r",ii,NN-1,FW[ii][NN-1]
221        endfor
222//
223//      special case: I=J=N
224        FW[NN-1][NN-1] = R_wave[NN-1]
225//
226        SetDataFolder root:
227        Return (0)
228End
229
230////
231Function Integrand_dsm(u,qi,m)
232        Variable u,qi,m
233       
234        Variable integrand
235        Integrand = (u^2+qi^2)^(m/2)
236        return integrand
237end
238
239///
240Function DU(ii,jj)
241        Variable ii,jj
242       
243        Wave Q_exp = root:DSM:Q_exp
244        Variable DU
245
246//      Wave DU_save=root:DSM:DU_save
247        If (ii == jj)
248          DU = sqrt(Q_exp[jj+1]^2-Q_exp[ii]^2)
249        Else
250          DU = sqrt(Q_exp[jj+1]^2-Q_exp[ii]^2) - sqrt(Q_exp[jj]^2-Q_exp[ii]^2)
251        EndIf
252       
253//      DU_save[ii][jj] = DU
254        Return DU
255End
256
257Function IG(ii,jj)
258        Variable ii,jj
259       
260        Wave Q_exp=root:DSM:Q_exp
261        Variable IG,UL,UU
262
263//      WAVE IG_save = root:DSM:IG_save
264        If (ii == jj)
265          UL=0.0
266        Else
267          UL=sqrt(Q_exp[jj]^2-Q_exp[ii]^2)
268        EndIf
269        UU=sqrt(Q_exp[jj+1]^2-Q_exp[ii]^2)
270       
271        //in FORTRAN log = natural log....     
272        IG = UU*Q_exp[jj+1]+Q_exp[ii]^2*ln(UU+Q_exp[jj+1])
273        IG -= UL*Q_exp[jj]
274        IG -= Q_exp[ii]^2*ln(UL+Q_exp[jj])
275        IG *= 0.5
276//      IG_save[ii][jj] = IG
277       
278        Return IG
279End
280
281//
282Function FF(ii,jj)
283        Variable ii,jj
284
285        Variable FF
286        NVAR dqv = root:DSM:gDqv
287
288        FF = (1.0/dqv)*(0.5+HH(ii,jj))
289        Return FF
290End
291
292//
293Function GG(ii,jj)
294        Variable ii,jj
295
296        Variable GG
297        NVAR dqv = root:DSM:gDqv
298       
299        GG = (1.0/dqv)*(0.5-HH(ii,jj))
300        Return  GG
301End
302//
303Function HH(ii,jj)
304        Variable ii,jj
305
306        Wave Q_exp=root:DSM:Q_exp
307        Variable HH
308       
309        HH = 0.5*(Q_exp[jj+1]+Q_exp[jj])/(Q_exp[jj+1]-Q_exp[jj])
310        HH -= (1.0/(Q_exp[jj+1]-Q_exp[jj]))*(CC(ii,jj+1)-CC(ii,jj))
311        return HH
312End
313//
314//
315Function CC(ii,jj)
316        Variable ii,jj
317       
318        wave Q_exp = root:DSM:Q_exp
319        Variable CC
320       
321        If (ii == jj)
322          CC = 0.0
323        Else
324          CC = (Q_exp[jj]-Q_exp[ii])^0.5
325        EndIf
326        Return CC
327End
328
329// QROMO is a gerneric NR routine that takes function arguments
330// Call Qromo(Integrand,lower,dqv,Q_exp(N),m,SS,Midpnt)
331// -- here, it is always called with Integrand and Midpnt as the functions
332// -- so rewrite in a simpler way....
333//
334// SS is the returned value?
335// H_wave, S_wave (original names H,S)
336//
337// modified using c-version in NR (pg 143)
338Function QROMO(A,B,qi,m)
339        Variable A,B,qi,m
340       
341        Variable EPS,JMAX,JMAXP,KM,K
342        EPS=1.E-6
343        JMAX=14
344        KM=4
345        K=KM+1
346       
347        Make/O/D/N=(JMAX) root:DSM:S_wave
348        Make/O/D/N=(JMAX+1) root:DSM:H_wave
349        wave S_wave=root:DSM:S_wave
350        wave H_wave=root:DSM:H_wave
351        S_wave=0
352        H_wave=0
353       
354        H_Wave[0] = 1
355       
356        variable jj,SS,DSS
357       
358        for(jj=0;jj<jmax;jj+=1)
359                S_wave[jj] = Midpnt(A,B,jj,qi,m)                //remove FUNC, always call Integrand from within Midpnt
360                IF (jj>=KM)             //after 1st 5 points calculated
361                    POLINT(H_wave,S_wave,jj-KM,KM,0.0,SS,DSS)           //ss, dss returned
362                    IF (ABS(DSS) < EPS*ABS(SS))
363                        RETURN ss
364                    endif
365                ENDIF
366//        S_wave[jj+1]=S_wave[jj]
367          H_wave[jj+1]=H_wave[jj]/9.
368        endfor
369       
370        DoAlert 0,"Too many steps in QROMO"
371        return 1                        //error if you get here
372END
373
374//
375// see NR pg 109
376// - Given input arrays xa[1..n] and ya[1..n], and a given value x
377// return a value y and an error estimate dy. if P(x) is the polynomial of
378// degree N-1 such that P(xai) = yai, then the returned value y=P(x)
379//
380// arrays XA[] and YA[] are passed in with an offset of the index
381// of where to start the interpolation
382Function POLINT(XA,YA,offset,N,X,Y,DY)
383        Wave XA,YA
384        Variable offset,N,X,&Y,&DY
385       
386        Variable ii,mm,nmax,ns,dif,den,ho,hp,wi,dift
387       
388        NMAX=10
389       
390        Make/O/D/N=(NMAX) root:DSM:Ci,root:DSM:Di
391        wave Ci = root:DSM:Ci
392        wave Di = root:DSM:Di
393
394        NS=1
395        DIF=ABS(X-XA[0])
396        for(ii=0;ii<n;ii+=1)
397                DIFT=ABS(X-XA[ii+offset])
398                IF (DIFT < DIF)
399                        NS=ii
400                        DIF=DIFT
401                ENDIF
402                Ci[ii]=YA[ii+offset]
403                Di[ii]=YA[ii+offset]
404        endfor
405       
406        Y=YA[NS+offset]
407        NS=NS-1
408        for(mm=1;mm<n-1;mm+=1)
409                for(ii=0;ii<(n-mm);ii+=1)
410                        HO=XA[ii+offset]-X
411                        HP=XA[ii+offset+mm]-X
412                        Wi=Ci[ii+1]-Di[ii]
413                        DEN=HO-HP
414                        IF(DEN == 0.)
415                                print "den == 0 in POLINT - ERROR!!!"
416                        endif
417                        DEN=Wi/DEN
418                        Di[ii]=HP*DEN
419                        Ci[ii]=HO*DEN
420                endfor  //ii
421                IF (2*NS < (N-mm) )
422                        DY=Ci[NS+1]
423                ELSE
424                        DY=Di[NS]
425                        NS=NS-1
426                ENDIF
427                Y=Y+DY
428        endfor  //mm
429        RETURN (0)              //y and dy are returned as pass-by-reference
430END
431
432//
433// FUNC is always Integrand()
434// again, see the c-version, NR pg 142
435Function MIDPNT(A,B,N,qi,m)
436        Variable A,B,N,qi,m
437       
438        Variable it,tnm,del,ddel,x,summ,jj
439       
440        NVAR S_ret = root:DSM:gS
441       
442        IF (N == 0)
443                S_ret=(B-A)*Integrand_dsm(0.5*(A+B),qi,m)
444//              Print "N==0, S_ret = ",s_ret
445                return(S_ret)
446                        //IT=1
447        ELSE
448                it = 1
449                for(jj=1;jj<n-1;jj+=1)
450                        it *= 3
451                endfor
452                TNM=IT
453                DEL=(B-A)/(3.*TNM)
454                DDEL=DEL+DEL
455                X=A+0.5*DEL
456                SUMm=0.
457                for(jj=1;jj<=it;jj+=1)
458                        SUMm=SUMm+Integrand_dsm(X,qi,m)
459                        X=X+DDEL
460                        SUMm=SUMm+Integrand_dsm(X,qi,m)
461                        X=X+DEL
462                endfor
463                S_ret=(S_ret+(B-A)*SUMm/TNM)/3.
464                IT=3*IT
465        ENDIF
466        return S_ret
467END
468
469//////////// end of Lake desmearing routines
470
471// - Almost everything below here is utility routines for data handling, panel
472// controls, and all the fluff that makes this useable. Note that this is
473// a lot more code than the actual guts of the Lake method!
474// (the guts of the iteration are in the DemsearButtonProc)
475
476
477// add three "fake" resolution columns to the output file so that it
478// looks like a regular SANS data set, but make the sigQ column 100x smaller than Q
479// so that there is effectively no resolution smearing. Set Qbar = Q, and fs = 1
480//
481// SRK 06 FEB 06
482//
483Function WriteUSANSDesmeared(fullpath,lo,hi,dialog)
484        String fullpath
485        Variable lo,hi,dialog           //=1 will present dialog for name
486       
487        String termStr="\r\n"
488        String destStr = "root:DSM:"
489        String formatStr = "%15.6g %15.6g %15.6g %15.6g %15.6g %15.6g"+termStr
490       
491        Variable refNum,integer,realval
492       
493        //*****these waves MUST EXIST, or IGOR Pro will crash, with a type 2 error****
494        WAVE Q_dsm =$(destStr + "Q_dsm")
495        WAVE I_dsm=$(destStr + "I_dsm")
496        WAVE S_dsm=$(destStr + "S_dsm")
497       
498        //check each wave
499        If(!(WaveExists(Q_dsm)))
500                Abort "Q_dsm DNExist in WriteUSANSDesmeared()"
501        Endif
502        If(!(WaveExists(I_dsm)))
503                Abort "I_dsm DNExist in WriteUSANSDesmeared()"
504        Endif
505        If(!(WaveExists(S_dsm)))
506                Abort "S_dsm DNExist in WriteUSANSDesmeared()"
507        Endif
508       
509        // 06 FEB 06 SRK
510        // make dummy waves to hold the "fake" resolution, and write it as the last 3 columns
511        //
512        Duplicate/O Q_dsm,res1,res2,res3
513        res3 = 1                // "fake" beamstop shadowing
514        res1 /= 100             //make the sigmaQ so small that there is no smearing
515       
516        if(dialog)
517                Open/D refnum as fullpath+".txt"                //won't actually open the file
518                If(cmpstr(S_filename,"")==0)
519                        //user cancel, don't write out a file
520                        Close/A
521                        Abort "no data file was written"
522                Endif
523                fullpath = S_filename
524        Endif
525       
526        //write out partial set?
527        Duplicate/O Q_dsm,tq,ti,te
528        ti=I_dsm
529        te=S_dsm
530        if( (lo!=hi) && (lo<hi))
531                redimension/N=(hi-lo+1) tq,ti,te,dumWave                //lo to hi, inclusive
532                tq=Q_dsm[p+lo]
533                ti=I_dsm[p+lo]
534                te=S_dsm[p+lo]
535        endif
536       
537        //tailor the output given the type of data written out...
538        String samStr="",dateStr="",str1,str2
539       
540        NVAR m = root:DSM:gPowerM                                       // power law exponent
541        NVAR chiFinal = root:DSM:gChi2Final             //chi^2 final
542        NVAR iter = root:DSM:gIterations                //total number of iterations
543       
544        //get the number of spline passes from the wave note
545        String noteStr
546        Variable boxPass,SplinePass
547        noteStr=note(I_dsm)
548        BoxPass = NumberByKey("BOX", noteStr, "=", ";")
549        splinePass = NumberByKey("SPLINE", noteStr, "=", ";")
550       
551        samStr = fullpath
552        dateStr="CREATED: "+date()+" at  "+time()
553        sprintf str1,"Chi^2 = %g   PowerLaw m = %4.2f   Iterations = %d",chiFinal,m,iter
554        sprintf str2,"%d box smooth passes and %d smoothing spline passes",boxPass,splinePass
555       
556        //actually open the file
557        Open refNum as fullpath
558       
559        fprintf refnum,"%s"+termStr,samStr
560        fprintf refnum,"%s"+termStr,str1
561        fprintf refnum,"%s"+termStr,str2
562        fprintf refnum,"%s"+termStr,dateStr
563       
564        wfprintf refnum, formatStr, tq,ti,te,res1,res2,res3
565       
566        Close refnum
567       
568        Killwaves/Z ti,tq,te,res1,res2,res3
569       
570        Return(0)
571End
572
573/// procedures to do the extrapolation to high Q
574// very similar to the procedures in the Invariant package
575//
576//create the wave to extrapolate
577// w is the input q-values
578Function DSM_SetExtrWaves(w)
579        Wave w
580
581        Variable num_extr=25
582       
583        SetDataFolder root:DSM
584       
585        Make/O/D/N=(num_extr) extr_hqq,extr_hqi
586        extr_hqi=1              //default values
587
588        //set the q-range
589        Variable qmax,num
590//      qmax=0.03
591
592        num=numpnts(w)
593        qmax=6*w[num-1]
594       
595        extr_hqq = w[num-1] + x * (qmax-w[num-1])/num_extr
596       
597        SetDataFolder root:
598        return(0)
599End
600
601//creates I_ext,Q_ext,S_ext with extended q-values
602//
603// if num_extr == 0 , the input waves are returned as _ext
604//
605// !! uses simple linear extrapolation at low Q - just need something
606// reasonable in the waves to keep from getting smoothing artifacts
607// - uses power law at high q
608Function ExtendToSmooth(qw,iw,sw,nbeg,nend,num_extr)
609        Wave qw,iw,sw
610        Variable nbeg,nend,num_extr
611       
612        Variable/G V_FitMaxIters=300
613        Variable/G V_fitOptions=4               //suppress the iteration window
614        Variable num_new,qmax,num,qmin
615        num=numpnts(qw)
616       
617        num_new = num + 2*num_extr
618       
619//      Print "num,num_new",num,num_new
620       
621        SetDataFolder root:DSM
622        Make/O/D/N=(num_new) Q_ext,I_ext,S_ext
623       
624        if(num_extr == 0)               //the extended waves are the input, get out
625                Q_ext = qw
626                I_ext = iw
627                S_ext = sw
628                setDatafolder root:
629                return (0)
630        endif
631       
632        //make the extensions
633        Make/O/D/N=(num_extr) hqq,hqi,lqq,lqi
634        hqi=1           //default values
635        lqi=0
636        //set the q-range based on the high/low values of the data set
637//      qmin=1e-6
638        qmin= qw[0]*0.5
639        qmax = qw[num-1]*2
640       
641        hqq = qw[num-1] + x * (qmax-qw[num-1])/num_extr
642        lqq = qmin + x * (qw[0]-qmin)/num_extr
643       
644        //do the fits
645        Duplicate/O iw dummy                    //use this as the destination to suppress fit_ from being appended
646       
647        //Use simple linear fits        line: y = K0+K1*x
648        CurveFit/Q line iw[0,(nbeg-1)] /X=qw /W=sw /D=dummy
649        Wave W_coef=W_coef
650        lqi=W_coef[0]+W_coef[1]*lqq
651//     
652// Guinier or Power-law fits
653//             
654//      Make/O/D G_coef={100,-100}              //input
655//      FuncFit DSM_Guinier_Fit G_coef iw[0,(nbeg-1)] /X=qw /W=sw /D
656//      lqi= DSM_Guinier_Fit(G_coef,lqq)
657       
658//      Printf "I(q=0) = %g (1/cm)\r",G_coef[0]
659//      Printf "Rg = %g ()\r",sqrt(-3*G_coef[1])
660               
661        Make/O/D P_coef={0,1,-4}                        //input  --- (set background to zero and hold fixed)
662        CurveFit/Q/N/H="100" Power kwCWave=P_coef  iw[(num-1-nend),(num-1)] /X=qw /W=sw /D=dummy
663        hqi=P_coef[0]+P_coef[1]*hqq^P_coef[2]
664       
665//      Printf "Power law exponent = %g\r",P_coef[2]
666//      Printf "Pre-exponential = %g\r",P_coef[1]
667       
668        // concatenate the extensions
669        Q_ext[0,(num_extr-1)] = lqq[p]
670        Q_ext[num_extr,(num_extr+num-1)] = qw[p-num_extr]
671        Q_ext[(num_extr+num),(2*num_extr+num-1)] = hqq[p-num_extr-num]
672        I_ext[0,(num_extr-1)] = lqi[p]
673        I_ext[num_extr,(num_extr+num-1)] = iw[p-num_extr]
674        I_ext[(num_extr+num),(2*num_extr+num-1)] = hqi[p-num_extr-num]
675        S_ext[0,(num_extr-1)] = sw[0]
676        S_ext[num_extr,(num_extr+num-1)] = sw[p-num_extr]
677        S_ext[(num_extr+num),(2*num_extr+num-1)] = sw[num-1]
678       
679        killwaves/z dummy
680        SetDataFolder root:
681        return(0)
682End
683
684// pass in the (smeared) experimental data and find the power-law extrapolation
685// 10 points is usually good enough, unless the data is crummy
686// returns the prediction for the exponent
687//
688Function DSM_DoExtrapolate(qw,iw,sw,nend)
689        Wave qw,iw,sw
690        Variable nend
691       
692        Setdatafolder root:DSM
693       
694//      Wave extr_lqi=extr_lqi
695//      Wave extr_lqq=extr_lqq
696        Wave extr_hqi=extr_hqi
697        Wave extr_hqq=extr_hqq
698        Variable/G V_FitMaxIters=300
699        Variable/G V_fitOptions=4               //suppress the iteration window
700        Variable num=numpnts(iw),retVal
701       
702        Make/O/D P_coef={0,1,-4}                        //input
703//      Make/O/T Constr={"K2<0","K2 > -20"}
704        //(set background to zero and hold fixed)
705        CurveFit/H="100" Power kwCWave=P_coef  iw[(num-1-nend),(num-1)] /X=qw /W=sw /D
706        extr_hqi=P_coef[0]+P_coef[1]*extr_hqq^P_coef[2]
707       
708        // for the case of data with a background
709//      Make/O/D P_coef={0,1,-4}                        //input
710//      //(set background to iw[num-1], let it be a free parameter
711//      P_coef[0] = iw[num-1]
712//      CurveFit Power kwCWave=P_coef  iw[(num-1-nend),(num-1)] /X=qw /W=sw /D
713//      extr_hqi=P_coef[0]+P_coef[1]*extr_hqq^P_coef[2]
714       
715       
716        Printf "Smeared Power law exponent = %g\r",P_coef[2]
717        Printf "**** For Desmearing, use a Power law exponent of %5.1f\r",P_coef[2]-1
718       
719        retVal = P_coef[2]-1
720        SetDataFolder root:
721        return(retVal)
722End
723
724Function DSM_Guinier_Fit(w,x) //: FitFunc
725        Wave w
726        Variable x
727       
728        //fit data to I(q) = A*exp(B*q^2)
729        // (B will be negative)
730        //two parameters
731        Variable a,b,ans
732        a=w[0]
733        b=w[1]
734        ans = a*exp(b*x*x)
735        return(ans)
736End
737
738Proc Desmear_Graph()
739        PauseUpdate; Silent 1           // building window...
740        Display /W=(5,44,408,558) /K=1
741        ModifyGraph cbRGB=(51664,44236,58982)
742        DoWindow/C Desmear_Graph
743        ControlBar 160
744        // break into tabs
745        TabControl DSM_Tab,pos={5,3},size={392,128},proc=DSM_TabProc
746        TabControl DSM_Tab,labelBack=(49151,49152,65535),tabLabel(0)="Load"
747        TabControl DSM_Tab,tabLabel(1)="Mask",tabLabel(2)="Extrapolate"
748        TabControl DSM_Tab,tabLabel(3)="Smooth",tabLabel(4)="Desmear",value= 0
749       
750        //always visible - revert and save
751        //maybe the wrong place here?
752        Button DSMControlA,pos={225,135},size={80,20},proc=DSM_RevertButtonProc,title="Revert"
753        Button DSMControlA,help={"Revert the smeared data to its original state and start over"}
754        Button DSMControlB,pos={325,135},size={50,20},proc=DSM_SaveButtonProc,title="Save"
755        Button DSMControlB,help={"Save the desmeared data set"}
756        Button DSMControlC,pos={25,135},size={50,20},proc=DSM_HelpButtonProc,title="Help"
757        Button DSMControlC,help={"Show the help file for desmearing"}
758       
759        // add the controls to each tab ---- all names start with "DSMControl_"
760
761        //tab(0) Load - initially visible
762        Button DSMControl_0a,pos={23,39},size={80,20},proc=DSM_LoadButtonProc,title="Load Data"
763        Button DSMControl_0a,help={"Load slit-smeared USANS data = \".cor\" files"}
764        CheckBox DSMControl_0b,pos={26,74},size={80,14},proc=DSM_LoadCheckProc,title="Log Axes?"
765        CheckBox DSMControl_0b,help={"Toggle Log/Lin Q display"},value= 1
766        TitleBox DSMControl_0c,pos={120,37},size={104,19},font="Courier",fSize=10
767        TitleBox DSMControl_0c,variable= root:DSM:gStr1
768        //second message string not used currently
769//      TitleBox DSMControl_0d,pos={120,57},size={104,19},font="Courier",fSize=10
770//      TitleBox DSMControl_0d,variable= root:DSM:gStr2
771       
772        //tab(1) Mask
773        Button DSMControl_1a,pos={20,35},size={90,20},proc=DSM_MyMaskProc,title="Mask Point"            //bMask
774        Button DSMControl_1a,help={"Toggles the masking of the selected data point"}
775        Button DSMControl_1a,disable=1
776        Button DSMControl_1b,pos={20,65},size={140,20},proc=DSM_MaskGTCursor,title="Mask Q >= Cursor"           //bMask
777        Button DSMControl_1b,help={"Toggles the masking of all q-values GREATER than the current cursor location"}
778        Button DSMControl_1b,disable=1
779        Button DSMControl_1c,pos={20,95},size={140,20},proc=DSM_MaskLTCursor,title="Mask Q <= Cursor"           //bMask
780        Button DSMControl_1c,help={"Toggles the masking of all q-values LESS than the current cursor location"}
781        Button DSMControl_1c,disable=1
782        Button DSMControl_1d,pos={180,35},size={90,20},proc=DSM_ClearMaskProc,title="Clear Mask"                //bMask
783        Button DSMControl_1d,help={"Clears all mask points"}
784        Button DSMControl_1d,disable=1
785//      Button DSMControl_1b,pos={144,66},size={110,20},proc=DSM_MaskDoneButton,title="Done Masking"
786//      Button DSMControl_1b,disable=1
787       
788        //tab(2) Extrapolate
789        Button DSMControl_2a,pos={31,42},size={90,20},proc=DSM_ExtrapolateButtonProc,title="Extrapolate"
790        Button DSMControl_2a,help={"Extrapolate the high-q region with a power-law"}
791        Button DSMControl_2a,disable=1
792        SetVariable DSMControl_2b,pos={31,70},size={100,15},title="# of points"
793        SetVariable DSMControl_2b,help={"Set the number of points for the power-law extrapolation"}
794        SetVariable DSMControl_2b,limits={5,100,1},value=  root:DSM:gNptsExtrap
795        SetVariable DSMControl_2b,disable=1
796        CheckBox DSMControl_2c,pos={157,45},size={105,14},proc=DSM_ExtrapolationCheckProc,title="Show Extrapolation"
797        CheckBox DSMControl_2c,help={"Show or hide the high q extrapolation"},value= 1
798        CheckBox DSMControl_2c,disable=1
799        SetVariable DSMControl_2d,pos={31,96},size={150,15},title="Power Law Exponent"
800        SetVariable DSMControl_2d,help={"Power Law exponent from the fit = the DESMEARED slope - override as needed"}
801        SetVariable DSMControl_2d format="%5.2f"
802        SetVariable DSMControl_2d,limits={-inf,inf,0},value= root:DSM:gPowerM
803        SetVariable DSMControl_2d,disable=1
804       
805        //tab(3) Smooth
806        Button DSMControl_3a,pos={34,97},size={70,20},proc=DSM_SmoothButtonProc,title="Smooth"
807        Button DSMControl_3a,disable=1
808                //BoxCheck
809        CheckBox DSMControl_3b,pos={34,39},size={35,14},title="Box",value= 1
810        CheckBox DSMControl_3b,help={"Use a single pass of 3-point box smoothing"}
811        CheckBox DSMControl_3b,disable=1
812                //SSCheck
813        CheckBox DSMControl_3c,pos={34,60},size={45,14},title="Spline",value= 0
814        CheckBox DSMControl_3c,help={"Use a single pass of a smoothing spline"}
815        CheckBox DSMControl_3c,disable=1
816                //extendCheck
817        CheckBox DSMControl_3d,pos={268,60},size={71,14},title="Extend Data"
818        CheckBox DSMControl_3d,help={"extends the data at both low q and high q to avoid end effects in smoothing"}
819        CheckBox DSMControl_3d,value= 0
820        CheckBox DSMControl_3d,disable=1
821        Button DSMControl_3e,pos={125,97},size={90,20},proc=DSM_SmoothUndoButtonProc,title="Start Over"
822        Button DSMControl_3e,help={"Start the smoothing over again without needing to re-mask the data set"}
823        Button DSMControl_3e,disable=1
824        SetVariable DSMControl_3f,pos={94,60},size={150,15},title="Smoothing factor"
825        SetVariable DSMControl_3f,help={"Smoothing factor for the smoothing spline"}
826        SetVariable DSMControl_3f format="%5.4f"
827        SetVariable DSMControl_3f,limits={0.01,2,0.01},value= root:DSM:gSmoothFac
828        SetVariable DSMControl_3f,disable=1
829        CheckBox DSMControl_3g,pos={268,39},size={90,14},title="Log-scale smoothing?"
830        CheckBox DSMControl_3g,help={"Use log-scaled intensity during smoothing (reverts to linear if negative intensity points found)"}
831        CheckBox DSMControl_3g,value=0
832        CheckBox DSMControl_3g,disable=1
833        ValDisplay DSMControl_3h pos={235,97},title="Chi^2",size={80,20},value=root:DSM:gChi2Smooth
834        ValDisplay DSMControl_3h,help={"This is the Chi^2 value for the smoothed data vs experimental data"}
835        ValDisplay DSMControl_3h,disable=1
836       
837        //tab(4) Desmear
838        Button DSMControl_4a,pos={35,93},size={80,20},proc=DSM_DesmearButtonProc,title="Desmear"
839        Button DSMControl_4a,help={"Do the desmearing - the result is in I_dsm"}
840        Button DSMControl_4a,disable=1
841        SetVariable DSMControl_4b,pos={35,63},size={120,15},title="Chi^2 target"
842        SetVariable DSMControl_4b,help={"Set the targetchi^2 for convergence (recommend chi^2=1)"}
843        SetVariable DSMControl_4b,limits={0,inf,0.1},value= root:DSM:gChi2Target
844        SetVariable DSMControl_4b,disable=1
845        SetVariable DSMControl_4c,pos={35,35},size={80,15},title="dQv"
846        SetVariable DSMControl_4c,help={"Slit height as read in from the data file. 0.117 is the NIST value, override if necessary"}
847        SetVariable DSMControl_4c,limits={-inf,inf,0},value= root:DSM:gDqv
848        SetVariable DSMControl_4c,disable=1
849        TitleBox DSMControl_4d,pos={160,37},size={104,19},font="Courier",fSize=10
850        TitleBox DSMControl_4d,variable= root:DSM:gIterStr
851        TitleBox DSMControl_4d,disable=1
852       
853       
854        SetDataFolder root:
855EndMacro
856
857// function to control the drawing of buttons in the TabControl on the main panel
858// Naming scheme for the buttons MUST be strictly adhered to... else buttons will
859// appear in odd places...
860// all buttons are named DSMControl_NA where N is the tab number and A is the letter denoting
861// the button's position on that particular tab.
862// in this way, buttons will always be drawn correctly :-)
863//
864Function DSM_TabProc(ctrlName,tab) //: TabControl
865        String ctrlName
866        Variable tab
867
868        String ctrlList = ControlNameList("",";"),item="",nameStr=""
869        Variable num = ItemsinList(ctrlList,";"),ii,onTab
870        for(ii=0;ii<num;ii+=1)
871                //items all start w/"DSMControl_"               //11 characters
872                item=StringFromList(ii, ctrlList ,";")
873                nameStr=item[0,10]
874                if(cmpstr(nameStr,"DSMControl_")==0)
875                        onTab = str2num(item[11])                       //12th is a number
876                        ControlInfo $item
877                        switch(abs(V_flag))     
878                                case 1:
879                                        Button $item,disable=(tab!=onTab)
880                                        break
881                                case 2:
882                                        CheckBox $item,disable=(tab!=onTab)
883                                        break
884                                case 5:
885                                        SetVariable $item,disable=(tab!=onTab)
886                                        break
887                                case 10:       
888                                        TitleBox $item,disable=(tab!=onTab)
889                                        break
890                                case 4:
891                                        ValDisplay $item,disable=(tab!=onTab)
892                                        break
893                                // add more items to the switch if different control types are used
894                        endswitch
895                       
896                endif
897        endfor
898       
899        if(tab==1)
900                DSM_MyMaskProc("")              //start maksing if you click on the tab
901        else
902                DSM_MaskDoneButton("")          //masking is done if you click off the tab
903        endif
904       
905        return 0
906End
907
908Proc AppendSmeared()
909        SetDataFolder root:DSM
910//      if( strsearch(TraceNameList("Desmear_Graph", "", 1),"I_exp_orig",0,2) == -1)            //Igor 5
911        if( strsearch(TraceNameList("Desmear_Graph", "", 1),"I_exp_orig",0) == -1)
912                AppendToGraph/W=Desmear_Graph  I_exp_orig vs Q_exp_orig
913                ModifyGraph mode=4,marker=19            //3 is markers, 4 is markers and lines
914                ModifyGraph rgb(I_exp_orig)=(0,0,0)
915                ModifyGraph msize=2,grid=1,log=1,mirror=2,standoff=0,tickunit=1
916                ErrorBars I_exp_orig Y,wave=(S_exp_orig,S_exp_orig)
917                Legend/N=text0/J "\\F'Courier'\\s(I_exp_orig) I_exp_orig"
918                Label left "Intensity"
919                Label bottom "Q (1/A)"
920        endif
921        //always update the textbox - kill the old one first
922        TextBox/K/N=text1
923//      TextBox/C/N=text1/F=0/A=MT/E=2/X=5.50/Y=0.00 root:DSM:gCurFile                  //Igor 5
924        TextBox/C/N=text1/F=0/A=MT/E=1/X=5.50/Y=0.00 root:DSM:gCurFile
925End
926
927Proc AppendMask()
928//      if( strsearch(TraceNameList("Desmear_Graph", "", 1),"MaskData",0,2) == -1)                      //Igor 5
929        if( strsearch(TraceNameList("Desmear_Graph", "", 1),"MaskData",0) == -1)
930                SetDataFolder root:DSM:
931                AppendToGraph/W=Desmear_Graph MaskData vs Q_exp_orig
932                ModifyGraph mode(MaskData)=3,marker(MaskData)=8,msize(MaskData)=2.5,opaque(MaskData)=1
933                ModifyGraph rgb(MaskData)=(65535,16385,16385)
934                setdatafolder root:
935        endif
936end
937
938Proc AppendSmoothed()
939//      if( strsearch(TraceNameList("Desmear_Graph", "", 1),"I_smth",0,2) == -1)                        //Igor 5
940        if( strsearch(TraceNameList("Desmear_Graph", "", 1),"I_smth",0) == -1)
941                SetDataFolder root:DSM:
942                AppendToGraph/W=Desmear_Graph I_smth vs Q_smth
943                ModifyGraph/W=Desmear_Graph rgb(I_smth)=(3,52428,1),lsize(I_smth)=2
944                setdatafolder root:
945        endif
946end
947
948Function RemoveSmoothed()
949        SetDataFolder root:DSM:
950        RemoveFromGraph/W=Desmear_Graph/Z I_smth
951        setdatafolder root:
952end
953
954Function RemoveMask()
955        SetDataFolder root:DSM:
956        RemoveFromGraph/W=Desmear_Graph/Z MaskData
957        setdatafolder root:
958end
959
960Proc AppendDesmeared()
961//      if( strsearch(TraceNameList("Desmear_Graph", "", 1),"I_dsm",0,2) == -1)         //Igor 5
962        if( strsearch(TraceNameList("Desmear_Graph", "", 1),"I_dsm",0) == -1)
963                SetDataFolder root:DSM:
964                AppendToGraph/W=Desmear_Graph I_dsm vs Q_dsm
965                ModifyGraph mode(I_dsm)=3,marker(I_dsm)=19
966                ModifyGraph rgb(I_dsm)=(1,16019,65535),msize(I_dsm)=2
967                ErrorBars I_dsm Y,wave=(S_dsm,S_dsm)
968                setdatafolder root:
969        endif
970end
971
972Function RemoveDesmeared()
973        SetDataFolder root:DSM:
974        RemoveFromGraph/W=Desmear_Graph/Z I_dsm
975        setdatafolder root:
976end
977
978Function AppendExtrapolation()
979//      if( strsearch(TraceNameList("Desmear_Graph", "", 1),"extr_hqi",0,2) == -1)              //Igor 5
980        if( strsearch(TraceNameList("Desmear_Graph", "", 1),"extr_hqi",0) == -1)
981                SetDataFolder root:DSM:
982                AppendToGraph/W=Desmear_Graph extr_hqi vs extr_hqq
983                ModifyGraph/W=Desmear_Graph lSize(extr_hqi)=2
984                setdatafolder root:
985        endif
986end
987
988Function RemoveExtrapolation()
989        SetDataFolder root:DSM:
990        RemoveFromGraph/W=Desmear_Graph/Z extr_hqi
991        setdatafolder root:
992end
993
994// step (1) - read in the data, and plot it
995// clear out all of the "old" waves, remove them from the graph first
996// reads in a fresh copy of the data
997//
998// produces Q_exp, I_exp, S_exp waves (and originals "_orig")
999// add a dummy wave note that can be changed on later steps
1000//
1001Function DSM_LoadButtonProc(ctrlName) : ButtonControl
1002        String ctrlName
1003
1004        String qStr,iStr,sStr,sqStr
1005        Variable nq,dqv,numBad,val
1006       
1007        // remove any of the old traces on the graph and delete the waves
1008        CleanUpJunk()
1009       
1010        SetDataFolder root:
1011        Execute "U_LoadOneDDataWithName(\"\")"
1012        //define the waves that the smoothing will be looking for...
1013        SVAR fname = root:myGlobals:gLastFileName                       //this changes as any data is loaded
1014        SVAR curFile = root:DSM:gCurFile                                                //keep this for title, save
1015        curFile = fname
1016       
1017        qStr = CleanupName((fName + "_q"),0)            //the q-wave
1018        iStr = CleanupName((fName + "_i"),0)            //the i-wave
1019        sStr = CleanupName((fName + "_s"),0)            //the s-wave
1020        sqStr = CleanupName((fName + "sq"),0)           //the sq-wave, which should have -dQv as the elements
1021
1022        Duplicate/O $qStr root:DSM:Q_exp                        //root:DSM:Q_exp_orig
1023        Duplicate/O $iStr root:DSM:I_exp                        //root:DSM:I_exp_orig
1024        Duplicate/O $sStr root:DSM:S_exp                //root:DSM:S_exp_orig
1025        wave Q_exp = root:DSM:Q_exp
1026        Wave I_exp = root:DSM:I_exp
1027        Wave S_exp = root:DSM:S_exp
1028        Wave/Z sigQ = $sqStr
1029       
1030        // remove any negative q-values (and q=0 values!)(and report this)
1031        // ? and trim the low q to be >= 3.0e-5 (1/A), below this USANS is not reliable.
1032        NumBad = RemoveBadQPoints(Q_exp,I_exp,S_exp,0)
1033        SVAR str1 = root:DSM:gStr1
1034        sprintf str1,"%d negative q-values were removed",numBad
1035       
1036// don't trim off any positive q-values
1037//      val = 3e-5              //lowest "good" q-value from USANS
1038//      NumBad = RemoveBadQPoints(Q_exp,I_exp,S_exp,val-1e-8)
1039//      SVAR str2 = root:DSM:gStr2
1040//      sprintf str2,"%d q-values below q = %g were removed",numBad,val
1041       
1042        Duplicate/O root:DSM:Q_exp root:DSM:Q_exp_orig
1043        Duplicate/O root:DSM:I_exp root:DSM:I_exp_orig
1044        Duplicate/O root:DSM:S_exp root:DSM:S_exp_orig
1045        wave I_exp_orig = root:DSM:I_exp_orig
1046       
1047        nq = numpnts(root:DSM:Q_exp)
1048       
1049        if(WaveExists(sigQ))                    //try to read dQv
1050                dqv = -sigQ[0]
1051//              DoAlert 0,"Found dQv value of " + num2str(dqv)
1052        else
1053                dqv = 0.117
1054        //      dqv = 0.037             //old value
1055                DoAlert 0,"Could not find dQv in the data file - using " + num2str(dqv)
1056        endif
1057        NVAR gDqv = root:DSM:gDqv                               //needs to be global for Weights_L()
1058        NVAR gNq = root:DSM:gNq
1059        //reset the globals
1060        gDqv = dqv
1061        gNq = nq
1062       
1063        // append the (blank) wave note to the intensity wave
1064        Note I_exp,"BOX=0;SPLINE=0;"
1065        Note I_exp_orig,"BOX=0;SPLINE=0;"
1066       
1067        //draw the graph
1068        Execute "AppendSmeared()"
1069       
1070        SetDataFolder root:
1071End
1072
1073// remove any q-values <= val
1074Function RemoveBadQPoints(qw,iw,sw,val)
1075        Wave qw,iw,sw
1076        Variable val
1077       
1078        Variable ii,num,numBad,qval
1079        num = numpnts(qw)
1080       
1081        ii=0
1082        numBad=0
1083        do
1084                qval = qw[ii]
1085                if(qval <= val)
1086                        numBad += 1
1087                else            //keep the points
1088                        qw[ii-numBad] = qval
1089                        iw[ii-numBad] = iw[ii]
1090                        sw[ii-numBad] = sw[ii]
1091                endif
1092                ii += 1
1093        while(ii<num)
1094        //trim the end of the waves
1095        DeletePoints num-numBad, numBad, qw,iw,sw
1096        return(numBad)
1097end
1098
1099// if mw = Nan, keep the point, if a numerical value, delete it
1100Function RemoveMaskedPoints(mw,qw,iw,sw)
1101        Wave mw,qw,iw,sw
1102       
1103        Variable ii,num,numBad,mask
1104        num = numpnts(qw)
1105       
1106        ii=0
1107        numBad=0
1108        do
1109                mask = mw[ii]
1110                if(numtype(mask) != 2)          //if not NaN
1111                        numBad += 1
1112                else            //keep the points that are NaN
1113                        qw[ii-numBad] = qw[ii]
1114                        iw[ii-numBad] = iw[ii]
1115                        sw[ii-numBad] = sw[ii]
1116                endif
1117                ii += 1
1118        while(ii<num)
1119        //trim the end of the waves
1120        DeletePoints num-numBad, numBad, qw,iw,sw
1121        return(numBad)
1122end
1123
1124// produces the _msk waves that have the new number of data points
1125//
1126Function DSM_MaskDoneButton(ctrlName) : ButtonControl
1127        String ctrlName
1128
1129//      Variable aExists= strlen(CsrInfo(A)) > 0                        //Igor 5
1130        Variable aExists= strlen(CsrWave(A)) > 0                        //Igor 4
1131        if(!aExists)
1132                return(1)               //possibly reverted data, no cursor, no Mask wave
1133        endif
1134       
1135        Duplicate/O root:DSM:Q_exp_orig,root:DSM:Q_msk
1136        Duplicate/O root:DSM:I_exp_orig,root:DSM:I_msk
1137        Duplicate/O root:DSM:S_exp_orig,root:DSM:S_msk
1138        Wave Q_msk=root:DSM:Q_msk
1139        Wave I_msk=root:DSM:I_msk
1140        Wave S_msk=root:DSM:S_msk
1141       
1142        //finish up - trim the data sets and reassign the working set
1143        Wave MaskData=root:DSM:MaskData
1144       
1145        RemoveMaskedPoints(MaskData,Q_msk,I_msk,S_msk)
1146
1147        //reset the number of points
1148        NVAR gNq = root:DSM:gNq
1149        gNq = numpnts(Q_msk)
1150       
1151        Cursor/K A
1152        HideInfo
1153       
1154        return(0)
1155End
1156
1157
1158// not quite the same as revert
1159Function DSM_ClearMaskProc(ctrlName) : ButtonControl
1160        String ctrlName
1161       
1162        Wave MaskData=root:DSM:MaskData
1163        MaskData = NaN
1164       
1165        return(0)
1166end
1167
1168// when the mask button is pressed, A must be on the graph
1169// Displays MaskData wave on the graph
1170//
1171Function DSM_MyMaskProc(ctrlName) : ButtonControl
1172        String ctrlName
1173       
1174        Wave data=root:DSM:I_exp_orig
1175       
1176//      Variable aExists= strlen(CsrInfo(A)) > 0                        //Igor 5
1177        Variable aExists= strlen(CsrWave(A)) > 0                        //Igor 4
1178       
1179// need to get rid of old smoothed data if data is re-masked
1180        Execute "RemoveSmoothed()"
1181        SetDataFolder root:DSM
1182        Killwaves/Z I_smth,Q_smth,S_smth
1183        SetDataFolder root:
1184               
1185        if(aExists)             //mask the selected point
1186                // toggle NaN (keep) or Data value (= masked)
1187                Wave MaskData=root:DSM:MaskData
1188                MaskData[pcsr(A)] = (numType(MaskData[pcsr(A)])==0) ? NaN : data[pcsr(A)]               //if NaN, doesn't plot
1189        else
1190                Cursor /A=1/H=1/L=1/P/W=Desmear_Graph A I_exp_orig leftx(I_exp_orig)
1191                ShowInfo
1192                //if the mask wave does not exist, make one
1193                if(exists("root:DSM:MaskData") != 1)
1194                        Duplicate/O root:DSM:Q_exp_orig root:DSM:MaskData
1195                        Wave MaskData=root:DSM:MaskData
1196                        MaskData = NaN          //use all data
1197                endif
1198                Execute "AppendMask()" 
1199        endif
1200       
1201        return(0)
1202End
1203
1204// when the mask button is pressed, A must be on the graph
1205// Displays MaskData wave on the graph
1206//
1207Function DSM_MaskLTCursor(ctrlName) : ButtonControl
1208        String ctrlName
1209       
1210        Wave data=root:DSM:I_exp_orig
1211       
1212//      Variable aExists= strlen(CsrInfo(A)) > 0                        //Igor 5
1213        Variable aExists= strlen(CsrWave(A)) > 0                        //Igor 4
1214       
1215        if(!aExists)
1216                return(1)
1217        endif
1218// need to get rid of old smoothed data if data is re-masked
1219        Execute "RemoveSmoothed()"
1220        SetDataFolder root:DSM
1221        Killwaves/Z I_smth,Q_smth,S_smth
1222        SetDataFolder root:
1223
1224        Wave MaskData=root:DSM:MaskData
1225        Variable pt,ii
1226        pt = pcsr(A)
1227        for(ii=pt;ii>=0;ii-=1)
1228                // toggle NaN (keep) or Data value (= masked)
1229                MaskData[ii] = (numType(MaskData[ii])==0) ? NaN : data[ii]              //if NaN, doesn't plot
1230        endfor
1231        return(0)
1232End
1233
1234// when the mask button is pressed, A must be on the graph
1235// Displays MaskData wave on the graph
1236//
1237Function DSM_MaskGTCursor(ctrlName) : ButtonControl
1238        String ctrlName
1239       
1240        Wave data=root:DSM:I_exp_orig
1241       
1242//      Variable aExists= strlen(CsrInfo(A)) > 0                        //Igor 5
1243        Variable aExists= strlen(CsrWave(A)) > 0                        //Igor 4
1244       
1245        if(!aExists)
1246                return(1)
1247        endif
1248// need to get rid of old smoothed data if data is re-masked
1249        Execute "RemoveSmoothed()"
1250        SetDataFolder root:DSM
1251        Killwaves/Z I_smth,Q_smth,S_smth
1252        SetDataFolder root:
1253
1254        Wave MaskData=root:DSM:MaskData
1255        Variable pt,ii,endPt
1256        endPt=numpnts(MaskData)
1257        pt = pcsr(A)
1258        for(ii=pt;ii<endPt;ii+=1)
1259                // toggle NaN (keep) or Data value (= masked)
1260                MaskData[ii] = (numType(MaskData[ii])==0) ? NaN : data[ii]              //if NaN, doesn't plot
1261        endfor
1262        return(0)
1263End
1264
1265Function CleanUpJunk()
1266        // clean up the old junk on the graph, /Z for no error
1267        Execute "RemoveExtrapolation()"
1268        Execute "RemoveDesmeared()"
1269        Execute "RemoveSmoothed()"
1270        Execute "RemoveMask()"
1271       
1272        //remove the cursor
1273        Cursor/K A
1274       
1275        //always initialize these
1276        String/G root:DSM:gStr1 = ""
1277        String/G root:DSM:gStr2 = ""
1278        String/G root:DSM:gIterStr = ""
1279       
1280        // clean up the old waves from smoothing and desmearing steps
1281        SetDataFolder root:DSM
1282        Killwaves/Z I_smth,I_dsm,I_dsm_sm,Q_smth,Q_dsm,S_smth,S_dsm,Yi_SS,Yq_SS
1283        Killwaves/Z Weights,FW,R_wave,S_wave,H_wave,Di,Ci
1284        Killwaves/Z Is_old,I_old,err
1285        Killwaves/Z Q_ext,I_ext,S_ext,hqq,hqi,lqq,lqi
1286        Killwaves/Z MaskData,Q_msk,I_msk,S_msk
1287        Killwaves/Z Q_work,I_work,S_work                                //working waves for desmearing step
1288        SetDataFolder root:
1289End
1290
1291// does not alter the data sets - reports a power law
1292// exponent and makes it global so it will automatically
1293// be used during the desmearing
1294//
1295// generates extr_hqi vs extr_hqq that are Appended to the graph
1296Function DSM_ExtrapolateButtonProc(ctrlName) : ButtonControl
1297        String ctrlName
1298
1299        NVAR nend = root:DSM:gNptsExtrap
1300        NVAR m_pred = root:DSM:gPowerM
1301       
1302        SetDataFolder root:DSM
1303//use masked data if it exists
1304        if(exists("I_msk")==1 && exists("Q_msk")==1 && exists("S_msk")==1)
1305                wave Qw = root:DSM:Q_msk
1306                Wave Iw = root:DSM:I_msk
1307                Wave Sw = root:DSM:S_msk
1308        else
1309                //start from the "_exp" waves
1310                if(exists("I_exp")==1 && exists("Q_exp")==1 && exists("S_exp")==1)
1311                        wave Qw = root:DSM:Q_exp
1312                        Wave Iw = root:DSM:I_exp
1313                        Wave Sw = root:DSM:S_exp
1314                endif
1315        endif
1316        SetDataFolder root:
1317       
1318        DSM_SetExtrWaves(Qw)
1319        m_pred = DSM_DoExtrapolate(Qw,Iw,Sw,nend)
1320        AppendExtrapolation()
1321       
1322        return(0)
1323End
1324
1325//smooths the data in steps as requested...
1326//
1327// typically do a simple Box smooth first,
1328// then do a smoothing spline, keeping the same number of data points
1329//
1330// chi-squared is reported - so you can see how "bad" the smoothing is
1331// smoothing of the data is done on a log(I) scale if requested
1332// setting doLog variable to 0 will return to linear smoothing
1333// (I see little difference)
1334//
1335// start from the "_smth" waves if they exist
1336// otheriwse start from the working waves
1337//
1338// create Q_smth, I_smth, S_smth
1339// keep track of the smoothing types and passes in the wave note
1340//
1341Function DSM_SmoothButtonProc(ctrlName) : ButtonControl
1342        String ctrlName
1343
1344        SetDataFolder root:DSM 
1345
1346        Variable ii,new_n,pass,nq_ext,offset,doLog=1
1347        String noteStr=""
1348       
1349        // want log scaling of intensity during smoothing?
1350        ControlInfo DSMControl_3g
1351        doLog = V_value
1352                       
1353////    if(exists("I_smth")==1 && exists("Q_smth")==1 && exists("S_smth")==1 && ! freshMask)
1354        if(exists("I_smth")==1 && exists("Q_smth")==1 && exists("S_smth")==1)
1355                //start from the smoothed waves --- just need the wave references
1356        else
1357                //start from the "msk", creating smth waves
1358                if(exists("I_msk")==1 && exists("Q_msk")==1 && exists("S_msk")==1)
1359                        wave Q_msk = root:DSM:Q_msk
1360                        Wave I_msk = root:DSM:I_msk
1361                        Wave S_msk = root:DSM:S_msk
1362                        Duplicate/O I_msk,I_smth,Q_smth,S_smth
1363                        I_smth = I_msk
1364                        Q_smth = Q_msk
1365                        S_smth = S_msk
1366                else
1367                        //start from the "_exp" waves
1368                        if(exists("I_exp")==1 && exists("Q_exp")==1 && exists("S_exp")==1)
1369                                wave Q_exp = root:DSM:Q_exp
1370                                Wave I_exp = root:DSM:I_exp
1371                                Wave S_exp = root:DSM:S_exp
1372                                Duplicate/O I_exp,I_smth,Q_smth,S_smth
1373                                I_smth = I_exp
1374                                Q_smth = Q_exp
1375                                S_smth = S_exp
1376                        endif
1377                endif
1378        endif
1379       
1380        wave Q_smth = root:DSM:Q_smth
1381        Wave I_smth = root:DSM:I_smth
1382        Wave S_smth = root:DSM:S_smth
1383
1384        // extend the data to avoid end effects
1385        //creates I_ext,Q_ext,S_ext with extended q-values
1386        // fit 15 pts at each end, typically add 10 pts to each end
1387        // does nothing if num_extr = 0
1388        ControlInfo/W=Desmear_Graph DSMControl_3d               //extendCheck
1389        if(V_value == 1)
1390                ExtendToSmooth(Q_smth,I_smth,S_smth,15,15,10)
1391        else
1392                ExtendToSmooth(Q_smth,I_smth,S_smth,15,15,0)            //don't extend, just use the original data
1393        endif
1394       
1395        //whether extending or not, the working data is "_ext", set by ExtendToSmooth()
1396        setDataFolder root:DSM
1397        wave Q_ext = root:DSM:Q_ext
1398        Wave I_ext = root:DSM:I_ext
1399        Wave S_ext = root:DSM:S_ext
1400       
1401        noteStr=note(I_smth)
1402        Note I_ext , noteStr
1403       
1404        WaveStats/Q I_ext
1405        if(V_min<=0)
1406                Print "negative itensity found, using linear scale for smoothing"
1407                doLog = 0
1408        endif
1409       
1410        if(doLog)
1411                //convert to log scale
1412                Duplicate/O I_ext root:DSM:I_log,root:DSM:I_log_err
1413                Wave I_log = root:DSM:I_log
1414                Wave I_log_err = root:DSM:I_log_err
1415                I_log = log(I_ext)
1416                WaveStats/Q I_log
1417                offset = 2*(-V_min)             
1418                I_log += offset
1419                I_log_err = S_ext/(2.30*I_ext)
1420                I_ext = I_log
1421        endif   
1422        //
1423       
1424        ControlInfo/W=Desmear_Graph DSMControl_3b               //BoxCheck
1425        if(V_value == 1)
1426                //do a simple Box smooth first - number of points does not change
1427                // fills ends with neighboring value (E=3) (n=3 points in smoothing window)     
1428                Smooth/E=3/B 3, I_ext
1429               
1430                noteStr=note(I_ext)
1431                pass = NumberByKey("BOX", noteStr, "=", ";")
1432                noteStr = ReplaceNumberByKey("BOX", noteStr, pass+1, "=", ";")
1433                Note/K I_ext
1434                Note I_ext , noteStr
1435//              Print "box = ",noteStr
1436        endif
1437
1438        NVAR sParam = root:DSM:gSmoothFac
1439       
1440        ControlInfo/W=Desmear_Graph DSMControl_3c               //SSCheck
1441        if(V_value == 1)
1442//              Interpolate2/T=3/N=(nq)/I=2/F=1/SWAV=S_ext/Y=Yi_SS/X=Yq_SS Q_ext, I_ext         //Igor 5
1443//              Interpolate/T=3/N=(nq)/I=2/F=1/S=S_ext/Y=Yi_SS/X=Yq_SS I_ext /X=Q_ext                   //Igor 4
1444        //Igor 4
1445                String str=""
1446                nq_ext = numpnts(Q_ext)
1447                str = "Interpolate/T=3/N=("+num2str(nq_ext)+")/I=1/F=("+num2str(sParam)+")/Y=Yi_SS/X=Yq_SS I_ext /X=Q_ext"             
1448                Execute str
1449//              Print Str
1450//              Interpolate2/T=3/N=(nq_ext)/I=1/F=(sParam)/Y=Yi_SS/X=Yq_SS Q_ext, I_ext
1451        // end Igor 4   
1452                wave yi_ss = root:DSM:yi_ss
1453                wave yq_ss = root:DSM:yq_ss
1454                // reassign the "working" waves with the result of interpolate, which has the same I/Q values
1455                I_ext = yi_ss
1456                Q_ext = yq_ss
1457               
1458                noteStr=note(I_ext)
1459                pass = NumberByKey("SPLINE", noteStr, "=", ";")
1460                noteStr = ReplaceNumberByKey("SPLINE", noteStr, pass+1, "=", ";")
1461                Note/K I_ext
1462                Note I_ext , noteStr
1463//              Print "spline = ",noteStr
1464        endif
1465       
1466        //undo the scaling
1467        If(doLog)
1468                I_ext -= offset
1469                I_ext = 10^I_ext
1470        endif
1471       
1472        // at this point, I_ext has too many points - and we need to just return the
1473        // center chunk that is the original q-values
1474        // as assign this to the working set for the desmear step
1475        // and to the _smth set in case another smoothng pass is desired
1476        // Q_smth has not changed, S_smth has not changed
1477        I_smth = interp(Q_smth[p], Q_ext, I_ext )
1478       
1479        Note/K I_smth
1480        Note I_smth , noteStr
1481//      Print "end of smoothed, note = ",note(I_smth)
1482       
1483        Variable nq = numpnts(root:DSM:Q_smth)
1484//      Print "nq after smoothing = ",nq
1485
1486        //reset the global
1487        NVAR gNq = root:DSM:gNq
1488        gNq = nq
1489        //report the chi^2 difference between the smoothed curve and the experimental data
1490        NVAR chi2 = root:DSM:gChi2Smooth
1491        chi2 = SmoothedChi2(I_smth)
1492       
1493        Execute "AppendSmoothed()"
1494       
1495        setDataFolder root:
1496        return(0)
1497End
1498
1499// undo the smoothing and start over - useful if you've smoothed
1500// too aggressively and need to back off
1501Function DSM_SmoothUndoButtonProc(ctrlName) : ButtonControl
1502        String ctrlName
1503       
1504        Execute "RemoveSmoothed()"
1505        SetDataFolder root:DSM
1506        Killwaves/Z I_smth,Q_smth,S_smth,Q_ext,I_ext,S_ext,Yi_SS,Yq_SS
1507        SetDataFolder root:
1508        return (0)
1509end
1510
1511//calculate the chi^2 value for the smoothed data
1512Function SmoothedChi2(I_smth)
1513        Wave I_smth
1514       
1515        //start from the "msk", if they exist
1516        if(exists("I_msk")==1 && exists("Q_msk")==1 && exists("S_msk")==1)
1517                Wave iw = root:DSM:I_msk
1518                Wave sw = root:DSM:S_msk
1519        else
1520                //start from the "_exp" waves
1521                if(exists("I_exp")==1 && exists("Q_exp")==1 && exists("S_exp")==1)
1522                        Wave iw = root:DSM:I_exp
1523                        Wave sw = root:DSM:S_exp
1524                endif
1525        endif
1526       
1527        Variable ii,chi2,num
1528        chi2=0
1529        num = numpnts(iw)
1530        for(ii=0;ii<num;ii+=1)
1531                chi2 += (iw[ii] - I_smth[ii])^2 / (sw[ii])^2
1532        endfor
1533        Chi2 /= (num-1)
1534        return (chi2)
1535end
1536
1537// I_dsm is the desmeared data
1538//
1539// step (7) - desmearing is done, write out the result
1540//
1541Function DSM_SaveButtonProc(ctrlName) : ButtonControl
1542        String ctrlName
1543
1544        String saveStr
1545        SVAR curFile = root:DSM:gCurFile
1546        saveStr = CleanupName((curFile + "_dsm"),0)             //the output filename
1547        //
1548        WriteUSANSDesmeared(saveStr,0,0,1)                      //use the full set (lo=hi=0) and present a dialog
1549       
1550        SetDataFolder root:
1551        return(0)
1552End
1553
1554Function DSM_HelpButtonProc(ctrlName) : ButtonControl
1555        String ctrlName
1556
1557        DisplayHelpTopic/K=1 "Desmearing USANS Data"
1558        return(0)
1559End
1560
1561
1562//toggles the log/lin display of the loaded data set
1563Function DSM_LoadCheckProc(ctrlName,checked) : CheckBoxControl
1564        String ctrlName
1565        Variable checked
1566
1567        ModifyGraph log=(checked)
1568        return(0)
1569End
1570
1571
1572Function DSM_ExtrapolationCheckProc(ctrlName,checked) : CheckBoxControl
1573        String ctrlName
1574        Variable checked
1575
1576        if(checked)
1577                AppendExtrapolation()
1578        else
1579                RemoveExtrapolation()
1580        endif
1581        return(0)
1582End
1583
1584
1585// takes as input the "working waves"
1586//
1587// creates intermediate waves to work on
1588//
1589// output of Q_dsm,I_dsm,S_dsm (and I_dsm_sm, the result of smearing I_dsm)
1590//
1591Function DSM_DesmearButtonProc(ctrlName) : ButtonControl
1592        String ctrlName
1593       
1594        SetDataFolder root:DSM
1595        if(exists("I_smth")==1 && exists("Q_smth")==1 && exists("S_smth")==1)
1596                wave Q_smth = root:DSM:Q_smth
1597                Wave I_smth = root:DSM:I_smth
1598                Wave S_smth = root:DSM:S_smth
1599                Duplicate/O I_smth,I_work,Q_work,S_work
1600                I_work = I_smth
1601                Q_work = Q_smth
1602                S_work = S_smth
1603        else
1604                //start from the "msk", creating work waves
1605                if(exists("I_msk")==1 && exists("Q_msk")==1 && exists("S_msk")==1)
1606                        wave Q_msk = root:DSM:Q_msk
1607                        Wave I_msk = root:DSM:I_msk
1608                        Wave S_msk = root:DSM:S_msk
1609                        Duplicate/O I_msk,I_work,Q_work,S_work
1610                        I_work = I_msk
1611                        Q_work = Q_msk
1612                        S_work = S_msk
1613                else
1614                        //start from the "_exp" waves
1615                        if(exists("I_exp")==1 && exists("Q_exp")==1 && exists("S_exp")==1)
1616                                wave Q_exp = root:DSM:Q_exp
1617                                Wave I_exp = root:DSM:I_exp
1618                                Wave S_exp = root:DSM:S_exp
1619                                Duplicate/O I_exp,I_work,Q_work,S_work
1620                                I_work = I_exp
1621                                Q_work = Q_exp
1622                                S_work = S_exp
1623                        endif
1624                endif
1625        endif
1626        SetDataFolder root:
1627       
1628        NVAR nq = root:DSM:gNq
1629        NVAR m = root:DSM:gPowerM
1630        NVAR chi2_target = root:DSM:gChi2Target
1631        NVAR maxFastIter = root:DSM:gMaxFastIter
1632        NVAR maxSlowIter = root:DSM:gMaxSlowIter
1633       
1634        //      SET WEIGHTING OF EXPERIMENTAL DATA.
1635        Duplicate/O Q_work root:DSM:weights
1636        Wave weights = root:DSM:weights
1637        weights = 1/S_work^2
1638
1639//      calculate weighting array for smearing of data
1640        Make/O/D/N=(nq,nq) root:DSM:FW
1641        Wave FW = root:DSM:FW
1642        Weights_L(m,FW,Q_work)
1643
1644//      ^^^^   Iterative desmearing   ^^^^*
1645        Variable chi2_old,chi2_new,done,iter
1646//      FOR 0TH ITERATION, EXPERIMENTAL DATA IS USED FOR Y_OLD, create ys_old for result
1647//      y_old = I_old, y_new = I_dsm, I_dsm_sm = ys_new,
1648// duplicate preserves the wave note!
1649        Duplicate/O I_work root:DSM:I_old,root:DSM:Is_old,root:DSM:I_dsm,root:DSM:I_dsm_sm
1650        Duplicate/O Q_work root:DSM:Q_dsm,root:DSM:S_dsm                //sets Q_dsm correctly
1651        wave S_dsm = root:DSM:S_dsm             //set correctly at end of this function
1652        wave I_old = root:DSM:I_old
1653        wave Is_old = root:DSM:Is_old
1654        wave I_dsm = root:DSM:I_dsm
1655        wave I_dsm_sm = root:DSM:I_dsm_sm
1656        I_old = I_work
1657        Is_old = 0
1658        I_dsm = 0
1659        I_dsm_sm = 0
1660       
1661//      Smear 0TH iter guess Y --> YS
1662        DSM_Smear(I_old,Is_old,NQ,FW)
1663        chi2_OLD = chi2(Is_old,I_work,weights,NQ)
1664       
1665        printf "starting chi^2 = %g\r",chi2_old
1666        Print "Starting fast convergence..... "
1667
1668        done = 0                //false
1669        iter = 0
1670        Do                                      // while not done, use Fast convergence
1671                I_dsm = I_work * I_old / Is_old //        Calculate corrected guess (no need for do-loop)
1672
1673//        Smear iter guess I_dsm --> I_dsm_sm and see how well I_dsm_sm matches experimental data
1674                DSM_Smear(I_dsm,I_dsm_sm,NQ,FW)
1675                chi2_new = chi2(I_dsm_sm,I_work,weights,NQ)
1676
1677//        Stop iteration if fit from new iteration has worse fit...
1678                If (chi2_new > chi2_old)
1679                        Done = 1
1680                Endif
1681               
1682//        Stop iteration if fit is better than target value...
1683                If (chi2_new < chi2_target)
1684                        Done = 1
1685                Else
1686                        Iter += 1
1687                        Printf "Iteration %d, Chi^2(new) = %g\r", Iter,chi2_new
1688                        I_old = I_dsm
1689                        Is_old = I_dsm_sm
1690                        if(iter>maxFastIter)
1691                                break
1692                        endif
1693                        CHI2_OLD = CHI2_NEW
1694                EndIf
1695        while( !done )
1696       
1697        // append I_dsm,I_exp to the graph
1698        Execute "AppendDesmeared()"
1699        DoUpdate
1700
1701// step (6) - refine the desmearing using slow convergence
1702        Print "Starting slow convergence..... "
1703        done = 0                //  ! reset flag for slow convergence
1704        Do
1705                I_dsm = I_old + (I_work - Is_old)       //        Calculate corrected guess
1706
1707//        Smear iter guess Y --> YS
1708                DSM_Smear(I_dsm,I_dsm_sm,NQ,FW)
1709                chi2_new = chi2(I_dsm_sm,I_work,weights,NQ)
1710
1711//        Stop iteration if fit from new iteration has worse fit...
1712                If (chi2_new > chi2_old)
1713                        Done = 1
1714                Endif
1715               
1716//        Stop iteration if fit is better than target value...
1717                If (chi2_new < chi2_target)
1718                        Done = 1
1719                Else
1720                        Iter += 1
1721                       
1722                        if(mod(iter, 50 ) ==0 )
1723                                Printf "Iteration %d, Chi^2(new) = %g\r", Iter,chi2_new
1724                                DoUpdate
1725                        endif
1726                        I_old = I_dsm
1727                        Is_old = I_dsm_sm
1728
1729                        CHI2_OLD = CHI2_NEW
1730                EndIf
1731                if(iter>maxSlowIter)
1732                        break
1733                endif
1734        While ( !done )
1735
1736        Printf "Iteration %d, Chi^2(new) = %g\r", Iter,chi2_new
1737
1738        // adjust the error
1739        SetDataFolder root:DSM
1740        Duplicate/O S_work err
1741        S_dsm = err/I_work*I_dsm                //proportional error
1742       
1743        NVAR gChi = root:DSM:gChi2Final         //chi^2 final
1744        NVAR gIter = root:DSM:gIterations               //total number of iterations
1745        gChi = chi2_new
1746        gIter = Iter
1747       
1748        SVAR str = root:DSM:gIterStr
1749        sprintf str,"%d iterations required",iter
1750       
1751        SetDataFolder root:
1752        return(0)
1753End
1754
1755Function DSM_RevertButtonProc(ctrlName) : ButtonControl
1756        String ctrlName
1757
1758        CleanUpJunk()
1759         
1760        //reset the working waves to the original
1761        wave Q_exp_orig = root:DSM:Q_exp_orig
1762        wave I_exp_orig = root:DSM:I_exp_orig
1763        wave S_exp_orig = root:DSM:S_exp_orig
1764
1765        Duplicate/O Q_exp_orig root:DSM:Q_exp
1766        Duplicate/O I_exp_orig root:DSM:I_exp
1767        Duplicate/O S_exp_orig root:DSM:S_exp
1768        // kill the data folder
1769        // re-initialize?
1770
1771        return(0)
1772End
Note: See TracBrowser for help on using the repository browser.