; ------------------------------------------------------------------------------------------------------- ; This code plots the spatial distribution plots of TEMPO derived true and derived tropospheric ; columns and WRF VMR and emissions of NO2 ; ; User can provide date (Aug 11, 12), averaging time ; ;- Procedures Used - ; read_TEMPO_hourly.pro ; scatter_tempo_wrf_density.pro ; This plots the scatter plot ; spl_trop_col_wrf.pro ; This plots the spatial distribution ; con_lc.pro ; -------------------------------------------------------------------------------------------------------- @read_TEMPO_hourly.pro @spl_trop_col_wrf.pro @scatter_tempo_wrf_density.pro !quiet=1 ; ----------------------------------------------------------------------------------- ; Print,'1: CODE will generate 5 maps for 5 WRF simulations ' ; --------------------------- Paths and relevant Info --------------------------- species = 'NO2' spc_ttl = 'NO!I2!N' tmp_path = '/home/pfister/TEMPO/FRAPPE_TEMPO/' wrf_path = '/home/pfister/TEMPO/FRAPPE_WRFchem/TEMPO/' plot_path= '/home/pfister/TEMPO/Piyush/plots/' wrf_stup = ['WRFV4_S05_d02_MOZAIC_final_adj_EPA2014_newWRF_2xOG/', $ 'WRFV4_S05_d02_MOZAIC_final_adj_EPA2014_newWRF_phys/', $ 'WRFV4_S05_d02_MOZAIC_final_adj_EPA2014_newWRF_RAQMS/',$ 'WRFV4_S05_d02_MOZAIC_final_NEI2011adj_newWRF/', $ 'WRFV4_S05_d02_MOZAIC_final_S05_newWRF/' ] tem_prfx = 'l2_synth/TEMPO_201408' wrf_prfx = 'wrfout_d02_2014-08-' scenario = ['adj_EPA2014_newWRF_2xOG','adj_EPA2014_newWRF_phys', $ 'adj_EPA2014_newWRF_RAQMS','NEI2011adj_newWRF','S05_newWRF'] ; ----------------------------------------------------------------------------------- ;goto, skip1 ; Specify date of interest; Read input from the terminal: dd = '' READ, dd, PROMPT='Enter one date of interest (month is Aug) e.g., 11 or 12: ' ; Specify single or multiple hour input; Read input from the terminal: avg = '' READ, avg, PROMPT='Enter hourly/averaged output e.g., hour or avg: ' if avg eq 'hour' then begin hh = '' READ, hh, PROMPT='Enter one hour between [12:23] e.g.,13 : ' endif if avg eq 'avg' then begin hh1 = ' ' hh2 = ' ' READ, hh1, PROMPT='Enter hourly range between [12:23] (e.g.,13,15 for 3 hours avg.) Enter start hour: ' READ, hh2, PROMPT='Enter hourly range between [12:23] (e.g.,13,15 for 3 hours avg.) Enter end hour: ' endif skip1: ; ----------------------------------------------------------------------------------- ; Used the section below for testing -- ;dd='11' ;avg='hour' ;hh=13 ;for hh1 = 12, 23 do begin ; use this loop for bulk hourly plots ;hh = strcompress(string(hh1),/remove_all) ; use this loop for bulk hourly plots ;avg='avg' ;hh1=13 ;hh2=16 ; ----------------------------------------------------------------------------------- for n = 0, n_elements(wrf_stup)-1 do begin ;for n = 0, 0 do begin ; Finding number of TEMPO files for selected day and hour(s) file_nm = strcompress(tmp_path+wrf_stup[n]+tem_prfx+dd+'_*h_CO_'+species $ +'.nc',/remove_all) files1 = file_search(file_nm, count=numfiles1) avl_hr1 = strmid(files1, strlen(tmp_path+wrf_stup[n]+tem_prfx+dd+'_'), 2) ; Finding number of WRF output files for selected day and hour(s) file_n1 = strcompress(wrf_path+wrf_stup[n]+wrf_prfx+dd+'_*', /remove_all) files2 = file_search(file_n1, count=numfiles2) avl_hr2 = strmid(files2, strlen(wrf_path+wrf_stup[n]+wrf_prfx+dd+'_'), 2) ; ======================================================================================== ; Part-I : If single hour is chosen, this section of the code works ; ---------------------------------------------------------------------------------------- if avg eq 'hour' then begin ;- Reading TEMPO data file_ind = where(avl_hr1 eq hh) if file_ind[0] ne -1 then begin filename = files1[file_ind] ;print, filename read_TEMPO_hourly, filename, output lat_all = output.lat lon_all = output.lon trop_column_all = output.trop_column true_trop_column_all = output.true_trop_column trop_column_all[where(trop_column_all lt -1.0e29)] = !Values.F_NAN true_trop_column_all[where(true_trop_column_all lt -1.0e29)] = !Values.F_NAN endif ;- Reading WRF-Chem output file_ind = where(avl_hr2 eq hh) if file_ind[0] ne -1 then begin filename = files2[file_ind] ;print, filename id = ncdf_open(filename, /nowrite) ncdf_varget, id, 'E_NO2', emi_no2 ncdf_varget, id, 'no2', vmr_no2 ncdf_varget, id, 'XLONG', wrf_lon ncdf_varget, id, 'XLAT', wrf_lat ncdf_close, id wrf_no2_vmr = reform(vmr_no2[*,*,0])*1e3 wrf_no2_emi = reform(emi_no2[*,*,0]) ;help, wrf_no2_vmr, wrf_no2_emi endif ; ----------------------------------------------------------------------------------------- ;- ploting section outfl = strcompress(plot_path+'WRF_TEMPO_trop_colum_'+species+'_'+scenario[n]+ $ '_'+dd+'aug_'+string(hh)+'UTC.eps' ,/remove_all) spl_trop_col_wrf, outfl, trop_column_all, true_trop_column_all, lat_all, lon_all, $ wrf_no2_vmr, wrf_no2_emi, wrf_lat, wrf_lon, dd, string(hh), scenario[n], spc_ttl ; ----------------------------------------------------------------------------------------- outfile = strcompress(plot_path+'WRF_TEMPO_scatter_'+species+'_'+scenario[n]+ $ '_'+dd+'aug_'+string(hh)+'UTC.eps' ,/remove_all) scatter_tempo_wrf_density, outfile, trop_column_all, true_trop_column_all, lat_all, lon_all, $ wrf_no2_vmr, wrf_no2_emi, wrf_lat, wrf_lon, dd, string(hh), scenario[n], spc_ttl ; ----------------------------------------------------------------------------------------- endif ; ======================================================================================== ; Part-II : If hourly average is chosed this section of the code works ; ---------------------------------------------------------------------------------------- if avg eq 'avg' then begin nfiles = long(hh2)-long(hh1)+1 ;- Reading TEMPO output ;- Declare arrays that would be used [column is [mirrorstep, xtrack]] mstep = 158 xtrack = 279 trop_column_all = fltarr(xtrack,mstep, nfiles) true_trop_column_all = fltarr(xtrack,mstep, nfiles) lon_all = fltarr(xtrack,mstep, nfiles) lat_all = fltarr(xtrack,mstep, nfiles) str_in = uint(where(avl_hr1 eq hh1)) end_in = uint(where(avl_hr1 eq hh2)) for nf = str_in[0], end_in[0] do begin ;print, nf,' ',files1[nf],' ',nf-str_in[0],'/',nfiles read_TEMPO_hourly, files1[nf], output lat_all[*,*,nf-str_in[0]] = output.lat lon_all[*,*,nf-str_in[0]] = output.lon trop_column_all[*,*,nf-str_in[0]] = output.trop_column true_trop_column_all[*,*,nf-str_in[0]] = output.true_trop_column endfor trop_column_all[where(trop_column_all lt -1.0e29)] = !Values.F_NAN true_trop_column_all[where(true_trop_column_all lt -1.0e29)] = !Values.F_NAN ; calculate and hourly averaged values lat_avg = fltarr(xtrack, mstep) lon_avg = fltarr(xtrack, mstep) trop_col_avg = fltarr(xtrack, mstep) true_trop_col_avg = fltarr(xtrack, mstep) for ii = 0, xtrack-1 do begin for jj = 0, mstep-1 do begin ; print, ii, jj lat_avg[ii,jj] = mean(reform(lat_all[ii,jj,*])) lon_avg[ii,jj] = mean(reform(lon_all[ii,jj,*])) xx = reform(trop_column_all[ii,jj,*]) xx1 = where(finite(xx) eq 1) if xx1[0] ne -1 then trop_col_avg[ii,jj] = mean(reform(xx),/NAN) yy = reform(true_trop_column_all[ii,jj,*]) yy1 = where(finite(yy) eq 1) if yy1[0] ne -1 then true_trop_col_avg[ii,jj] = mean(reform(yy),/NAN) ;trop_col_avg[ii,jj] = mean(reform(trop_column_all[ii,jj,*]),/NAN) ;true_trop_col_avg[ii,jj] = mean(reform(true_trop_column_all[ii,jj,*]),/NAN) endfor endfor ; ----------------------------------------------- ;- Reading WRF-Chem files ;- getting input (lat/lon) dims and array declaration id = ncdf_open(files2[0], /nowrite) ncdf_varget, id, 'XLONG', xlong ncdf_close, id dims = size(xlong) wlon = dims[1] wlat = dims[2] wrf_no2_vmr = fltarr(wlon, wlat, nfiles) wrf_no2_emi = fltarr(wlon, wlat, nfiles) str_in = uint(where(avl_hr2 eq hh1)) end_in = uint(where(avl_hr2 eq hh2)) ;print,'4: ',str_in, end_in, nfiles for nf = str_in[0], end_in[0] do begin ;print, nf,' ',files2[nf],' ',nf-str_in[0],'/',nfiles id = ncdf_open(files2[nf], /nowrite) ncdf_varget, id, 'E_NO2', emi_no2 ncdf_varget, id, 'no2', vmr_no2 ncdf_varget, id, 'XLONG', wrf_lon ncdf_varget, id, 'XLAT', wrf_lat ncdf_close, id wrf_no2_vmr[*, *, nf-str_in[0]] = reform(vmr_no2[*,*,0])*1e3 wrf_no2_emi[*, *, nf-str_in[0]] = reform(emi_no2[*,*,0]) ;print, minmax(reform(wrf_no2_vmr[*, *, nf-str_in[0]]),/NAN) endfor ;- Averaging WRF output wno2_emi_avg = fltarr(wlon, wlat) wno2_vmr_avg = fltarr(wlon, wlat) for ii = 0, wlon-1 do begin for jj = 0, wlat-1 do begin wno2_emi_avg[ii,jj] = mean(reform(wrf_no2_emi[ii,jj,*]),/NAN) wno2_vmr_avg[ii,jj] = mean(reform(wrf_no2_vmr[ii,jj,*]),/NAN) endfor endfor ;print, ' ------------------------------------------ ' ; ----------------------------------------------------------------------------------------- ;- ploting section int = strcompress(string(hh1)+'-'+string(hh2), /remove_all) outfl = strcompress(plot_path+'WRF_TEMPO_trop_colum_'+species+'_'+scenario[n]+$ '_'+dd+'aug_avg_'+int+'UTC.eps' ,/remove_all) spl_trop_col_wrf, outfl, trop_col_avg, true_trop_col_avg, lat_all, lon_all, $ wno2_vmr_avg, wno2_emi_avg, wrf_lat, wrf_lon, dd, int, scenario[n], spc_ttl ; ----------------------------------------------------------------------------------------- outfile = strcompress(plot_path+'WRF_TEMPO_scatter_'+species+'_'+scenario[n]+ $ '_'+dd+'aug_avg_'+int+'UTC.eps' ,/remove_all) scatter_tempo_wrf_density, outfile, trop_col_avg, true_trop_col_avg, lat_all, lon_all, $ wno2_vmr_avg, wno2_emi_avg, wrf_lat, wrf_lon, dd, int, scenario[n], spc_ttl ; ----------------------------------------------------------------------------------------- endif ; ------------------------------------------------------------------------------------------- endfor ; num WRF simulatons ;endfor ; num hours ; use this loop for bulk hourly plots ; -------------------------------------------------------------------------------------------------------- end