;
;  This program plots a daily summary of the cloud forcing
;
pro plot_average_daily_1d_cforce2,site,date,avg_int,ret_type

;*************************************
;*************************************
;  INPUTS
;*************************************
;*************************************

;  Site is the following 'sgp','nsa','twpc1','twpc2','twpc3'
;site='sgp'

;  Date in the form yyyymmdd 
;date='19970101'

;  Averaging interval
;avg_int='5min';'1hr';'5min'

;  Retrieval version
;ret_type='ciret4'

;**************************************
;  Commons
;**************************************

;  Color common
common colors,r_orig,g_orig,b_orig,r_curr,g_curr,b_curr

;*****************************************
;  Constants
;*****************************************

;  Indicate which files you want to use, the monthly or daily files.
m_or_d='day'  ;day,month

;  Missing data value
missing=-8888

;  Seconds in a day
sec_per_day=double(24L*60L*60L)

;  Get IDL julian day of day 1
day1=julday(1,1,1970,0,0,0)

;  Arm data directory
arm_dir='/uufs/chpc.utah.edu/common/home/mace_grp/data/arm/'+site+'/'

;*****************************************
;  Read in the monthly averaged data
;*****************************************

;  Make sure date is a string
date=strcompress(date,/remove_all)

;  Set the fill flag - 1d data doesn't need to be filled
fill_flag='no'

;  List of variables to read in  
varname=['base_time','time_offset',$
         'net_toa_cld_forcing',$
         'net_toa_cld_forcing_obs',$
         'longwave_toa_cld_forcing',$
         'longwave_toa_cld_forcing_obs',$
         'shortwave_toa_cld_forcing',$
         'shortwave_toa_cld_forcing_obs',$
         'net_atm_cld_forcing',$
         'net_atm_cld_forcing_obs',$
         'longwave_atm_cld_forcing',$
         'longwave_atm_cld_forcing_obs',$
         'shortwave_atm_cld_forcing',$
         'shortwave_atm_cld_forcing_obs',$
         'net_sfc_cld_forcing',$
         'net_sfc_cld_forcing_obs',$
         'longwave_sfc_cld_forcing',$
         'longwave_sfc_cld_forcing_obs',$
         'shortwave_sfc_cld_forcing',$
         'shortwave_sfc_cld_forcing_obs']

;  Get the variables in the list from the file
get_average_data_list,varname,data,$
  site,date,avg_int,ret_type,m_or_d,fill_flag

;  Pull out the individual data arrays
for v=0,n_elements(varname)-1 do begin
  xstr=varname[v]+'=data.'+varname[v]
  result=execute(xstr)
endfor

;  Get the size of the arrays
numtimes=n_elements(time_offset)

;*******************************************
;  Convert time to julian day for plotting
;*******************************************

;  Calculate julian_day
julian_day=double(day1+((base_time+time_offset)/sec_per_day))

;  Calculate the standard date arrays
caldat,julian_day,mm,dd,yy,hh,mi,ss

;  File start and end time
syear=yy[0]
smonth=mm[0]
sday=dd[0]
shour=hh[0]
smin=mi[0]
ssec=ss[0]
;print,'file start',syear,smonth,sday,shour,smin,ssec

eyear=yy[numtimes-1]
emonth=mm[numtimes-1]
eday=dd[numtimes-1]
ehour=hh[numtimes-1]
emin=mi[numtimes-1]
esec=ss[numtimes-1]
;print,'file end',eyear,emonth,eday,ehour,emin,esec

;********************************
;  Set up the plot device
;********************************

;  Set up the device
set_plot,'z'

;  xsize of the plot window
xwin=800  

;  ysize of the plot window
ywin=800  

;  Size the window
device,set_resolution=[xwin,ywin]

;********************************
;  Set up the colors
;  0-253 rainbow
;  254   gray
;  255   white
;********************************

;  Load the color table
loadct,39   

;  Make the background white
!p.background=!d.n_colors-1  

;  Make gray for missing data in color 253
r_curr[254]=200
g_curr[254]=200
b_curr[254]=200

;********************************
;  Set up the plot
;********************************

;  Date label format
dummy=label_date(date_format=['%H:%I'])

;  Xticks point out from plot
!x.ticklen=-0.02
!x.style=1
!y.style=1

;  Char sizes
tsz=1.0
stsz=0.9
!p.charsize=0.8

;********************************
;  Position vectors
;********************************

xl=0.05  &  xr=0.99
yb=0.05  &  yt=0.97
sx=0.05        
sy=0.03 ;0.05        
numplots_x=2  ;number of plots in x direction
numplots_y=5  ;number of plots in y direction
position_plots,xl,xr,yb,yt,sx,sy,numplots_x,numplots_y,pos

;  String positioning
dr=0.09  ;string offset from the right xcoor
dl=0.025  ;string offset from the left xcoor
dt=0.025  ;string offset from the top ycoor

