function azam_a_image, a_file, kk $ , bkg=bkg, ref_frame=ref_frame, pix_deg=pix_deg ;+ ; ; function: azam_a_image ; ; purpose: return 2D image of ASP data file ; ; author: paul@ncar 4/95 ; ;============================================================================= ; ; Check number of parameters. ; if n_params() eq 0 then begin print print, "usage: array = azam_a_image( file [, kk] ) print print, " Return 2D array of ASP file. The data file must" print, " have the same format as imageable a_* files output" print, " by program 'bite' with -x option. print print, "arguments:" print, " file - string with path to input file." print, " kk - ASP structure of the type computed by" print, " 'azam_b_image.pro', 'azam_c_image.pro', 'e_image.pro'" print, " or 'azam_op.pro'. The structure is" print, " computed if necessary or reused if" print, " consistant." print, " (AZAM structures assumed already computed.)" print, "keywords:" print, " bkg - input background for unsolved points" print, " (def: 0.0)" print, " ref_frame - input string for reference frame." print, " 'instrument': returns data in instrument" print, " coordinates." print, " 'zenith': do neighbor interpolation to" print, " map data for observer at zenith to solar" print, " surface where the data was taken." print, " 'disk': do neighbor interpolation to" print, " map data for observer looking normal to" print, " disk center." print, " (Def: 'instrument' or use kk.rframe if" print, " kk structure useable as is.)" print, " (Note: only the first character of the" print, " string has to be correct.)" print, " pix_deg - great circle pixels per degree for scaling" print, " reference frames 'zenith' or 'disk';" print, " ignored for 'instrument' frame." print, " (def: 50. or use kk.pix_deg if" print, " kk structure useable as is)." print, " (Use about 45 for no magnification;" print, " > 70 may not map properly)." print return, 0 endif ;- ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;; 4/11/95 ;; ;; Function ~stokes/src/idl/a_image.pro is now available. The routine ;; selects and calls one of b_image.pro, c_image.pro, e_image.pro ;; or s_image.pro. ;; ;; It is hoped that a_image.pro is natural to use ('natural to use' ;; is not the same as 'easy to document'). ;; ;; With a_image.pro comes some upgraded terminology. 'Neighbor interpolated' ;; replaces 'stretched'. 'Zenith frame' replaces 'local frame'. ;; ;; Consider example sequence of statements: ;; ;; (1) b_cct = a_image( '/hilo/d/asp/data/red/92.06.19/o7/a__cct', b $ ;; (2) , ref_frame='instrument', bkg=0. ) ;; (3) b_fld = a_image( '/hilo/d/asp/data/red/92.06.19/o7/a_fld' , b $ ;; (4) , ref_frame='instrument', bkg=0. ) ;; (5) ;; (6) c_cct = a_image( '/hilo/d/asp/data/red/92.06.19/o7/a__cct', c $ ;; (7) , ref_frame='zenith', bkg=0., pix_deg=50. ) ;; (8) c_fld = a_image( '/hilo/d/asp/data/red/92.06.19/o7/a_fld' , c $ ;; (9) , ref_frame='zenith', bkg=0., pix_deg=50. ) ;; (10) ;; (11) e_cct = a_image( '/hilo/d/asp/data/red/92.06.19/o7/a__cct', e $ ;; (12) , ref_frame='disk', bkg=0., pix_deg=50. ) ;; (13) e_fld = a_image( '/hilo/d/asp/data/red/92.06.19/o7/a_fld' , e $ ;; (14) , ref_frame='disk', bkg=0., pix_deg=50. ) ;; ;; With some simplifications allowed by the code, lines (1 to 14) can be ;; replace by ( bkg=0. and pix_deg=50. are the defaults ): ;; ;; (15) dty = '/hilo/d/asp/data/red/92.06.19/o7/' ;; (16) ;; (17) b_cct = a_image( dty+'a__cct', b ) ;Default ref_frame='instrument' ;; (18) b_fld = a_image( dty+'a_fld' , b ) ;; (19) ;; (20) c_cct = a_image( dty+'a__cct', c , ref='z' ) ;; (21) c_fld = a_image( dty+'a_fld' , c ) ;; (19) ;; (20) e_cct = a_image( dty+'a__cct', e , ref='d' ) ;; (21) e_fld = a_image( dty+'a_fld' , e ) ;; ;; Equivalent code using b_image.pro, c_image.pro and e_image.pro would be: ;; ;; (17) dty = '/hilo/d/asp/data/red/92.06.19/o7/' ;; (18) ;; (19) b_cct = b_image( dty+'a__cct', b_str=b ) ;; (20) b_fld = b_image( dty+'a_fld' , b_str=b, /reuse ) ;; (21) ;; (22) c_cct = c_image( dty+'a__cct', c_str=c ) ;; (23) c_fld = c_image( dty+'a_fld' , c_str=c, /reuse ) ;; (21) ;; (22) e_cct = e_image( dty+'a__cct', e_str=e ) ;; (23) e_fld = e_image( dty+'a_fld' , e_str=e, /reuse ) ;; ;; ---------------------------------------------- ;; ;; Reference frame in a_image.pro is picked by keyword ref_frame. ;; Choices for reference frame are 'instrument', 'zenith' and 'disk'. ;; ;; ref_frame = 'instrument': ;; ;; Reference frame returned by b_image.pro. ;; It is a neighbor interpolation of the 'instrument' frame. ;; (x,y) coordinates correspond to the data as worked on by the Stokes ;; inversion code. ;; ;; ref_frame = 'zenith': ;; ;; Reference frame returned by c_image.pro. ;; It is a neighbor interpolation of 'instrument' frame. ;; (x,y) coordinates are for observer looking down the normal to solar ;; surface at the site where data was taken. ;; X increases to west along constant hiliocentric latitude. ;; Y is 90 degrees CCW from X. ;; ;; ref_frame = 'disk': ;; ;; Reference frame returned by e_image.pro. ;; It is a neighbor interpolation of 'instrument' frame. ;; (x,y) coordinates for an observer looking normal to solar disk center. ;; X increases to west and is parallel to the line drawn from observed ;; disk center to solar west. ;; Y is 90 degrees CCW from X and is parallel to central the meridian. ;; ;; --------------------------------------- ;; ;; Some things to note are (not all are in above examples): ;; ;; (A) a_image.pro returns the structure as a second (optional) argument. ;; The older programs used keywords. ;; ;; (B) a_image.pro does not have keyword 'reuse'. ;; ;; The structures has info about 'ref_frame' and 'pix_deg'. ;; The structures also has the directory path for the a_* file data set. ;; a_image.pro will reuse an input structure if keywords are not explicitly ;; reset. ;; ;; (C) Files corresponding to an a_* data set are sometimes written to other ;; directories. a_image.pro can read these files if it has input the ;; structure as input. ;; ;; (D) Reference frame will default to 'instrument' if not given by ;; keyword and the structure has to be computed or recomputed. ;; ;; (E) Only the first character for keyword 'ref_frame' is checked. ;; ref_frame='zenith', ref='z', ref='zzzz' are the same thing. ;; (Standard IDL allows abbreviation of keyword names to shortest ;; unambiguous case.) ;; ;; For more information see usage information in the source. ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; ;----------------------------------------------------------------------------- ;Let me explain the following a_image.pro zenith frame structure. ;(All examples have been tested form the text). ; ;----------------------------------------------------------------------------- ;cct = a_image( 'a__cct', c, r='z' ) ; ;a_image.pro: computing 'zenith' reference frame structure, ;file: a__cct ; ;help,cct ;CCT FLOAT = Array(78, 98) ; ;help, /str, c ;** Structure <27e000>, 28 tags, length=141820, refs=1: ;INDEX LONG 1 ;STRETCH LONG 1 ;NPOINTS LONG 3200 ;NSOLVED LONG 2801 ;NPXY LONG 7259 ;NSXY LONG 6358 ;X0 LONG 30 ;Y0 LONG 60 ;XPNT LONG Array(3200) ;YPNT LONG Array(3200) ;CCT_MIN FLOAT 12517.0 ;CCT_MAX FLOAT 39409.0 ;MXFLD LONG 4000 ;HEAD FLOAT Array(128) ;XDIM LONG 78 ;YDIM LONG 98 ;PXY LONG Array(7259) ;SXY LONG Array(6358) ;VEC_PXY LONG Array(7259) ;VEC_SXY LONG Array(6358) ;PBKG LONG Array(385) ;SBKG LONG Array(1286) ;PIX_DEG FLOAT 50.0000 ;MM_PER_DEG FLOAT 12.1480 ;DTY STRING '' ;RFRAME STRING 'zenith' ;FNAME STRING 'a__cct' ;FDIM LONG 3200 ; ;----------------------------------------------------------------------------- ;INDEX LONG 1 ; ; index = 0 for instrument frame structure. ; 1 for zenith frame structure. ; 2 for AZAM structure. ; 3 for AZAM ZOOM structure. ; 4 for disk frame structure. ; ; AZAM structures are not returned by a_image.pro but ; have most information. The 2 3 distiction is used internal ; to AZAM. ; ;----------------------------------------------------------------------------- ;STRETCH LONG 1 ; ; stretch = 0 if images are not neighbor interpolated (instrument frame). ; = 1 if images are neighbor interpolated (zenith or disk frame). ; ; Computer time can be saved if the code knows the images are not ; interpolated. Otherwise, instrument frame structures has one-to-one ; mapping arrays for trivial neighbor interpolating. ; ;----------------------------------------------------------------------------- ;NPOINTS LONG 3200 ; This is the total number of raster points output by inversion code. ; It is the length of a__* (double underscore) files. ; ;----------------------------------------------------------------------------- ;NSOLVED LONG 2801 ; This is the number of solved points output by inversion code. ; It is the length of a_* (single underscore) files. ; ;----------------------------------------------------------------------------- ;NPXY LONG 7259 ; This is the number of raster points after neighbor interpolating. ; In instrument frame npxy equals npoints. ; ;----------------------------------------------------------------------------- ;NSXY LONG 6358 ; This is the number of solved points after neighbor interpolating. ; In instrument frame nsxy equals nsolved. ; ;----------------------------------------------------------------------------- ;X0 LONG 30 ; Smallest scan step in the data set. ; ;----------------------------------------------------------------------------- ;Y0 LONG 60 ; Smallest pixel position along the slit in the data set. ; ;----------------------------------------------------------------------------- ;XPNT LONG Array(3200) ; Scan step numbers corresponding to raster points output ; by the inversion code. ; ;----------------------------------------------------------------------------- ;YPNT LONG Array(3200) ; Pixel along the slit corresponding to raster points output ; by the inversion code. ; ;----------------------------------------------------------------------------- ;CCT_MIN FLOAT 12517.0 ;CCT_MAX FLOAT 39409.0 ; Minimum & maximum value of the continuum (from a__cct file). ; ;----------------------------------------------------------------------------- ;MXFLD LONG 4000 ; MXFLD is included only to be compatable with AZAM structures; ; gives color scale maximum for ploting field strength. ; ;----------------------------------------------------------------------------- ;HEAD FLOAT Array(128) ; Contents of a___header file. ; Array elements are: ; ; ( 0) month ; ( 1) day ; ( 2) year ; ( 3) inversion code level ; ( 4) number of pixels along the slit (usually extended to 512). ; ( 5) number of multiplet lines. ; ( 6) First pixel in B0 average (rarely used inversion option) ; ( 7) Last pixel in B0 average (rarely used inversion option) ; ( 8) 1 if psi angle is fit. ; ( 9) 1 if field is fit. ; ( 10) 1 if Doppler width is fit. ; ( 11) 1 if 1st line of multiplet is fit. ; ( 12) 1 if damping is fit. ; ( 13) 1 if eta0 is fit. ; ( 14) 1 if azimuth is fit. ; ( 15) 1 if B1*mu is fit. ; ( 16) 1 if B0 is fit. ; ; ( 21) 1 if 2nd line of multiplet is fit. ; ; ( 23) 1 if eta0 if fit for 2nd multiplet line independent ; of eta0 for 1st line. Rarely used inversion option; ; strength ratio is usually used. ; ( 24) 1 if scattered light shift is fit. ; ( 25) 1 if scattered light amplitude is fit. ; ( 26) 0 if magneto optics is turned off (rarely used inversion option). ; ; ( 32) least squares weight used for I ; ( 33) least squares weight used for Q ; ( 34) least squares weight used for U ; ( 35) least squares weight used for V ; ; ( 36) initial guess center 1st line (FORTRAN wavelength pixel) ; ( 37) initial guess center 2nd line (FORTRAN wavelength pixel) ; ; ( 39) spectra dispersion, mA. ; ( 40) wavelength 1st spectra line, A. ; ( 41) wavelength 2nd spectra line, A. ; ; ( 43-44) band pass about 1st spectra line, mA. ; ( 45-46) band pass about 2nd spectra line, mA. ; ; ( 49) number of smearing function or laser profile used. ; ( 50) width of theoretical smearing function (rarely used). ; ; ( 57) number of points along the slit (usually extended to 512). ; This is used to assign raster point numbers. ; ; ( 96) = ephemeris, right ascension of disk center, degrees. ; ( 97) = ephemeris, declination of disk center, degrees. ; ( 98) = ephemeris, right ascension from central meridian, degrees. ; ( 99) = ephemeris, b0-angle, degrees. ; (100) = ephemeris, p-angle, degrees. ; (101) = ephemeris, carrington longitude disk center, degrees. ; (102) = ephemeris, solar radius, arcsecs. ; (103) = azimuth of disk center from data set center (zenith ; frame), degrees ccw from west. ; (104) = parallactic angle, degrees. ; ;----------------------------------------------------------------------------- ;XDIM LONG 78 ;YDIM LONG 98 ; x & y dimension of resulting 2D arrays. ; ;----------------------------------------------------------------------------- ;VEC_PXY LONG Array(7259) ;PXY LONG Array(7259) ; ; Where arrays used to map data from a__* (double underscore) ; files to 2D arrays. ; ; VEC_PXY is where array used to expand data output by the inversion ; code to the number of data points in the image arrays. ; In instrument frame VEC_PXY is a trivial one-to-one mapping. ; ; PXY is where array to map the above expansion into 2D arrays. ; ; ;example to form 2D array of a__pip file: ; ; ; ; read data file ; a__pip_data = read_floats( 'a__pip' ) ; ; ; ; initialize 2D array ; pip = fltarr(c.xdim,c.ydim) ; ; ; ; expand and move data into 2D array. ; ; In instrument frame the (c.vec_pxy) is not required. ; pip(c.pxy) = a__pip_data(c.vec_pxy) ; ; ;example to do the reverse of above: ; ; ; ; initialize vector array to hold data. ; a__pip_data = fltarr(c.npoints) ; ; ; ; move data from 2D image and compress. ; ; Note: it is assumed c.npoints is less than or equal c.npxy ; a__pip_data(c.vec_pxy) = pip(c.pxy) ; ; ; ; write the data file. ; error = write_floats( 'a__pip', a__pip_data ) ; ;----------------------------------------------------------------------------- ;VEC_SXY LONG Array(6358) ;SXY LONG Array(6358) ; ; Where arrays used to map data from a_* (single underscore) ; files to 2D arrays. ; ; VEC_SXY is where array used to expand data output by the inversion ; code to the number of data points in the image arrays. ; In instrument frame VEC_SXY is a trivial one-to-one mapping. ; ; SXY is where array to map the above expansion into 2D arrays. ; ; ;example to form 2D array of a_fld file: ; ; ; ; read data file ; a_fld_data = read_floats( 'a_fld' ) ; ; ; ; initialize 2D array ; fld = fltarr(c.xdim,c.ydim) ; ; ; ; expand and move data into 2D array. ; ; In instrument frame the (c.vec_sxy) is not required. ; fld(c.sxy) = a_fld_data(c.vec_sxy) ; ; ;example to do the reverse of above: ; ; ; ; initialize vector array to hold data. ; a_fld_data = fltarr(c.nsolved) ; ; ; ; move data from 2D image and compress. ; ; Note: it is assumed c.nsolved is less than or equal c.nsxy ; a_fld_data(c.vec_sxy) = fld(c.sxy) ; ; ; ; write the data file. ; error = write_floats( 'a_fld', a_fld_data ) ; ;----------------------------------------------------------------------------- ;PBKG LONG Array(385) ; Where array for the background locations in 2D images formed from ; a__* (double underscore) files. A null where array is -1. ; ; ;example to set background locations to .5 in pip array. ; if (c.pbkg)(0) ge 0 then pip(c.pbkg)=.5 else print, 'no background' ; ; Note: (c.pbkg)(0) is the first element of c.pbkg; ; the parentheses are required. ;----------------------------------------------------------------------------- ;SBKG LONG Array(1286) ; Where array for the background locations in 2D images formed from ; a_* (single underscore) files. A null where array is -1. ; ; ;example to set background locations to 500 in fld array. ; if (c.sbkg)(0) ge 0 then fld(c.sbkg)=500. else print, 'no background' ; ; Note: (c.sbkg)(0) is the first element of c.sbkg; ; the parentheses are required. ; ;----------------------------------------------------------------------------- ;PIX_DEG FLOAT 50.0000 ; Number of image pixels per degree of great circle on the sun. ; Does not apply in instrument frame. 50 is the default for ; for zenith and disk reference frames. This can be changed ; on keyword to 'a_image.pro'. ; ;----------------------------------------------------------------------------- ;MM_PER_DEG FLOAT 12.1480 ; Hardwired mega meters per degree of great circle on the sun. ; ;----------------------------------------------------------------------------- ;RFRAME STRING 'zenith' ; Reference frame as a string value. Other values could be ; 'instrument' or 'disk'. ; ;----------------------------------------------------------------------------- ;DTY STRING '' ; Path to the directory with the a_* files. In this example ; the path is the current working directory. ; ;----------------------------------------------------------------------------- ;FNAME STRING 'a__cct' ; Name of last file input (not often used variable). ; ; ; print full path from current working directory to last file read. ; print, c.dty+c.fname ; ;----------------------------------------------------------------------------- ;FDIM LONG 3200 ; Length of the last file input (not often used variable). ; FDIM should equal either c.npoints or c.nsolved. ; ;----------------------------------------------------------------------------- ;Set scale or default. if n_elements(pix_deg) ne 0 then begin pix_scale = pix_deg end else begin if n_elements(kk) ne 0 $ then pix_scale = kk.pix_deg $ else pix_scale = 50. end ;Set reference frame or default. if n_elements(ref_frame) ne 0 then begin rframe = ref_frame end else begin if n_elements(kk) ne 0 $ then rframe = kk.rframe $ else rframe = 'instrument' end ;Check first character for valid ;reference frame. if strmid(rframe,0,1) eq 'i' then rframe='instrument' if strmid(rframe,0,1) eq 'z' then rframe='zenith' if strmid(rframe,0,1) eq 'd' then rframe='disk' if rframe ne 'instrument' then $ if rframe ne 'zenith' then $ if rframe ne 'disk' then begin stop, rframe+': unknown reference frame' return, 0 end ;If structure isn't defined or reference ;frame changed, compute 2D array and ;structure simultaneously. if n_elements(kk) eq 0 then goto, recompute0 if rframe ne kk.rframe then goto, recompute0 ;Check if image scale has changed, ;compute 2D array and structure ;simultaneously. if (rframe eq 'zenith' or rframe eq 'disk') $ and pix_scale ne kk.pix_deg $ then goto, recompute0 ;Get file length (in status). on_ioerror, break001 openr, /get_lun, unit, a_file status = fstat( unit ) sz = status.size/4 free_lun, unit goto, break002 break001: print, 'azam_a_image.pro: '+!err_string stop return, -1 break002: ;Reuse structure if file length ok. if sz ne kk.npoints and sz ne kk.nsolved then begin ;Get directory path by finding last '/'. j = -1 i = 0 while i ne -1 do begin i = strpos( a_file, '/', i ) if i ne -1 then begin j = i i = i+1 end end dty = strmid( a_file, 0, j+1 ) ;Check if file is in another directory. if dty ne kk.dty then begin ;Key on a___header file to check if ;file is in another a_* directory. ;Try to open a___header file. on_ioerror, break000 openr, /get_lun, unit, dty+'a___header' free_lun, unit ;File is in another a_* directory. goto, recompute0 end ;Inform user of error. break000: print, 'azam_a_image.pro: file is wrong length,' stop, 'file: '+a_file return, -1 end ;Return 2D array using current already ;defined structure. return, s_image( a_file, kk, bkg=bkg ) ; ------------------------------------------- ;Have to recompute structure. recompute0: ;Set scale and reference frame. if n_elements(pix_deg) ne 0 then pix_scale=pix_deg else pix_scale=50. if n_elements(ref_frame) ne 0 then rframe=ref_frame else rframe='instrument' if strmid(rframe,0,1) eq 'i' then rframe='instrument' if strmid(rframe,0,1) eq 'z' then rframe='zenith' if strmid(rframe,0,1) eq 'd' then rframe='disk' print print, "azam_a_image.pro: computing '"+rframe+"' reference frame structure," print, "file: "+a_file print ;Simultaneous compute and return 2D array ;and structure. Return 2D array. case rframe of 'instrument': return, azam_b_image( a_file, b_str=kk, bkg=bkg ) 'zenith': return, azam_c_image( a_file, c_str=kk, bkg=bkg, pix_deg=pix_scale ) 'disk': return, azam_e_image( a_file, e_str=kk, bkg=bkg, pix_deg=pix_scale ) end end