Changes in D16 macro suggested by Stef to handle correctly the case where a Cd...

Changes in D16 macro suggested by Stef to handle correctly the case where a Cd is given, but not an empty cell measurement

lamp_mac/D16/d16_sheet_normTime.prox

@@ -15,11 +15,11 @@
 ;;; the same and they must correspond to the same 'gamma' acquisitions to be consistent.
 ;;; This macro will not check this and results will be incorrect if the numors given
 ;;; for the different scans are not consistent.
-;;; 
+;;;
 ;;; The macro assumes as default that the beam center is at the center of the Miland detector,
 ;;; i.e. at pixel position 160,160 (note that the position of the lower left pixel is taken as 0.5,0.5
 ;;; and the position of the upper right pixel is taken as 319.5,319.5). However if a
-;;; different center has been previously calculated it can be given using the Center cell. 
+;;; different center has been previously calculated it can be given using the Center cell.
 ;;;
 ;;; The numors will be read using the RDOPR command in Lamp, so if a calibration file is
 ;;; needed it has to be loaded into Lamp before calling the macro. Use the Calibration files
@@ -28,7 +28,7 @@
 ;;; The macro does the following operations:
 ;;;
 ;;;    1. Reads a scan as a 3D workspace, e.g. w1[320, 320, N_acquisitions], where N_acquisitions
-;;;       corresponds to the number of gamma positions measured. If a calibration file has been 
+;;;       corresponds to the number of gamma positions measured. If a calibration file has been
 ;;;       previously loaded, each numor is divided by the calibration file.
 ;;;
 ;;;    2. Combines scans corresponding to the same gamma position and sorts them by gamma value
@@ -46,7 +46,7 @@
 ;;;
 ;;;       If a normalization (water/vanadium) run has been given, its background is subtracted in
 ;;;       a similar way and the corrected intensity stored in w16.
-;;;       
+;;;
 ;;;       The needed transmissions must have been calculated previously and be given in
 ;;;       the corresponding T_S, T_BgrS, T_V, and T_BgrV cells.
 ;;;       Any transmission value not given is set to 1 by default.
@@ -57,7 +57,7 @@
 ;;;       Default values for the water thickness (0.1 cm) and cross section (0.87 cm-1) are set in the
 ;;;       variables d and c, respectively. The user can modify those values if needed, or set them to 1
 ;;;       together with the sample thickness to skip this scaling.
-;;;       
+;;;
 ;;;
 ;;;    6. Scale the result by multiplying each of the N_acquisitions by a factor computed as the
 ;;;       mean of the intensity of two numors in the overlapping region. The first numor is taken
@@ -69,7 +69,7 @@
 ;;;    8. Reduce the number of points and improve the statistics of each point by averaging
 ;;;       N points together, where N is the value given in the Grouping cell --> w21.
 ;;;
-;;; The macro will go through all the lines and for each line it will save the content of 
+;;; The macro will go through all the lines and for each line it will save the content of
 ;;; the following workspaces:
 ;;;
 ;;;    w16 = Sample intensity (integrated over 2-theta and corrected for background) for each gamma position
@@ -78,32 +78,26 @@
 ;;;    w20 = Sample normalized intensity vs 2-theta (all gamma positions combined)
 ;;;    w21 = Sample normalized intensity vs 2-theta (after grouping)
 ;;;
-
 ;sets workspace manipulation mode to raw --> operations on workspaces are not extended to associated variables
 setmanip, /raw
-
 ;Water thickness (in cm) and cross section (in cm-1)
-d = 0.1    
-c = 0.87 
-
+d = 0.1
+c = 0.87
 ;Initialize wokspaces when treating first line
 if (cell.Firstime eq 1) then begin & w_clear,w=1 & w_clear,w=2 & w_clear,w=3 & w_clear,w=4 & w_clear,w=5 & endif
 if (cell.Firstime eq 1) then begin & w_clear,w=12 & w_clear,w=12 & w_clear,w=13 & w_clear,w=14 & w_clear,w=15 & endif
 if (cell.Firstime eq 1) then begin & w_clear,w=16 & w_clear,w=17 & w_clear,w=18 & w_clear,w=19 & w_clear,w=20 & w_clear,w=21 & endif
-
 ;Set coordinates for center
 if (cell.Center gt '')  then a = float(strsplit(cell.Center, ',', /extract))
 if (n_elements(a) ne 2) then a = [160.0,160.0]
-
 ;Sample thickness: Set to 1 if empty
 if (cell.d_S gt '') then e = float(cell.d_S)
 if n_elements(e) eq 0 then e = 1.0
-
 ;Default transmissions = 1