;*********************************************
;  Plot 0 
;*********************************************

;  Plot location
pnum=0

alldata=[net_toa_cld_forcing,net_toa_cld_forcing_obs]
result=where(alldata gt missing and alldata lt 1e36,count)
if count gt 0 then begin
  dmin=min(alldata[where(alldata gt missing and alldata lt 1e36)])
  dmax=max(alldata[where(alldata gt missing and alldata lt 1e36)])
  plot,julian_day,alldata,color=0,$
    yrange=[dmin,dmax],position=pos[pnum,*],$
    xtickunits='time',xtickformat='label_date',$
    title='net_toa_cld_forcing dmin='+strcompress(dmin,/remove_all) +$
    ' dmax='+strcompress(dmax,/remove_all)
  oplot,julian_day,net_toa_cld_forcing_obs,color=250,psym=4
endif
;*********************************************
;  Plot 1 
;*********************************************

;  Plot location
pnum=1

alldata=[longwave_toa_cld_forcing,longwave_toa_cld_forcing_obs]
result=where(alldata gt missing and alldata lt 1e36,count)
if count gt 0 then begin
  dmin=min(alldata[where(alldata gt missing)])
  dmax=max(alldata[where(alldata gt missing)])
  plot,julian_day,longwave_toa_cld_forcing,color=0,/noerase,$
    yrange=[dmin,dmax],position=pos[1,*],$
    xtickunits='time',xtickformat='label_date',$
    title='longwave_toa_cld_forcing dmin='+strcompress(dmin,/remove_all) +$
    ' dmax='+strcompress(dmax,/remove_all)
  oplot,julian_day,longwave_toa_cld_forcing_obs,color=250,psym=4
endif

;**********************************
;  Plot 2 - ir_flux_up_clear
;**********************************

;  Plot location
pnum=2

alldata=[shortwave_toa_cld_forcing,shortwave_toa_cld_forcing_obs]
result=where(alldata gt missing and alldata lt 1e36,count)
if count gt 0 then begin
  dmin=min(alldata[where(alldata gt missing)])
  dmax=max(alldata[where(alldata gt missing)])
  plot,julian_day,shortwave_toa_cld_forcing,color=0,/noerase,$
    yrange=[dmin,dmax],position=pos[2,*],$
    xtickunits='time',xtickformat='(A1)',$
    title='shortwave_toa_cld_forcing dmin='+strcompress(dmin,/remove_all) +$
    ' dmax='+strcompress(dmax,/remove_all)
  oplot,julian_day,shortwave_toa_cld_forcing_obs,color=250,psym=4
endif


;**********************************
;  Plot 3 
;**********************************

;  Plot location
pnum=3

alldata=[net_atm_cld_forcing,net_atm_cld_forcing_obs]
result=where(alldata gt missing and alldata lt 1e36,count)
if count gt 0 then begin
  dmin=min(alldata[where(alldata gt missing)])
  dmax=max(alldata[where(alldata gt missing)])
  plot,julian_day,net_atm_cld_forcing,color=0,/noerase,$
    yrange=[dmin,dmax],position=pos[3,*],$
    xtickunits='time',xtickformat='(A1)',$
    title='net_atm_cld_forcing dmin='+strcompress(dmin,/remove_all) +$
    ' dmax='+strcompress(dmax,/remove_all)
  oplot,julian_day,net_atm_cld_forcing_obs,color=250,psym=4
endif

;**********************************************
;  Plot 4 - downwelling_solar_diffuse_clear
;**********************************************

;  Plot location
pnum=4

alldata=[longwave_atm_cld_forcing,longwave_atm_cld_forcing_obs]
result=where(alldata gt missing and alldata lt 1e36,count)
if count gt 0 then begin
  dmin=min(alldata[where(alldata gt missing)])
  dmax=max(alldata[where(alldata gt missing)])
  plot,julian_day,longwave_atm_cld_forcing,color=0,/noerase,$
    yrange=[dmin,dmax],position=pos[4,*],$
    xtickunits='time',xtickformat='(A1)',$
    title='longwave_atm_cld_forcing dmin='+strcompress(dmin,/remove_all) +$
    ' dmax='+strcompress(dmax,/remove_all)
  oplot,julian_day,longwave_atm_cld_forcing_obs,color=250,psym=4
endif

;*********************************************
;  Plot 5 - downwelling_solar_direct_clear
;*********************************************

;  Plot location
pnum=5

alldata=[shortwave_atm_cld_forcing,shortwave_atm_cld_forcing_obs]
result=where(alldata gt missing and alldata lt 1e36,count)
if count gt 0 then begin
  dmin=min(alldata[where(alldata gt missing)])
  dmax=max(alldata[where(alldata gt missing)])
  plot,julian_day,shortwave_atm_cld_forcing,color=0,/noerase,$
    yrange=[dmin,dmax],position=pos[5,*],$
    xtickunits='time',xtickformat='(A1)',$
    title='shortwave_atm_cld_forcing dmin='+strcompress(dmin,/remove_all) +$
    ' dmax='+strcompress(dmax,/remove_all)
  oplot,julian_day,shortwave_atm_cld_forcing_obs,color=250,psym=4
