source: trunk/src/flux_evaluation_tpr_timeseries.pro

;+
;
; .. _flux_evaluation_tpr_timeseries.pro:
;
; ==================================
; flux_evaluation_tpr_timeseries.pro
; ==================================
;
; DESCRIPTION
; ===========
;
; SEE ALSO
; ========
;
; :ref:`project_profile.sh`
; :ref:`project_init.pro`
; :ref:`cm_project.pro`
;
; TODO
; ====
;
; make it work : missing READ_TPR_NETFLUX
;
; EXAMPLES
; ========
;
; .. code-block:: idl
;
;    date1 = 20000101L
;    date2 = 20091231L
;    flux_evaluation_tpr_timeseries, date1, date2
;
; EVOLUTIONS
; ==========
;
; $Id$
;
; $URL$
;
; - fplod 20111130T141341Z 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_timeseries, $
;                    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 = 0
openps, FILENAME = project_od_env+'flux_evaluation_tpr_timeseries.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

;------------------------------------------------------------
;; Before running this program, you have to compile the following subroutines
;; .r extract_tpr_location
;; .r read_tpr_netflux

;; 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

fi=project_id_env + 'flux_stat.txt'
openw,1,fi

printf,1, 'x     y      cor    bias     std ratio     rmsd'
;------------------------------------------------------------------------------------------------------------------------
;; First, this program reads the full TropFlux data and later extract it at specific TPR locations

file = project_id_env + "shf_tropflux_1d_1989_2010.nc"
initncdf, file
var=-1*read_ncdf("shf", date1, date2, file=file,/nostr)
help, var

;------------------------------------------------------------------------------------------------------------------------
nn=n_elements(sitelist)
date1=date1
date2=date2
nsmooth=1. ;; this will return daily TPR values
mooring=0 & product=0

for n=0, nn-1 do begin

    ;; reading data from mooring

    site=sitelist(n) & csite=site
    print, csite
    x=x_site_location(site)
    y=y_site_location(site)
    if (y ge 0. and y le 30.) then y=y+360.
    dx=0.5 & dy=0.5 & box=[y-dy, y+dy, x-dx, x+dx]

    read_tpr_netflux, csite,date1,date2,nsmooth, $
            sw,lw,sh,lh

    ;; select the appropriate variables for evaluation (trp = sw or lw or sh or lh)

    tpr=sh

    ind=where(finite(tpr)) & no_valid=n_elements(ind)

    if (no_valid ge min_obs) then begin
        extract_tpr_location,var,box, $
        var_tpr
        var_tpr=reform(var_tpr)

        stats_5d, tpr,var_tpr, $  ;; tpr=TPR observation and var_tpr=gridded product extracted at TPR location
        cor, bias, std, rmsd

        printf, 1, x, y, cor, bias, std, rmsd, format='(f6.2, 3x, f6.2, 3x, f4.2, 3x, f7.2, 3x, f4.2, 3x, f5.2)'
        cstat=string(cor, bias, std, rmsd, format='(f4.2,3x,f7.2,3x,f4.2,3x,f5.2)')
        print, cstat

        ;;     PLOTTING THE TIME-SERIES
        array=[tpr, var_tpr] & mi=min(array,/nan) & ma=max(array,/nan) & int=(ma-mi)/3.
        pltt, ts_smooth(tpr,5,/nan), "t",/rempl, small=[1,3,1], lct=65, $
        title='Five day stats are shown below.  TPR (black) and Product (red) at'+csite+' ', charsize=1., $
        subtitle=cstat
        ind=where(finite(tpr,/nan)) & var_tpr(ind)=!Values.f_nan
        pltt, ts_smooth(var_tpr,5,/nan), "t",/ov1d, color=250
        erase
        mooring=[mooring,tpr] & product=[product,var_tpr]
    endif
endfor

close,/all
closeps

end