source: trunk/src/flux_evaluation_tpr_map.pro

;+
;
; .. _flux_evaluation_tpr_map.pro:
;
; ===========================
; flux_evaluation_tpr_map.pro
; ===========================
;
; DESCRIPTION
; ===========
;
; SEE ALSO
; ========
;
; :ref:`project_profile.sh`
; :ref:`project_init.pro`
; :ref:`cm_project.pro`
;
; TODO
; ====
;
; make it work : need output of flux_evaluation_tpr_timeseries.pro
;
; uncomment var
;
; coding rules
;
; EXAMPLES
; ========
;
; .. code-block:: idl
;
;    date1 = 20000101L
;    date2 = 20091231L
;    flux_evaluation_tpr_map, date1, date2
;
; EVOLUTIONS
; ==========
;
; $Id$
;
; $URL$
;
; - fplod 20111130T142224Z cratos (Linux)
;
;   * try to make it work on my account on cratos
;   * remove x_site_location and x_site_location because already exist
;
;-
pro flux_evaluation_tpr_map, $
;                    var,    $ ;;  flux variable (swr, lwr, lhf, shf) to calculate the statistics
date1,  $ ;;  start date (in Julian date. eg. 20000101)
date2     ;;  end date (in Julian date eg. 20091231)

@common
@cm_project
reinitplt, /z,/invert
key_portrait = 1
openps, FILENAME = project_od_env+'flux_evaluation_tpr_map.ps'
;
;; part to change

min_obs=10.  ;; this will allow to calculate statistics at locations with more than 180. valid observation

;; choose the appropriate min and max values for the following.

bias_mi=-5     &  bias_ma=5    & bias_int=0.5
std_mi=0.5      &  std_ma=1.5    & std_int=0.05
rmsd_mi=0       &  rmsd_ma=15    & rmsd_int=1.5
cor_mi=0.5      &  cor_ma=1.     & cor_int=0.025

;; TPR locations.  This needs to be updated with time since more locations are added to the array.

sitelist=['5n165e','8s67e','12s55e', '8s55e', '8s80.5e', '1.5s80.5e', '0n80.5e', '1.5n80.5e', '1.5s90e', $
           '0n90e', '1.5n90e', '4n90e','8n90e','12n90e', '15n90e', '5s95e', $
           '8s165e', '8s180w',  '8s155w', '8s125w', '8s110w', '8s95w',  '5s156e', '5s165e', '5s180w', '5s170w', $
          '5s155w', '5s140w', '5s125w', '5s110w', '5s95w', '2s156e', '2s165e', '2s180w', '2s170w', '2s155w', '2s140w', $
          '2s125w', '2s110w', '2s95w', '0n147e', '0n156e', '0n165e', '0n180w', '0n170w', '0n155w', '0n140w', '0n125w', $
          '0n110w', '0n95w', '2n147e', '2n156e', '2n165e', '2n180w', '2n170w', '2n155w', '2n140w', '2n125w', '2n110w', $
          '2n95w', '5n147e', '5n156e', '5n170w', '5n155w', '5n140w', '5n125w', '5n110w', '5n95w', $
          '8n156e', '8n165e', '8n180w', '8n170w', '9n140w', '8n125w', '8n110w', '8n95w', $
          '0n0e', '0n10w', '0n23w', '0n35w', '10s10w', '12n23w', '12n38w', '14s32w', '15n38w', '19s34w', '20n38w', $
          '21n23w', '4n23w', '4n38w', '6s10w', '8n38w', '8s30w']

;;   This program will create the following text files with statistics of respective variables

close,/all

erase
;; PLOTTING THE MAPS

fi = project_od_env + 'flux_stat.txt'
res=read_ascii(fi,data_start=1)
ff=res.field1

lat=reform(ff(0,*))
lon=reform(ff(1,*))
cor=reform(ff(2,*))
bias=reform(ff(3,*))
std=reform(ff(4,*))
rmsd=reform(ff(5,*))

corr_t=total(cor)/n_elements(cor)
bias_t=total(bias)/n_elements(bias)
std_t=total(std)/n_elements(std)
rmsd_t=total(rmsd)/n_elements(rmsd)
cstat=string(corr_t, bias_t, std_t, rmsd_t, format='(f4.2,3x,f7.2,3x,f4.2,3x,f5.2)')

file= project_id_env + 'mask_oaflux_30N30S.nc'
initncdf, file
domdef, 30,390,-30,30
msk=ncdf_lec(file,var='msk')
marge1=[0,0,-5,0]

plt, msk,realcont=2,/nocont,/nofill,xminor=1,yminor=1,lct=64,cor_mi, cor_ma,int=cor_int, $
title='1) Correlation ', subtitle='', small=[1,4,1],/rempl,marge=marge1

NN=n_elements(lat)