endif

;*********************************************
;  Plot 6 
;*********************************************

;  Plot location
pnum=6

alldata=[net_sfc_cld_forcing,net_sfc_cld_forcing_obs]
result=where(alldata gt missing and alldata lt 1e36,count)
if count gt 0 then begin
dmin=min(alldata[where(alldata gt missing)])
dmax=max(alldata[where(alldata gt missing)])
plot,julian_day,net_sfc_cld_forcing,color=0,/noerase,$
  yrange=[dmin,dmax],position=pos[6,*],$
  xtickunits='time',xtickformat='(A1)',$
  title='net_sfc_cld_forcing dmin='+strcompress(dmin,/remove_all) +$
  ' dmax='+strcompress(dmax,/remove_all)
oplot,julian_day,net_sfc_cld_forcing_obs,color=250,psym=4,symsize=0.5
endif


;*********************************************
;  Plot 7 
;*********************************************

;  Plot location
pnum=7

alldata=[longwave_sfc_cld_forcing,longwave_sfc_cld_forcing_obs]
result=where(alldata gt missing and alldata lt 1e36,count)
if count gt 0 then begin
dmin=min(alldata[where(alldata gt missing)])
dmax=max(alldata[where(alldata gt missing)])
plot,julian_day,longwave_sfc_cld_forcing,color=0,/noerase,$
  yrange=[dmin,dmax],position=pos[7,*],$
  xtickunits='time',xtickformat='(A1)',$
  title='longwave_sfc_cld_forcing dmin='+strcompress(dmin,/remove_all) +$
  ' dmax='+strcompress(dmax,/remove_all)
oplot,julian_day,longwave_sfc_cld_forcing_obs,color=250,psym=4,symsize=0.5
endif

;*********************************************
;  Plot 8 
;*********************************************

;  Plot location
pnum=8

alldata=[shortwave_sfc_cld_forcing,shortwave_sfc_cld_forcing_obs]
result=where(alldata gt missing and alldata lt 1e36,count)
if count gt 0 then begin
dmin=min(alldata[where(alldata gt missing)])
dmax=max(alldata[where(alldata gt missing)])
plot,julian_day,shortwave_sfc_cld_forcing,color=0,/noerase,$
  yrange=[dmin,dmax],position=pos[8,*],$
  xtickunits='time',xtickformat='(A1)',$
  title='shortwave_sfc_cld_forcing dmin='+strcompress(dmin,/remove_all) +$
  ' dmax='+strcompress(dmax,/remove_all)
oplot,julian_day,shortwave_sfc_cld_forcing_obs,color=250,psym=4,symsize=0.5
endif

;************************************
;  Plot title 
;************************************

;  Add a title
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'
endif else if site eq 'nsa' then begin
  capsite='NSA'
endif else if site eq 'twpc1' then begin
  capsite='TWPC1'
endif else if site eq 'twpc2' then begin
  capsite='TWPC2'
endif else if site eq 'twpc3' then begin
  capsite='TWPC3'
endif

str_year=string(syear,format='(I4)')
str_month=string(smonth,format='(I02)')
str_day=string(sday,format='(I02)')

title_str=capsite+' '+avg_int+' '+ret_type+$
  '!C'+str_day+' '+smon+' '+str_year
xyouts,0.90,0.97,align=0.5,charsize=1.7,/normal,color=0,$
title_str

comment_str='The red values are from !C observations'
xyouts,0.80,0.92,align=0.5,charsize=1.7,/normal,color=0,$
comment_str

;**********************************
;  Create gif file
;**********************************

if avg_int eq '1hr' or avg_int eq '6hr' then begin
  imagedir=arm_dir+site+'.average.'+avg_int+'.'+ret_type+$
           '/images_daily/'+str_year+'/'
endif else if avg_int eq '5min' and ret_type eq 'mine' then begin
  imagedir=arm_dir+site+'.average.5min/images_daily/'+str_year+'/'
endif else if avg_int eq '5min' and ret_type ne 'mine' then begin
  imagedir=arm_dir+site+'.average.5min.'+ret_type+'/images_daily/'+$
  str_year+'/'
endif else begin
  imagedir=''
endelse

if ret_type eq 'mine' then begin
  imagename=site+'.average.'+avg_int+'.'+date+'.1d_cforce'
endif else begin
  imagename=site+'.average.'+avg_int+'.'+ret_type+'.'+date+'.1d_cforce'
endelse

print,imagedir+imagename
write_gif,imagedir+imagename+'.gif',tvrd( )

end