-if (cell.T_S gt '')    then s = float(cell.T_S)    
-if (cell.T_BgrS gt '') then t = float(cell.T_BgrS) 
-if (cell.T_V gt '')    then v = float(cell.T_V)    
-if (cell.T_BgrV gt '') then w = float(cell.T_BgrV) 
+if (cell.T_S gt '')    then s = float(cell.T_S)
+if (cell.T_BgrS gt '') then t = float(cell.T_BgrS)
+if (cell.T_V gt '')    then v = float(cell.T_V)
+if (cell.T_BgrV gt '') then w = float(cell.T_BgrV)
 if n_elements(s) eq 0 then s = 1.0
 if n_elements(t) eq 0 then t = 1.0
 if n_elements(v) eq 0 then v = 1.0
@@ -112,77 +106,68 @@ if ((s lt 0.0) or (s gt 1.0)) then s = 1.0
 if ((t lt 0.0) or (t gt 1.0)) then t = 1.0
 if ((v lt 0.0) or (v gt 1.0)) then v = 1.0
 if ((w lt 0.0) or (w gt 1.0)) then w = 1.0
-
 ;print, 'Sample ', cell.Sn, ' center = ', a, ' transm = ', s, t, v, w
-
-;Sample 
-w1 = rdopr(cell.Sample) & w1 = sort_d16(w1) & w11 = straight_d16(w1, center=a, norm='time') 
+;Sample
+w1 = rdopr(cell.Sample) & w1 = sort_d16(w1) & w11 = straight_d16(w1, center=a, norm='time')
 y_tit[11] = 'Gamma scan'
-
 ;Background for sample
-if (cell.Bgr_S gt '') then w2 = rdopr(cell.Bgr_S) 
-if (cell.Bgr_S gt '') then w2 = sort_d16(w2) 
+if (cell.Bgr_S gt '') then w2 = rdopr(cell.Bgr_S)
+if (cell.Bgr_S gt '') then w2 = sort_d16(w2)
 if (cell.Bgr_S gt '') then w12 = straight_d16(w2, center=a, norm='time')
 if (cell.Bgr_S gt '') then y_tit[12] = 'Gamma scan'
-
 ;Vanadium
-if (cell.Vana gt '') then w3 = rdopr(cell.Vana)  
-if (cell.Vana gt '') then w3 = sort_d16(w3) 
-if (cell.Vana gt '') then w13 = straight_d16(w3, center=a, norm='time') 
+if (cell.Vana gt '') then w3 = rdopr(cell.Vana)
+if (cell.Vana gt '') then w3 = sort_d16(w3)
+if (cell.Vana gt '') then w13 = straight_d16(w3, center=a, norm='time')
 if (cell.Vana gt '') then y_tit[13] = 'Gamma scan'
-
 ;Background for vanadium
-if (cell.Bgr_V gt '') then w4 = rdopr(cell.Bgr_V) 
-if (cell.Bgr_V gt '') then w4 = sort_d16(w4) 
-if (cell.Bgr_V gt '') then w14 = straight_d16(w4, center=a, norm='time') 
+if (cell.Bgr_V gt '') then w4 = rdopr(cell.Bgr_V)
+if (cell.Bgr_V gt '') then w4 = sort_d16(w4)
+if (cell.Bgr_V gt '') then w14 = straight_d16(w4, center=a, norm='time')
 if (cell.Bgr_V gt '') then y_tit[14] = 'Gamma scan'
-
 ;Cadmium
-if (cell.Cd gt '') then w5 = rdopr(cell.Cd)    
-if (cell.Cd gt '') then w5 = sort_d16(w5) 
-if (cell.Cd gt '') then w15 = straight_d16(w5, center=a, norm='time') 
+if (cell.Cd gt '') then w5 = rdopr(cell.Cd)
+if (cell.Cd gt '') then w5 = sort_d16(w5)
+if (cell.Cd gt '') then w15 = straight_d16(w5, center=a, norm='time')
 if (cell.Cd gt '') then y_tit[15] = 'Gamma scan'
-
 ;Compare dimensions of ancillary measurements with those of sample & set to zero any wks with different dimensions
-
 lamp_journal,'-------------------------------------------------------------------------------------'
-
 if n_elements(w12) ne n_elements(w11) then w12 = w11 * 0.0
 if n_elements(w12) ne n_elements(w11) then e12 = w11 * 0.0
-if total(w12) eq 0 then lamp_journal, 'D16 Prox_sheet for Sample ' + cell.Sn + ': No background measurement.' 
-
+if total(w12) eq 0 then lamp_journal, 'D16 Prox_sheet for Sample ' + cell.Sn + ': No background measurement.'
+if total(w12) eq 0 then print, 'D16 Prox_sheet for Sample ' + cell.Sn + ': No background measurement.'
 if n_elements(w13) ne n_elements(w11) then w13 = 0
 if n_elements(w13) ne n_elements(w11) then e13 = 0
 if total(w13) eq 0 then lamp_journal, 'D16 Prox_sheet for Sample ' + cell.Sn + ': No vanadium measurement.'
