;+ ; Project : HESI ; ; Name : GD ; ; Purpose : This procedure either reads or creates grid parameters. ; It reads them from the specified grid parameter file, gdfile, ; or takes a few input values, number of grids, grid spacing, ; angular resolution log-spacing, etc, to create a grid ; parameter file. ; ; Category : HESI ; ; Explanation : ; If an gdfile is specified, call read_gdfile_gd to read in ; values needed for grid dimensions, instead of calculating them! ; Otherwise, the program will calculate it's own values based on ; formulae given by Brian Dennis, GSFC ; These original formulae are based on a Tungsten Grid. ; ; Use : ; ; Inputs : ; ; Opt. Inputs : None ; ; Calling Sequence: ; gd, xx, grid_sep, pitch, e_diff_10, e_diff_50, e_diff_fl, e_peak,$ ; thickness, gridnm, density,$ ; e_diff_manual, man_rmod_amp, res, subs, ngrids, fine_res, coarse_res, subarr ; Outputs : ; XX - cross_section arrays vs energy ; Grid_sep- grid separation distance in mm ; Pitch - Array of grid pitches in mm ; E_diff_10 - energies corresponding to 1st 3 harmonic numbers and 90% modulation. ; E_diff_50 - energies corresponding to 1st 3 harmonic numbers and 50% modulation. ; E_diff_fl - energies corresponding to 1st 3 harmonic numbers and 0 modulation. ; E_peak - Energy (keV) of peak of diffraction amp. curve ; Thickness - Array of grid thicknesses in mm ; Gridnm - collimator number, starts at 1 ; Density - Array of Density of grid materials in gm/cm3 ; E_diff_manual- energies corresponding to 1st 3 harmonic numbers and Man_rmod_amp. ; Man_rmod_amp - If set then input relative modulation amplitude, ; using Rmod_amp as a percentage ; Res - angular resolution of each grid pair, unused. ; Subs- Grid materials, string ; ngrids- total number of grids ; Fine_res- finest resolution in arcseconds ; Coarse_res- coarsest resolution in arcseconds ; Subarr- Grid materials, string ; ; Opt. Outputs: None ; ; Keywords : GDFILE- Text file to be read by read_gdfile_gd.pro with ; grid parameters. ; OUT - Text file to be created with grid parameters. ; Rmod_amp - Find the energies corresponding this percentage ; of the diffraction amplitude curve and 1st-3rd harmonics. ; Calls : ; ; Common : None ; ; Restrictions: ; ; Side effects: None. ; ; Prev. Hist : Written by Eric Carzon, Hughes/STX, May 1993 ; ; Modified : Version 2 RAS, 16-apr-1997 ; Modified : Version 2 RAS, 13-oct-1999. ;- Pro gd, xx, grid_sep, pitch, e_diff_10, e_diff_50, e_diff_fl, e_peak,$ thickness, gridnm, density,$ e_diff_manual, man_rmod_amp, res, subs, ngrids, fine_res, coarse_res, subarr,$ cgrid_sep=cgrid_sep,rmod_amp=rmod_amp,gdfile=gdfile,out=out man_rmod_amp = 0 ; initial value of allowing user specified ; reduction in modulation ; amplititue, program is set to do 10,50, and 90% redu. if keyword_set(out) then outfile = out else outfile ='' ; WRITE HESI STRAWMAN GRID CHARACTERISTICS rim = 10 ; if (keyword_set(gdfile)) then gdfile = gdfile else gdfile = '' ; ; If using a file to specify strawman grid characteristics, read it ; and compare density and substance of file with that of density and ; substance retrieved from XCOM.DAT file or from XCOM program ; if (gdfile ne '') then begin read_gdfile_gd,gdfile,grid_sep,db,dt,dfov,frac,subarr,gdensity,$ gridnm,res,pitch,sw,thickness,fov,ngrids,fine_res,coarse_res,tfset XX = fltarr(550,8,ngrids) denarr = fltarr(ngrids) process_xcom_gd,subarr(0),X,density,subs xlim = size(X) xx(0:xlim(1)-1,0:*,0) = X(0:*,0:*) denarr(0) = density for i=1,ngrids-1 do begin if (strupcase(subarr(i-1)) ne strupcase(subarr(i))) then $ process_xcom_gd,subarr(i),X,density,subs xlim = size(X) xx(0:xlim(1)-1,0:*,i) = X(0:*,0:*) denarr(i) = density endfor density = denarr text = strarr(ngrids) for i=0,ngrids-1 do text(i)=string(/print,i,subarr(i),density(i)) more,text goto,gweights endif else begin ; ; ; This block of the code creates grid parameters if no grid ; parameter file was input. ; IF (KEYWORD_SET(cgrid_sep)) THEN grid_sep = cgrid_sep else $ grid_sep = fcheck(grid_sep, 1750.) ;collimator length (mm) ;separation betw. top and bott. colls. ; ngrids = 12 ;number of grids (collimators) gridnm = indgen(ngrids) + 1 ;array of collimator numbers db = 71 ;set bottom coll. diameter (mm) dfov=.3 ;set angular overlap (degrees) ; diam of top coll. (mm) = ; diam. of bottom + 2 * coll. length * arctangent(angul overlap * pi/180) dt = db + 2.0 * grid_sep * atan(dfov * !pi/180) ;diameter of top coll. (mm) rim = 10 ;rim thickness around grids (mm?) ; testing a new slit width/pitch ratio on 1/24/95 ; .6 for 1-10, .55 for 11,12 frac = fltarr(12) & frac(0:9) = .6 & frac(10:11) = .55 ;frac = 0.6 ; original slit width/pitch ; test out new amin 1/24/95 amin = fltarr(12) & amin(0) = 38.66 amin(1:8) = 1.0 & amin(9:11) = [1.11,4.41,6.39] ;amin = 1. ;original minimum collimation angle in degrees rmin = 2.0 ;minimum angular resolution in arcseconds rmax = 120.0 ;maximum angular resolution in arcseconds ; log step size = (nat. log(max. ang. res.) - nat. log(min. ang. res.) / ; number of grids - 1 lrstep=(alog(rmax) - alog(rmin))/(ngrids-1) ;Logarithmic step size ; log resolution = nat. log(min. ang. res.) + [array 0,1...11] * step size lres = alog(rmin) + findgen(ngrids) * lrstep ;Logarithmic resolutions ; resolution in arcsecconds = exponential(log. resolution) res = exp(lres) ;Resolutions in arcseconds ; pitch = 2 * coll. lengt * tangent( resolution*pi/(3600.*180.)) pitch = 2.0 * grid_sep * tan(res*!pi/(3600.*180.)) ;Pitch ; slit width = pitch * slit width/pitch ratio sw = pitch * frac ;slit width ; calcu. thick = slit width/tangent(min. collimation angle*pi/180) cth = sw/tan(amin * !pi/180) ;calculated thickness to give 1 degree FOV if (grid_sep eq 5000.) then mth = 40 else mth = 10 ;maximum thickness of grids thickness = cth < mth ; allow only grids smaller than max. thickness ;field of view = arctangent(slit width/thickness) * 180/pi FOV = atan(sw/thickness) * 180/!pi ;FOV in degrees endelse gweights: ; Write an output file if keyword out is set ; Subs and Dens will follow that of XCOM Program! if (outfile ne '') then write_gd_data,outfile,grid_sep,db,$ dt,dfov,frac,subs,density,$ gridnm,res,pitch,sw,thickness,fov,fine_res,coarse_res,tfset,ngrids ; CALCULATE GRID WEIGHTS IN KG. ; wttop (top grids) = thickness*pi/4* ; ( (top diameter + rim thick)^2 - slit/width ratio * top d^2) ; * density/1e6 ; wtbot (bottom grids) = same formula, except top diam. = bottom diameter ; note: diameter refers to collimator diameter wttop = thickness * !pi/4. * ((dt + rim)^2. - frac * dt^2.) * density/1000000. wttoptot = total(wttop) ; sum of the 12 weights wtbot = thickness * !pi/4. * ((db + rim)^2. - frac * db^2.) * density/1000000. wtbottot = total(wtbot) ; sum of the 12 weights ; Calculate new Energy resulting from reduction in modulation amplitude ; and Energy at modulation peak (10%,50%, and 100% reductions) ; Look for the manual keyword, if the user has specified a reduction in ; modulation amplitude, then use that value in energy calucations ; Otherwise calculate red. in mod. amplitude effects at 10,50, &100% red. ; If keyword_set(rmod_amp) then begin man_rmod_amp = 1 goto, manual_rmodamp endif ; ; initialize arrays for 10,50, 100%, 0% reduction in mod. amplitude from ; the peak value of 1 for 12 collimators at 3 harmonics. ; gdiff = grid_diffraction([0.90,0.50,0.0,1.0],3,pitch,grid_sep, $ elim=0.1) E_diff_10= gdiff(*,*,0) E_diff_50= gdiff(*,*,1) E_diff_fl= gdiff(*,*,2) e_peak = gdiff(*,*,3) goto, skip_manual ; passes by the manual reduction (user defined) ;in modulation amplitude, rmodamp:reduction in mod. amp. ; ; User input value for reduction in modulation amplitude from diffraction. ; manual_rmodamp: c1 = 1.-(rmod_amp/100.) ; calculate diffraction limit (1-reduction) gdiff = grid_diffraction([c1, 1.0],3,pitch,grid_sep,elim=0.1) e_diff_manual = gdiff(*,*,0) e_peak = gdiff(*,*,1) skip_manual: ; Use Attenuation Coefficients Array (X) retrieved by subprogram ; process_xcom_gd to interpolate energies for GD.PRO ; User now has attenuation coefficient ; data in the array called X. The energies are in X(*,0) and the various ; attenuation coefficients are in X(*,1-7). ftran = 0.5 ;Fraction transmitted at highest energy covered. ; Calculate total attenuation coefficient that would give a transmision factor ; of ftran for the given thickness and density of the grids. ; newy = -natural_log(fraction transmitted at hightest en) / thickness/10. * density newy = -ALOG(ftran) / ( (thickness/10.) * density) newe = fltarr(ngrids) for i=0,ngrids-1 do begin newe(i) = exp(interpol(alog(xx(*,0,i)),alog(xx(*,6,i)),alog(newy(i)))) if (newe(i) eq 0) then newe(i) = 10. endfor ;print,'newy:',newy ;print,'newe',newe*1000. ; OPEN & WRITE GRID_DIMENSIONS.DATA With ALL PROGRAM INFORMATION ng = ngrids openw,lun,'grid_dimensions.data',/get_lun Printf,lun,'HEISPEC Grid Dimensions, Compiled: '+strmid(atime(sys2ut(0)),0,15) printf,lun,grid_sep,format= '(/"Collimator Length (mm) ",f10.0)' Printf,lun,db,format= '("Bottom collimator diameter (mm)", f10.1)' Printf,lun,dt,format= '("Top collimator diameter (mm) ", f10.1)' Printf,lun,dfov,format= '("Angular overlap (degrees) ",f10.3)' printf,lun,gridnm,format= '(/"Collimator ",13i7)' printf,lun,res,format= '("FWHM Resolution (arcsec)",13f7.2)' printf,lun,pitch,format='("Pitch (mm) ",13f7.3)' printf,lun,sw,format= '("Slit width (mm) ",13f7.3)' printf,lun,thickness,format= '("Thickness (mm) ",13f7.2)' printf,lun,fov,format= '("Field of View (deg.) ",13f7.2)' Printf,lun,rim,format= '(/"Grid rim thickness (mm) ",f7.2)' Printf,lun,density,format= '("Grid density (g/cm^3) ",f7.2)' printf,lun,wttop,format= '(/"Top grid weight (kg) ",13f7.2)' printf,lun,wttoptot,format= '("Total (kg) ",f7.2)' printf,lun,wtbot,format= '(/"Bottom grid weight (kg) ",13f7.2)' printf,lun,wtbottot,format= '("Total (kg) ",f7.2)' printf,lun,ftran*100.,format='(/"Energy for ", F5.1, "% transmission (10000 = > 10 MeV) ")' printf,lun,newe*1000.,format= '("Energy (keV) ",13i7/)' if (man_rmod_amp eq 1) then begin printf,lun,'Diffraction Energy limit (keV) with 10%, 50%, and 100% reduction in modulation amplitude' printf,lun,' (0.0 = < .1 keV) ' printf,lun,rmod_amp,e_diff_manual(*,0),format='(f5.1,"% red, @ fundamental ",13f7.1)' printf,lun,rmod_amp,e_diff_manual(*,0),format='(f5.1,"% red, @ 2nd harmon. ",13f7.1)' printf,lun,rmod_amp,e_diff_manual(*,0),format='(f5.1,"% red, @ 3rd harmon. ",13f7.1)' endif else begin printf,lun,'Diffraction Energy limit (keV) with 10%, 50%, and 100% reduction in modulation amplitude' printf,lun,' (0.0 = < .1 keV) ' printf,lun,e_diff_10(*,0),format='("10% red, @ fundamental ",13f7.1)' printf,lun,e_diff_10(*,1),format='("10% red. @ 2nd harmon. ",13f7.1)' printf,lun,e_diff_10(*,2),format='("10% red. @ 3rd harmon. ",13f7.1)' printf,lun,e_diff_50(*,0),format='("50% red. @ fundamental ",13f7.1)' printf,lun,e_diff_50(*,1),format='("50% red. @ 2nd harmon. ",13f7.1)' printf,lun,e_diff_50(*,2),format='("50% red. @ 3rd harmon. ",13f7.1)' printf,lun,e_diff_fl(*,0),format='("100% red. @ fundamental ",13f7.1)' printf,lun,e_diff_fl(*,1),format='("100% red. @ 2nd harmon. ",13f7.1)' printf,lun,e_diff_fl(*,2),format='("100% red. @ 3rd harmon. ",13f7.1)' endelse printf,lun,format='(/"First Diffraction Peak in Modulation Amplitude (keV)")' printf,lun,e_peak(*,0),format='("1st pk, fundamental ",13f7.1)' printf,lun,e_peak(*,1),format='("1st pk, 2nd harmonic ",13f7.1)' printf,lun,e_peak(*,2),format='("1st pk, 3rd harmonic ",13f7.1)' close,lun free_lun,lun ;spawn,'print grid_dimensions.data /setup=(r612,land)' goto,get_out read_error: Print,'Reading error in your input file, please see sys$user:[eric.gd]' Print,'File gd_data.dat for an example of the correct format' get_out: return end