; ; Plot a monthly summary ; pro histogram_average_cfrac_height;,site,searchdate,avg_int site='sgp' searchdate='200004' avg_int='1hr' ; ; common containing the monthly averaged data ; common monthly_averaged_data, base_time,time_offset,height,reflectivity,$ reflectp, reflectc, cbase, ctop, velocity, clutter,maxclutter,spectrum,$ lwp,vceil,numtimes,numheights,allclutter,goodclutter,temp,press,$ vap,dir,dsol,sfpres,numlayers,twsflux,zenith,solarratio,sfcalbedo,$ column_cfrac ; ; common containing the monthly averaged data ; common monthly_averaged_data_rad, base_time_rad,time_offset_rad,height_rad,$ numtimes_rad,numheights_rad,$ IWC_best_estimate,re_ice_best_estimate,ice_micro_error_estimate, $ iwc_source,ice_effective_radius_source,ice_tau_solar, $ ice_omega_solar,ice_g_solar,ice_tau_IR,ice_omega_IR,ice_g_IR, $ ice_RadiationParameterization_source_solar,ice_RadiationParameterization_source_ir, $ lwc_best_estimate,re_liquid_best_estimate,liquid_micro_error_estimate, $ lwc_source,liquid_effective_radius_source,liquid_tau_solar,liquid_omega_solar, $ liquid_g_solar,liquid_tau_IR,liquid_omega_IR,liquid_g_IR, $ liquid_RadiationParameterization_source_solar,$ liquid_RadiationParameterization_source_ir,$ ice_fraction_CCM3,solar_flux_down_clear,solar_flux_up_clear, $ IR_flux_up_clear,IR_flux_down_clear,downwelling_solar_diffuse_clear, $ downwelling_solar_direct_clear,TOA_solar_up_clear,TOA_solar_down_clear, $ TOA_IR_up_clear,solar_flux_down,solar_flux_up,IR_flux_up,IR_flux_down, $ downwelling_solar_diffuse,downwelling_solar_direct,TOA_solar_up,TOA_solar_down, $ TOA_IR_up ; ; Read in the monthly averaged data ; get_average_data,site,searchdate,avg_int get_average_data_rad,site,searchdate,avg_int ; ; Convert heights to km ; height=height/1000.0 ; ; Create a 2d array of heights to match cfrac ; heights2d=make_array(numtimes,numheights,/float,value=0) for i=0,numtimes-1 do begin heights2d[i,*]=height endfor ;********************************************************************** ; Hour of the day vs column cloud fraction ;********************************************************************** ; ; Strings ; dataname='Cloud Fraction vs. Height' datatag='cfrac_height' cb_title='counts' ;************************** ; Calculate cloud fraction ;************************** ; ; First if statement for twp,nsa ; if site eq 'nsa' or site eq 'twpc1' or site eq 'twpc2' then begin ; Percentage of points that are cloudy cfrac=float(reflectc)/float(reflectp)*100. ; Where possible points is 0, there is missing data. ; Where possible points is -8888, there is missing data ; Give these points the missing data value,-8888,gray result=where(reflectp le 0.0 or reflectp eq -8888,count) print,'reflectp is 0 for ',count,' points' if count ne 0 then cfrac[result]=-8888 ; ; Second if statement for sgp ; endif else if site eq 'sgp' or site eq 'crystal' then begin ; Percentage of points that are cloudy cfrac=float(goodclutter)/float(allclutter)*100. ; Where possible points is 0, there is missing data. ; Where possible points is -8888, there is missing data. ; Give these points the missing data value,-8888,gray result=where(allclutter eq 0.0 or allclutter eq -8888,count) print,'allclutter is 0 for ',count,' points' if count ne 0 then cfrac[result]=-8888 endif ;************************** ; End of cloud fraction ;************************** ; ; Set up cfrac bins with mid points 5,15,25,35,45,55,65 % ; cfrac_bins=5.0 dcfrac_bins=5.0 end_cfrac_bins=105.0 while cfrac_bins[n_elements(cfrac_bins)-1] lt end_cfrac_bins $ do cfrac_bins=[cfrac_bins,cfrac_bins[n_elements(cfrac_bins)-1]+(2.*dcfrac_bins)] ; ; Set up height bins with mid points 25,75,125 km ; height_bins=1.0 dheight_bins=0.5 end_height_bins=15.0 while height_bins[n_elements(height_bins)-1] lt end_height_bins $ do height_bins=[height_bins,height_bins[n_elements(height_bins)-1]+(2.*dheight_bins)] ; ; Set up the 2D histogram arrays ; cfrac_height_freq= $ make_array(n_elements(cfrac_bins),n_elements(height_bins),/float,value=0) ; ; Build the 2D histogram for temp-thick ; for j=0,n_elements(cfrac_bins)-1 do begin for k=0,n_elements(height_bins)-1 do begin result=where(cfrac ge cfrac_bins[j]-dcfrac_bins and $ cfrac lt cfrac_bins[j]+dcfrac_bins and $ heights2d ge height_bins[k]-dheight_bins and $ heights2d lt height_bins[k]+dheight_bins,count) cfrac_height_freq[j,k]=count ;if hour_bins[j]-dhour_bins eq 5 then begin ;print,cfrac_bins[j]-dcfrac_bins,cfrac_bins[j]+dcfrac_bins,$ ; height_bins[k]-dheight_bins,height_bins[k]+dheight_bins,$ ; cfrac_height_freq[j,k],format='(5(F6.2,1X))' ;print,cfrac[result] ;print,heights2d[result] ;stop ;endif endfor ;end of loop through thick bins endfor ;end of loop through temp bins ; Change the counts array to a frequency array ;hour_cfrac_freq=float(hour_cfrac_freq)/float(total(hour_cfrac_freq)) ; Frequency array per hour ;for j=0,n_elements(hour_bins)-1 do begin ; hour_cfrac_freq[j,*]=float(hour_cfrac_freq[j,*])/float(total(hour_cfrac_freq[j,*])) ;endfor ;****************** ; Set up the plot ;****************** window,0,/pixmap loadct,39 !p.background=!d.n_colors-1 !p.charsize=1.5 !p.multi=[1,1,1] !x.style=1 !y.style=1 !x.ticklen=-0.02 !y.ticklen=-0.02 pos_vect=[0.15,0.15,0.7,0.8] cb1=[0.9,0.15,0.95,0.8] dmax=1000.0;max(cfrac_height_freq) dmin=min(cfrac_height_freq) contour,cfrac_height_freq,cfrac_bins,height_bins,/nodata,$ xstyle=4,ystyle=4,position=pos_vect px=!x.window*!d.x_vsize py=!y.window*!d.y_vsize sx=px[1]-px[0]+1 sy=py[1]-py[0]+1 ; Scale and create the image image=bytscl(cfrac_height_freq,max=dmax,min=dmin,top=254) result=where(cfrac_height_freq eq 0) image[result]=255 tv,congrid(image,sx,sy),px[0],py[0] ; Set up the axes contour,cfrac_height_freq,cfrac_bins,height_bins,/noerase,/nodata,$ color=0,xtitle='cfrac',ytitle='height',$ xstyle=1,ystyle=1,position=pos_vect,$ xrange=[cfrac_bins[0]-dcfrac_bins,cfrac_bins[n_elements(cfrac_bins)-1]+dcfrac_bins],$ yrange=[height_bins[0]-dheight_bins,height_bins[n_elements(height_bins)-1]+dheight_bins] colorbar_david_fanning,maxrange=dmax,minrange=dmin,$ title='counts',$ ncolors=254,format='(F6.2)',$ vertical=1,color=0,position=cb1,charsize=1.5,$ divisions=8 ; ; Add a title ; smonth=strmid(searchdate,4,2) syear=strmid(searchdate,0,4) sday='01' if smonth eq 1 then smon='Jan' if smonth eq 2 then smon='Feb' if smonth eq 3 then smon='Mar' if smonth eq 4 then smon='Apr' if smonth eq 5 then smon='May' if smonth eq 6 then smon='Jun' if smonth eq 7 then smon='Jul' if smonth eq 8 then smon='Aug' if smonth eq 9 then smon='Sep' if smonth eq 10 then smon='Oct' if smonth eq 11 then smon='Nov' if smonth eq 12 then smon='Dec' if site eq 'sgp' then begin capsite='SGP' datastream='average.'+avg_int endif else if site eq 'nsa' then begin capsite='NSA' datastream='average.'+avg_int endif else if site eq 'twpc1' then begin capsite='TWP Manus' datastream='average.'+avg_int endif else if site eq 'twpc2' then begin capsite='TWP Nauru' datastream='average.'+avg_int endif syr=strcompress(string(syear),/remove_all) xyouts,0.5,0.95,align=0.5,charsize=1.7,/normal,color=0,$ 'Histogram of Hourly Averages!C'+ $ capsite+' - '+smon+' '+syr + '!C' + dataname ; ; Create the gif file ; strmon=string(smonth) if smonth lt 10 then strmon='0'+strmon strday=string(sday) if sday lt 10 then strday='0'+strday curdate=strcompress(syr+strmon+strday,/remove_all) image_dir='/data/mace4/arm_data/temp_image/' image_name=site+'.'+datastream+'.'+curdate+'.'+ $ datatag+'.2dhist' print,image_dir+image_name write_gif,image_dir+image_name+'.gif',tvrd( ) ; ; Remove the .png file ; unix_command='rm '+image_dir+image_name+'.png' spawn,unix_command ;**************************************** result=where(cfrac gt -8888,count) dmin=min(cfrac[result]) dmax=max(cfrac[result]) contour,cfrac,time_offset,height,/nodata,$ xstyle=4,ystyle=4,position=pos_vect px=!x.window*!d.x_vsize py=!y.window*!d.y_vsize sx=px[1]-px[0]+1 sy=py[1]-py[0]+1 ; Scale and create the image image=bytscl(cfrac,max=dmax,min=dmin,top=254) result=where(cfrac eq 0) image[result]=255 tv,congrid(image,sx,sy),px[0],py[0] ; Set up the axes contour,cfrac,time_offset,height,/noerase,/nodata,$ color=0,xtitle='cfrac',ytitle='height',$ xstyle=1,ystyle=1,position=pos_vect ;colorbar_david_fanning,maxrange=dmax,minrange=dmin,$ ; title='counts',$ ; ncolors=254,format='(F6.2)',$ ; vertical=1,color=0,position=cb1,charsize=1.5,$ ; divisions=8 ; ; Add a title ; smonth=strmid(searchdate,4,2) syear=strmid(searchdate,0,4) sday='01' if smonth eq 1 then smon='Jan' if smonth eq 2 then smon='Feb' if smonth eq 3 then smon='Mar' if smonth eq 4 then smon='Apr' if smonth eq 5 then smon='May' if smonth eq 6 then smon='Jun' if smonth eq 7 then smon='Jul' if smonth eq 8 then smon='Aug' if smonth eq 9 then smon='Sep' if smonth eq 10 then smon='Oct' if smonth eq 11 then smon='Nov' if smonth eq 12 then smon='Dec' if site eq 'sgp' then begin capsite='SGP' datastream='average.'+avg_int endif else if site eq 'nsa' then begin capsite='NSA' datastream='average.'+avg_int endif else if site eq 'twpc1' then begin capsite='TWP Manus' datastream='average.'+avg_int endif else if site eq 'twpc2' then begin capsite='TWP Nauru' datastream='average.'+avg_int endif syr=strcompress(string(syear),/remove_all) xyouts,0.5,0.95,align=0.5,charsize=1.7,/normal,color=0,$ 'Histogram of Hourly Averages!C'+ $ capsite+' - '+smon+' '+syr + '!C' + dataname ; ; Create the gif file ; strmon=string(smonth) if smonth lt 10 then strmon='0'+strmon strday=string(sday) if sday lt 10 then strday='0'+strday curdate=strcompress(syr+strmon+strday,/remove_all) image_dir='/data/mace4/arm_data/temp_image/' image_name=site+'.'+datastream+'.'+curdate+'.'+ $ 'cfrac_only'+'.2dhist' print,image_dir+image_name write_gif,image_dir+image_name+'.gif',tvrd( ) ; ; Remove the .png file ; unix_command='rm '+image_dir+image_name+'.png' spawn,unix_command ;**************************************** dmin=min(heights2d) dmax=max(heights2d) contour,heights2d,time_offset,height,/nodata,$ xstyle=4,ystyle=4,position=pos_vect px=!x.window*!d.x_vsize py=!y.window*!d.y_vsize sx=px[1]-px[0]+1 sy=py[1]-py[0]+1 ; Scale and create the image image=bytscl(heights2d,max=dmax,min=dmin,top=254) ;result=where(heights2d eq 0) ;image[result]=255 tv,congrid(image,sx,sy),px[0],py[0] ; Set up the axes contour,heights2d,time_offset,height,/noerase,/nodata,$ color=0,xtitle='cfrac',ytitle='height',$ xstyle=1,ystyle=1,position=pos_vect colorbar_david_fanning,maxrange=dmax,minrange=dmin,$ title='counts',$ ncolors=254,format='(F6.2)',$ vertical=1,color=0,position=cb1,charsize=1.5,$ divisions=8 ; ; Add a title ; smonth=strmid(searchdate,4,2) syear=strmid(searchdate,0,4) sday='01' if smonth eq 1 then smon='Jan' if smonth eq 2 then smon='Feb' if smonth eq 3 then smon='Mar' if smonth eq 4 then smon='Apr' if smonth eq 5 then smon='May' if smonth eq 6 then smon='Jun' if smonth eq 7 then smon='Jul' if smonth eq 8 then smon='Aug' if smonth eq 9 then smon='Sep' if smonth eq 10 then smon='Oct' if smonth eq 11 then smon='Nov' if smonth eq 12 then smon='Dec' if site eq 'sgp' then begin capsite='SGP' datastream='average.'+avg_int endif else if site eq 'nsa' then begin capsite='NSA' datastream='average.'+avg_int endif else if site eq 'twpc1' then begin capsite='TWP Manus' datastream='average.'+avg_int endif else if site eq 'twpc2' then begin capsite='TWP Nauru' datastream='average.'+avg_int endif syr=strcompress(string(syear),/remove_all) xyouts,0.5,0.95,align=0.5,charsize=1.7,/normal,color=0,$ 'Histogram of Hourly Averages!C'+ $ capsite+' - '+smon+' '+syr + '!C' + dataname ; ; Create the gif file ; strmon=string(smonth) if smonth lt 10 then strmon='0'+strmon strday=string(sday) if sday lt 10 then strday='0'+strday curdate=strcompress(syr+strmon+strday,/remove_all) image_dir='/data/mace4/arm_data/temp_image/' image_name=site+'.'+datastream+'.'+curdate+'.'+ $ 'heights2d'+'.2dhist' print,image_dir+image_name write_gif,image_dir+image_name+'.gif',tvrd( ) ; ; Remove the .png file ; unix_command='rm '+image_dir+image_name+'.png' spawn,unix_command stop end