usersym, [-.5,.5,.5,-.5,-.5],[-.5,-.5,.5,.5,-.5],/fill

for n=0,NN-1 do begin
    x=lon(n)
    y=lat(n)
    c=cor(n)
    cmi=cor_mi & cma=cor_ma & dc=cma-cmi
    col=((10+244*(c-cmi)/dc) > 10) < 254
    plots, x,y,psym=8,symsize=1.5,color=col
endfor

usersym, [-.5,.5,.5,-.5,-.5],[-.5,-.5,.5,.5,-.5]

for n=0,NN-1 do begin
    x=lon(n)
    y=lat(n)
    c=cor(n)
    plots, x,y,psym=8,symsize=1.5,color=0
endfor

plt,msk,realcont=2,/nocont,/nofill,xminor=1,yminor=1,lct=64,bias_mi, bias_ma, int=bias_int,/noer, marge=marge1, $
title='2) Mean bias', subtitle='', small=[1,4,2],/rempl

NN=n_elements(lat)

usersym, [-.5,.5,.5,-.5,-.5],[-.5,-.5,.5,.5,-.5],/fill

for n=0,NN-1 do begin
    x=lon(n)
    y=lat(n)
    c=bias(n)
    cmi=bias_mi & cma=bias_ma & dc=cma-cmi
    col=((10+244*(c-cmi)/dc) > 10) < 254
    plots, x,y,psym=8,symsize=1.5,color=col
endfor

usersym, [-.5,.5,.5,-.5,-.5],[-.5,-.5,.5,.5,-.5]

for n=0,NN-1 do begin
    x=lon(n)
    y=lat(n)
    c=bias(n)
    plots, x,y,psym=8,symsize=1.5,color=0
endfor
plt,msk,realcont=2,/nocont,/nofill,xminor=1,yminor=1,lct=64,rmsd_mi, rmsd_ma, int=rmsd_int,/noer, $
title='3) RMSD ', subtitle='', small=[1,4,3],/rempl, marge=marge1

NN=n_elements(lat)

usersym, [-.5,.5,.5,-.5,-.5],[-.5,-.5,.5,.5,-.5],/fill

for n=0,NN-1 do begin
    x=lon(n)
    y=lat(n)
    c=rmsd(n)
    cmi=rmsd_mi & cma=rmsd_ma & dc=cma-cmi
    col=((10+244*(c-cmi)/dc) > 10) < 254
    plots, x,y,psym=8,symsize=1.5,color=col
endfor

usersym, [-.5,.5,.5,-.5,-.5],[-.5,-.5,.5,.5,-.5]
for n=0,NN-1 do begin
    x=lon(n)
    y=lat(n)
    c=rmsd(n)
    plots, x,y,psym=8,symsize=1.5,color=0
endfor
plt,msk,realcont=2,/nocont,/nofill,xminor=1,yminor=1,lct=64,std_mi, std_ma, int=std_int,/noer, $
title='4) STD ratio', subtitle='', small=[1,4,4],/rempl, marge=marge1

NN=n_elements(lat)

usersym, [-.5,.5,.5,-.5,-.5],[-.5,-.5,.5,.5,-.5],/fill

for n=0,NN-1 do begin
    x=lon(n)
    y=lat(n)
    c=std(n)
    cmi=std_mi & cma=std_ma & dc=cma-cmi
    col=((10+244*(c-cmi)/dc) > 10) < 254
    plots, x,y,psym=8,symsize=1.5,color=col
endfor

usersym, [-.5,.5,.5,-.5,-.5],[-.5,-.5,.5,.5,-.5]

for n=0,NN-1 do begin
    x=lon(n)
    y=lat(n)
    c=std(n)
    plots, x,y,psym=8,symsize=1.5,color=0
endfor

erase
;; computing the pdf

x=mooring & y=product
mio=-50 & mao=300 & inx=25 & nx=(mao-mio)/inx+1l & xx=mio+indgen(nx)*inx
mis=-50 & mas=300 & iny=25 & ny=(mas-mis)/iny+1l & yy=mis+indgen(ny)*iny
pdf=fltarr(nx,ny)
for i=0, nx-2 do begin
    for j=0,ny-2 do begin
        xmi=xx(i) & xma=xx(i+1)
        ymi=yy(j) & yma=yy(j+1)
        ind=where((x ge xmi) and (x le xma) and (y ge ymi) and (y le yma))
        pdf(i,j)=float(n_elements(ind))
        ;    print, xmi,xma,ymi,yma,n_elements(ind)
    endfor
endfor

scontour, pdf,xx,yy ,small=[2,1,1],/noer, charsize=1, nlevels=30,/fill, $
xtitle='TAO Flux', ytitle='Product', title='PDF'
xyouts, 0, -150, cstat, charsize=1.2

closeps

end