-
+if total(w13) eq 0 then print, 'D16 Prox_sheet for Sample ' + cell.Sn + ': No vanadium measurement.'
 if n_elements(w14) ne n_elements(w11) then w14 = w11 * 0.0
 if n_elements(w14) ne n_elements(w11) then e14 = w11 * 0.0
-if (total(w14) eq 0) and (total(w13) gt 0) then lamp_journal, 'D16 Prox_sheet for Sample ' + cell.Sn + ': No background for vanadium.' 
-
+if total(w14) eq 0 then lamp_journal, 'D16 Prox_sheet for Sample ' + cell.Sn + ': No background for vanadium.'
+if total(w14) eq 0 then print, 'D16 Prox_sheet for Sample ' + cell.Sn + ': No background for vanadium..'
 if n_elements(w15) ne n_elements(w11) then w15 = w11 * 0.0
 if n_elements(w15) ne n_elements(w11) then e15 = w11 * 0.0
 if total(w15) eq 0 then begin & lamp_journal,' ' & lamp_journal, 'D16 Prox_sheet for Sample ' + cell.Sn + ': No cadmium measurement.' & lamp_journal,' ' & endif
-
+if total(w15) eq 0 then print, 'D16 Prox_sheet for Sample ' + cell.Sn + ': No cadmium measurement.'
 lamp_journal,'-------------------------------------------------------------------------------------'
-
-;Subtract background and Cd (sample)
-w16 = ((w11-w15)/(s*t)) - ((w12-w15)/t) 
-e16 = sqrt( (e11/(s*t))^2 + (e15/(s*t))^2 + (e12/t)^2 + (e15/t)^2 )
-
-;Subtract background and Cd (vanadium, only if contains data)
-if n_elements(w13) gt 1 then w17 = ((w13-w15)/(v*w)) - ((w14-w15)/w) 
-if n_elements(w13) gt 1 then e17 = sqrt( (e13/(s*t))^2 + (e15/(s*t))^2 + (e14/t)^2 + (e15/t)^2 )
-if n_elements(w13) gt 1 then x17 = x13 
-
+;Subtract Cd to All
+w11 = w11 - w15 & e11 = sqrt(e11^2 + e15^2)
+if total(w12) gt 1e-9 then begin & w12 = w12 - w15 & e12 = sqrt(e12^2 + e15^2) & endif
+if total(w13) gt 1e-9 then begin & w13 = w13 - w15 & e13 = sqrt(e13^2 + e15^2) & endif
+if total(w14) gt 1e-9 then begin & w14 = w14 - w15 & e14 = sqrt(e14^2 + e15^2) & endif
+;Subtract background (sample)
+w16 = (w11/(s*t)) - (w12/t)
+e16 = sqrt( (e11/(s*t))^2 + (e12/t)^2)
+;Subtract background (vanadium, only if contains data)
+if n_elements(w13) gt 1 then w17 = (w13/(v*w)) - (w14/w)
+if n_elements(w13) gt 1 then e17 = sqrt( (e13/(s*t))^2 + (e14/t)^2)
+if n_elements(w13) gt 1 then x17 = x13
 ;Normalise by vanadium (only if contains data)
 w18 = w16
-if n_elements(w17) gt 1 then w18 = (c*d/e) * w16/w17 
+if n_elements(w17) gt 1 then w18 = (c*d/e) * w16/w17
 if n_elements(w17) gt 1 then e18 = sqrt( ((c*d/e)*e16/w17)^2 + ((c*d/e)*w16*e17/(w17^2))^2 )
 if n_elements(w17) gt 1 then x18 = x16
-
 ;Scale consecutive scans so that they overlap
 w19 = scale_d16(w18)
-
 ;Create a single vector and group points
 w20 = reform(w19, n_elements(w19))
 e20 = reform(e19, n_elements(e19))
@@ -191,18 +176,16 @@ w20 = w20[sort(x20)]
 e20 = e20[sort(x20)]
 x20 = x20[sort(x20)]
 y_tit[20] = 'Intensity'
-
 ;Number of theta points to group
-if (cell.Grouping gt '') then n = fix(cell.Grouping)    
+if (cell.Grouping gt '') then n = fix(cell.Grouping)
 if n_elements(n) eq 0 then n = 1
-
 w21 = groupx(w20,n) / float(n) & e21 = e21 / float(n)
-
-;Convert x-axis to Q
+;Convert x-axis to Q'
 l = p21[8] ;wavelength
 if l gt 3.0 and l lt 8.0 then x21 = (4*!pi/l)*sin(x21*!pi/360.0)
 if l gt 3.0 and l lt 8.0 then x_tit[21] = 'Q/AA-1'
-
-;Write grouped data into an ASCII file
+;print,'Write grouped data into an ASCII file'
 z = 'D16_' + strcompress(cell.SId, /remove_all) + '.dat'
-write_lamp, z, W=21, format = 'column'
\ No newline at end of file
+write_lamp, z, W=21, format = 'column'
+;
+print,'<----- END ----->'