;------------------------------------------------------------------------------ ;+ ; NAME: AUX_READ_GALAXEV ; ; RESTRICTIONS: This routine is not approved part of the ULySS ; package. It is not tested and was designed for a ; particular purpose. ; ; PURPOSE: ; Read a grid of GalaxEV models ; ; USAGE: grid = aux_read_galaxev( model_dir, $ ; WAVERANGE=waverange, VELSCALE=velscale, $ ; AGERANGE=agerange, METALRANGE=metalrange, $ ; MGFE=mgfe, QUIET=quiet ) ; ; ARGUMENTS: ; model_dir: string, directory where the SSP .txt files are located ; To read models with 2 Mg/Fe enhancements, give ; a string array with the two directories. ; KEYWORDS: ; WAVERANGE: Two element array giving the limits of the range to read, ; in Angstrom. For example WAVERANGE=[4000d,4786d]. ; By default the whole data range is read. ; ; VELSCALE: Size of the pixels in km/s. For example VELSCALE=30. ; If VELSCALE is not given, the models are not rebinned, ie. ; they are as stored in the disk file (in linear or ; log scale). ; ; AGERANGE/METALRANGE: Extracted region of the grid. Ages in ; log(Myr) and metallicities in dex. ; By default the whole data range is read. ; For example ; AGERANGE=alog([100,15000]),METALRANGE=[-1.4,+0.7] ; ; MGFE: Mg/Fe ratio. A value to be added in the header ; ; QUIET: Verbosity control ; ; RETURN: A structure containing the readed models with the ; associeted WCS (see also ULY_SSP_ALLOC) ; ; DESCRIPTION: This routine can be used to generate a new grid of models ; by reading the txt files produced by the program GalaxEv/galaxpl. ; ; It is useful when, for example, a new library, or a ; new version becomes available. ; ; The grid read by AUX_READ_GALAXEV can be saved to a file ; and afterward used to define a component with ULY_SSP. ; The set of GalaxEv grids of models provided with ; ULySS in the directory 'models/' were generated by ; this program. ; ; EXAMPLE: To read the latest version of GalaxEv/Miles and save it ; grid = aux_read_galaxev(uly_root+'/models/galaxev/bc2003_hrs_miles_0p1_ssp') ; uly_ssp_write, grid, uly_root+'/models/test.fits' ; ; TESTS: ; test_models, uly_root+'/models/test.fits' ; sp = uly_spect_read(uly_root+'/data/VazMiles_z-0.40t07.94.fits') ; ulyss, sp, uly_ssp(model=uly_root+'/models/test.fits'), /plot ; ; HISTORY: Mina Koelva 2011/11/11, modified version of aux_read_pegase.pro ;- ; CATEGORY: AUX ;------------------------------------------------------------------------------ function read_ssp_glxv, filename, $ WAVERANGE=waverange, VELSCALE=velscale, $ AGERANGE=agerange, METALRANGE=metalrange, $ ORIGSSP=rr0, CONVSSP=rr, $ OUTAGES=outages, $ ; log(age) QUIET=quiet, $ STATUS=status, MODEL_DIR=model_dir ; READ_SSP_GLXV (internal routine) ; one-by-one processing of each SSP file compile_opt idl2, hidden on_error, 2 status = 0 ; initialize the error status ;; if filename is not set, or if the pointed file does not exist display a ;; message and return with an error if (filename eq "") then begin message, "no file to read (1st arg is empty)" endif if (file_test(filename) eq 0) then begin message, "filename is not valid (1st arg="+filename+")" endif template_grid_glxv = uly_ssp_alloc() time0 = systime(1) ;; read the entire lines of the first nage = 221 ; number of age points (see the first column of the file) nwl = 8195 ; number of wl points (see the first column of the file) print message, /INFO, '!!! WE ASSUME THAT THE WL POINTS AND AGE POINTS ARE FIXED TO 8195 AND 221 !!!' print outages = replicate(0., nage+1) ; first column only shows the number of entries wavelength = replicate(!Values.F_NaN, nwl+1) trash = ' ' fluxes = replicate(!Values.F_NaN, nwl+1) rr = fltarr(nwl,nage) openr, lun, filename, /GET_LUN readf, lun, outages ; 1st readf, lun, trash readf, lun, trash readf, lun, trash readf, lun, trash readf, lun, trash ; rejection test upon ages range outages = outages[1:*] outages = alog(float(outages) / 1.e6) ;; ages if max(alog(agerange)) ge min(outages) and min(alog(agerange)) le max(outages) then begin ageind = where(outages ge alog(agerange[0]) and outages le alog(agerange[1])) outages = outages[ageind] endif else begin message, 'intersection between agerange and template too small: STOP' endelse ;; reading fluxes and wavelenghts readf, lun, wavelength wavelength = wavelength[1:*] nlines = 0 while nlines lt nage do begin readf, lun, fluxes rr[*,nlines] = fluxes[1:*] nlines ++ endwhile rr = rr[*,ageind] ; rejection test upon wavelength range good = where(wavelength ge waverange[0] and wavelength le waverange[1], cnt) if cnt eq 0 then message, $ 'No intersection between waverange and models: STOP' nummin = min(good, MAX=nummax) wavelength = wavelength[good] sampling = 2 slashpos = strpos(filename, '/', /REVERSE_SEARCH) filename0 = strsplit(strmid(filename, slashpos+1), '_', /EXTR) met = {names:['m22','m32','m42','m52','m62','m72'], values:['-2.30103','-1.69897','-0.698970','-0.397940','0.00000','0.397940']} nn = 0 mcnt = 0 while nn lt n_elements(met.names) and mcnt eq 0 do begin mets = where(met.names[nn] eq filename0, mcnt) nn ++ end if mcnt gt 0 then outmetal=met.values[nn-1] else begin print, filename message, "Metallicity is not of the form 'm22','m32','m42','m52','m62','m72'... STOP" endelse ; print, 'reading time1=', systime(1)-time0 ; time 1 is small time0 = systime(1) rr = rr[good,*] * 1e+07 bunit = '1e-07 LSUN/MSUN/0.1nm' ; The rebinning it taking 60% of the total model reading time GLXV_spectrum = uly_spect_alloc(DATA=rr, SAM=sampling, WAVELEN=wavelength) GLXV_log_spectrum = uly_spect_logrebin(GLXV_spectrum, velscale) *template_grid_glxv.data = *GLXV_log_spectrum.data template_grid_glxv.start = GLXV_log_spectrum.start template_grid_glxv.step = GLXV_log_spectrum.step template_grid_glxv.sampling = GLXV_log_spectrum.sampling uly_spect_free,GLXV_spectrum uly_spect_free,GLXV_log_spectrum ;print, 'reading time2=', systime(1)-time0 ;time0 = systime(1) *template_grid_glxv.o_age = outages *template_grid_glxv.o_metal = outmetal mkhdr, hhh, *template_grid_glxv.data, /EXTEND *template_grid_glxv.hdr = hhh return, template_grid_glxv end ;------------------------------------------------------------------------------ function aux_read_galaxev, model_dir, $ WAVERANGE=waverange, VELSCALE=velscale, $ AGERANGE=agerange, METALRANGE=metalrange,$ MGFE=mgfe, QUIET=quiet compile_opt idl2 on_error, 0 if not keyword_set(waverange) then waverange = [3000, 9000] ; 0.1nm if not keyword_set(agerange) then agerange = [1, 23000] ; Myr if not keyword_set(metalrange) then metalrange = [-5.0, +3.0] ; dex template_grid = uly_ssp_alloc() ; detect the fits files in the given directory ; print, model_dir[0]+'/'+'*_ASCII' filename = file_search(model_dir[0]+'/'+'*_ASCII', COUNT=nmetall) if (n_elements(filename) eq 1) then begin if (filename eq "") then print, "model not found ("+model_dir[0]+")" return endif ; make an array for the output metallicities ; outmetall = make_array(nmetall,/DOUBLE,/NOZERO) template_grid_glxv = read_ssp_glxv(filename[0], $ WAVERANGE=waverange, $ AGERANGE=agerange, $ METALRANGE=metalrange, $ VELSCALE=velscale, $ ORIGSSP=origssp, $ QUIET=quiet, $ MODEL_DIR=model_dir $ ) ; grasp general data from the 1st readed fits file sz1 = (size(*template_grid_glxv.data))[1] outages = *template_grid_glxv.o_age nages = n_elements(outages) if template_grid_glxv.step ne 0 then begin template_grid.start = template_grid_glxv.start template_grid.step = template_grid_glxv.step template_grid.sampling = template_grid_glxv.sampling *template_grid.hdr = *template_grid_glxv.hdr endif outmetall[0] = *template_grid_glxv.o_metal if(nmetall gt 1) then begin templates = dblarr(sz1,nages,nmetall) templates[*,*,0] = (*template_grid_glxv.data)[*,*] uly_ssp_free, template_grid_glxv for i = 1, nmetall-1 do begin template_grid_glxv = read_ssp_glxv(filename[i], $ WAVERANGE=waverange, $ AGERANGE=agerange, $ METALRANGE=metalrange, $ VELSCALE=velscale, $ ORIGSSP=origssp, $ QUIET=quiet, $ MODEL_DIR=model_dir $ ) if template_grid_glxv.step ne 0 then begin template_grid.start = template_grid_glxv.start template_grid.step = template_grid_glxv.step template_grid.sampling = template_grid_glxv.sampling endif outmetall[i] = *template_grid_glxv.o_metal templates[*,*,i] = (*template_grid_glxv.data)[*,*] uly_ssp_free, template_grid_glxv endfor endif *template_grid.data = templates *template_grid.o_age = outages *template_grid.o_metal = outmetall templates = 0 outages = 0 uly_ssp_free, template_grid_glxv if n_elements(model_dir) eq 2 then begin ; This is a grid with [Mg/Fe] templates = dblarr(sz1,nages,nmetall) filename2 = file_search(model_dir[1]+'/' + '*.fits') nmetall2 = n_elements(filename) ; nb of files to be readed if nmetall2 ne nmetall then $ message, 'The enhanced grid must have the same N of metallicities' for i = 0, nmetall-1 do begin template_grid_glxv = read_ssp_phr(filename2[i], $ WAVERANGE=waverange, $ AGERANGE=agerange, $ METALRANGE=metalrange, $ VELSCALE=velscale, $ ORIGSSP=origssp, $ QUIET=quiet, $ MODEL_DIR=model_dir $ ) outmetall[i] = *template_grid_glxv.o_metal outages = *template_grid_glxv.o_age templates[*,*,i] = (*template_grid_glxv.data)[*,*] uly_ssp_free,template_grid_glxv endfor a = where(outages - *template_grid.o_age,n) if n ne 0 then $ message,'ages from both grids must be the same' b = where(outages - *template_grid.o_age,m) if m ne 0 then $ message,'ages from both grids must be the same' *template_grid.o_mgfe = [0.0, 0.4] ; always these 2 enhancement levels form = [size(templates,/DIM), 2] tot = n_elements(templates) *template_grid.data = reform([reform(*template_grid.data, tot), reform(templates,tot)], form) endif template_grid.title = model_dir[0] *template_grid.waverange = waverange ; compute the 2nd derivative in case we want to pass this grid to ulyss if size(*template_grid.data, /N_DIM) eq 3 then begin ; normal grid *template_grid.d2t = derive_3d(TMPARR=*template_grid.data, GRID=*template_grid.o_age) endif else begin ; grid with Mg/Fe resolution form = size(*template_grid.data,/DIM) tot = n_elements((*template_grid.data)[*,*,*,0]) *template_grid.d2t = [reform(derive_3d(TMPARR=(*template_grid.data)[*,*,*,0], GRID=*template_grid.o_age), tot),$ reform(derive_3d(TMPARR=(*template_grid.data)[*,*,*,1], GRID=*template_grid.o_age),tot)] *template_grid.d2t = reform(*template_grid.d2t, form) endelse ;; fix the header sxaddpar, *template_grid.hdr, 'COMMENT', 'Grid distributed with ULySS http://ulyss.univ-lyon1.fr/' if n_elements(mgfe) then sxaddpar, *template_grid.hdr, 'MGFE', mgfe, 'Mg/Fe ratio' sxaddpar, *template_grid.hdr, 'COMMENT', 'Metal min='+string(min(*template_grid.o_metal),F='(f5.2)')$ +' max='+string(max(*template_grid.o_metal),F='(f5.2)')+' dex' sxaddpar, *template_grid.hdr, 'COMMENT', 'Age min='+strtrim(min(exp(*template_grid.o_age)),2)+$ ' max='+strtrim(max(exp(*template_grid.o_age)),2)+' Myr' sxaddpar, *template_grid.hdr, 'COMMENT', 'Evolutionary tracks Padova 1994' sxaddpar, *template_grid.hdr, 'ULY_TYPE', 'SSP' return, template_grid end ;-- end -----------------------------------------------------------------------