! cd /aosc/horse/carton/fluxes/
*
* compute monthly mean fields
*
!cdo monmean /aosc/greenland/forcing/era_i/era_i_rain_1dy_19790101_20181231.nc era_i_rain_1x1_1mn.nc
!cdo monmean /aosc/greenland/forcing/era_i/era_i_snowfall_1dy_19790101_20181231.nc era_i_tsnflt_1x1_1mn.nc
!cdo seldate,1980-01-01,2010-12-31 /aosc/greenland/forcing/CERA-20C/combined/CERA20C_rain_1901_2010.nc CERA20C_rain_1980_2010_dy.nc
!cdo monmean CERA20C_rain_1980_2010_dy.nc  CERA20C_rain_1980_2010_mn.nc
!cdo seldate,1980-01-01,2010-12-31 /aosc/greenland/forcing/CERA-20C/combined/CERA20C_sf_1901_2010.nc CERA20C_sf_1980_2010_dy.nc
!cdo monmean CERA20C_sf_1980_2010_dy.nc CERA20C_sf_1980_2010_mn.nc 
!rm CERA20C_rain_1980_2010_dy.nc CERA20C_sf_1980_2010_dy.nc
*
* Compare downwelling rain and snow from a new reanalysis with ERA-I.  
*
* The upper two panels compare the 
* LEFT:  difference in record means and 
* RIGHT:  difference in monthly variances.
* The lower panel compares a one-year long time series at a point
* (42N, 163E) in the central North Pacific
reinit
set display color white
c
* load the color table
/homes/metofac/carton/colors
*
* plot rain
*
sdfopen  CERA20C_rain_1980_2010_mn.nc
mean = ave(rain,time=jan1980,time=dec2018)
var = ave((rain-mean)*(rain-mean),time=jan1980,time=dec2018))
set lon 162
set lat 42
set t 1 last
ts = rain
close 1
sdfopen /aosc/horse/carton/fluxes/era_i/era_i_rain_1x1_1mn.nc
set lon 162
set lat 42
set t 1 last
tss = fpr
set  lon 0 360
set lat -90 90
set t 1
means = ave(fpr,time=jan1980,time=dec2018)
vars = ave((fpr-means)*(fpr-means),time=jan1980,time=dec2018))
means2 = (lterp(mean,means)-means)
vars2 = (lterp((var),vars)-vars)
set vpage 0 4.5 8.5 10.5
set xlopts 1 3 .15
set ylopts 1 3 .15
set grads off
set grid off
set cmark 0
set gxout shaded
set clevs -2.0 -1. -.5 -.2 .2 .5 1. 2.0 
set ccols 47 45 43 42 0 61 63 65 67
d means2*1.e5
basemap L 15 15
cbarn
draw title AVG_1-AVG_era (x10^5)
set ccolor 1
set vpage 4 8.5 8.5 10.5
set xlopts 1 3 .15
set ylopts 1 3 .15
set grads off
set grid off
set cmark 0
set clevs  -2. -1 -.5 -.1 .1 .5 1 2  
set ccols 47 45 43 42 0 61 63 65 67
d vars2*1.e9
basemap L 15 15
cbarn
draw title VAR_1-VAR_era (x10^9)
set vpage off
set strsiz .14
set string 1 c 4 0
draw string 4.25 10.8  RAIN (kg/(m2s))
set x 1
set y 1
set time jan1980 dec2019
set vpage 0 4.25 6.5 8.75
set cmark 0
set ccolor 1
set cthick 8
set grads off
set grid off
set axlim 00 15
set yaxis 00 15 5
set xlopts 1 3 .15
set ylopts 1 3 .15
set ylab on
set ylpos 0 l 
d ts*1.e5
set axlim 00 15
set yaxis 00 15 5
set xlopts 1 3 .15
set ylopts 1 3 .15
set ylab on
set ylpos 0 l 
set cthick 2
set ccolor 2
set cmark 0
d tss*1.e5
set vpage off
set strsiz .13
set string 1 c 4 0
draw string 2.45606 8.55556 MONTHLY 42N,162E
set string 2 l 3 0
set strsiz .09
draw string 2.4787 8.21605 -- ERA-I
set line 1 1 1
draw line 2 6.3 7 6.3
close 1
*
* plot snow
*
sdfopen  CERA20C_sf_1980_2010_mn.nc 
mean = ave(sf,time=jan1980,time=dec2018)
var = ave((sf-mean)*(sf-mean),time=jan1980,time=dec2018))
set lon 162
set lat 42
set t 1 last
ts = sf
close 1
sdfopen /aosc/horse/carton/fluxes/era_i/era_i_tsnflt_1x1_1mn.nc
set lon 162
set lat 42
set t 1 last
tss = sf
set  lon 0 360
set lat -90 90
set t 1
means = ave(sf,time=jan1980,time=dec2018)
vars = ave((sf-means)*(sf-means),time=jan1980,time=dec2018))
means2 = (lterp(mean,means)-means)
vars2 = (lterp((var),vars)-vars)
set vpage 0 4.5 3.5 5.5
set xlopts 1 3 .15
set ylopts 1 3 .15
set grads off
set grid off
set cmark 0
set gxout shaded
set clevs  -2  -1 -.5 .5 1 2  
set ccols  47 45 43  0  63 65 67
d means2*1.e5
basemap L 15 15
cbarn
draw title AVG_1-AVG_era (x1.e5)
set ccolor 1
set vpage 4 8.5 3.5 5.5
set xlopts 1 3 .15
set ylopts 1 3 .15
set grads off
set grid off
set cmark 0
set clevs  -4  -2 -1 1 2 4  
set ccols  47 45 43  0  63 65 67
d vars2*1.e11
basemap L 15 15
cbarn
draw title VAR_1-VAR_era (x1.e11)
set vpage off
set strsiz .14
set string 1 c 4 0
draw string 4.25 5.8  SNOW (kg/(m2s))
set x 1
set y 1
set time jan1980 dec2019
set vpage 0 4.25 1.5 3.75
set cmark 0
set grads off
set grid off
set axlim 0 3
set yaxis 0 3 1
set xlopts 1 3 .15
set ylopts 1 3 .15
set ylab on
set ylpos 0 l 
set cthick 8
d ts*1.e5
*set vpage 3.75 8 1.5 3.75
set cmark 0
set axlim 0 3
set yaxis 0 3 1
set xlopts 1 3 .15
set ylopts 1 3 .15
set ylab off
set ccolor 2
set cthick 2
set grads off
set grid off
d tss*1.e5
set vpage off
set strsiz .13
set string 1 c 4 0
draw string  2.45606 3.4856 MONTHLY 42N,162E
set string 2 l 3 0
set strsiz .09
draw string  1.08655 3.27058 -- ERA-I
gxprint /aosc/greenland/forcing/CERA-20C/combined/rain_snow_flux-comparison.pdf