[Dart-dev] [4534] DART/trunk/diagnostics/threed_sphere: Finished the documentation for the rank histogram. Made a tr

nancy at ucar.edu nancy at ucar.edu
Mon Oct 18 16:13:37 MDT 2010


Revision: 4534
Author:   thoar
Date:     2010-10-18 16:13:37 -0600 (Mon, 18 Oct 2010)
Log Message:
-----------
Finished the documentation for the rank histogram. Made a trivial change
to the fortran code so one diagnostic comment will not wrap on output.

Modified Paths:
--------------
    DART/trunk/diagnostics/threed_sphere/obs_diag.f90
    DART/trunk/diagnostics/threed_sphere/obs_diag.html

-------------- next part --------------
Modified: DART/trunk/diagnostics/threed_sphere/obs_diag.f90
===================================================================
--- DART/trunk/diagnostics/threed_sphere/obs_diag.f90	2010-10-18 18:11:04 UTC (rev 4533)
+++ DART/trunk/diagnostics/threed_sphere/obs_diag.f90	2010-10-18 22:13:37 UTC (rev 4534)
@@ -3628,8 +3628,8 @@
                 dimids=(/ RegionDimID, LevelDimID, RankDimID, TimeDimID /), &
                 varid=VarID2), 'WriteTLRV', 'rank_hist:def_var')
          else
-            write(logfileunit,*)'Variable '//string1//' has ',ndata,'"rank"able observations.'
-            write(     *     ,*)'Variable '//string1//' has ',ndata,'"rank"able observations.'
+            write(logfileunit,*)string1//' has ',ndata,'"rank"able observations.'
+            write(     *     ,*)string1//' has ',ndata,'"rank"able observations.'
          endif
 
       endif

Modified: DART/trunk/diagnostics/threed_sphere/obs_diag.html
===================================================================
--- DART/trunk/diagnostics/threed_sphere/obs_diag.html	2010-10-18 18:11:04 UTC (rev 4533)
+++ DART/trunk/diagnostics/threed_sphere/obs_diag.html	2010-10-18 22:13:37 UTC (rev 4534)
@@ -81,11 +81,11 @@
    information to calculate a rank histogram, while the ensemble means are still
    used for all the other calculations.
    <table width="100%"><tr>
-   <td rowspan=2><img src="../../doc/html/obs_diag_evolution_example.png" width="300">
-   <td rowspan=2><img src="../../doc/html/obs_diag_profile_example.png" width="300">
-   <td><img src="../../doc/html/RankHistogram_ncview.png" width="200"></td>
+   <td rowspan=2><a href="../../doc/html/obs_diag_evolution_example.png"><img src="../../doc/html/obs_diag_evolution_example.png" width="300"></a>
+   <td rowspan=2><a href="../../doc/html/obs_diag_profile_example.png"><img src="../../doc/html/obs_diag_profile_example.png" width="300"></a>
+   <td><a href="../../doc/html/RankHistogram_ncview.png"><img src="../../doc/html/RankHistogram_ncview.png" width="200"></a></td>
    </tr>
-   <tr><td><img src="../../doc/html/RankHistogram_matlab.png" width="200"></td>
+   <tr><td><a href="../../doc/html/RankHistogram_matlab.png"><img src="../../doc/html/RankHistogram_matlab.png" width="200"></a></td>
    </tr>
    </table>
    <br />
@@ -111,19 +111,19 @@
    time-averaged vertical profile (figure in the middle), and sometimes 
    3) in terms of a rank histogram - "Where does the actual observation rank 
    relative to the  rest of the ensemble?" (figures on the right).
-   The rank histogram information can easily be plotted with <a href="http://meteora.ucsd.edu/~pierce/ncview_home_page.html">ncview</a>,
-   a free third-party piece of software; the rest of the
-   figures were created by several Matlab&#174; scripts that query 
-   the <em class="file">obs_diag_output.nc</em> file: 
-   <em class="file">DART/diagnostics/matlab/</em>
-   <a href=https://proxy.subversion.ucar.edu/DAReS/DART/trunk/diagnostics/matlab/plot_evolution.m>plot_evolution.m</a> ,
-   <a href=https://proxy.subversion.ucar.edu/DAReS/DART/trunk/diagnostics/matlab/plot_profile.m>plot_profile.m</a>, and 
-   <a href=https://proxy.subversion.ucar.edu/DAReS/DART/trunk/diagnostics/matlab/plot_rank_histogram.m>plot_rank_histogram.m</a>. Each of these takes as input a 
-   file name and a 'quantity' to plot ('rmse','spread','totalspread', ...)
+   The figures on the left and center were created by several Matlab&#174; 
+   scripts that query the <em class="file">obs_diag_output.nc</em> file: 
+   <em class="file">DART/diagnostics/matlab/<a href="https://proxy.subversion.ucar.edu/DAReS/DART/trunk/diagnostics/matlab/plot_evolution.m">plot_evolution.m</a></em> and 
+   <a href="https://proxy.subversion.ucar.edu/DAReS/DART/trunk/diagnostics/matlab/plot_profile.m">plot_profile.m</a>.
+   Both of these takes as input a 
+   file name and a 'quantity' to plot ('rmse','spread','totalspread',&nbsp;...)
    and exhaustively plots the quantity (for every variable, every level, 
    every region) in a single matlab figure window  - and creates a series 
    of .ps files with multiple pages for each of the figures. 
    The directory gets cluttered with them.
+   The rank histogram information can easily be plotted with <a href="http://meteora.ucsd.edu/~pierce/ncview_home_page.html">ncview</a>,
+   a free third-party piece of software or with 
+   <a href="https://proxy.subversion.ucar.edu/DAReS/DART/trunk/diagnostics/matlab/plot_rank_histogram.m">plot_rank_histogram.m</a>.
    <br />
    <br />
    <em class="program">obs_diag</em>
@@ -144,27 +144,13 @@
    <br />
    <br />
    <em class=program>obs_diag</em> reads <em class=file>obs_seq.final</em>
-   files and calculates the following quantities for an arbitrary number of
-   regions and levels:
+   files and calculates the following quantities (in no particular order) 
+   for an arbitrary number of regions and levels. It is necessary to query
+   the <em class=code>CopyMetaData</em> variable to determine the storage order
+   ( i.e. "which copy is what?").
 </P>
 
    <table width="90%">
-   <tr><td valign="top"><b>rmse</b></td>
-       <td>The root-mean-squared error (the horizontal wind components 
-           are also used to calculate the vector wind velocity 
-           and its RMS error).</td></tr>
-   <tr><td valign="top"><b>bias</b></td>
-       <td>The simple sum of forecast - observation. The bias of the 
-           horizontal wind speed (not velocity) is also computed.</td></tr>
-   <tr><td valign="top"><b>spread</b></td>
-       <td>The standard deviation of the univariate obs.
-           DART does not exploit the bivariate nature of U,V winds
-           and so the spread of the horizontal wind is defined as
-           the sum of the spreads of the U and V components.</td></tr>
-   <tr><td valign="top"><b>totalspread&nbsp;&nbsp;&nbsp;</b></td>
-       <td>The total standard deviation of the estimate. We pool the
-           ensemble variance of the observation plus the observation error 
-           variance and take the square root.</td></tr>
    <tr><td valign="top"><b>Nposs</b></td>
        <td>The number of observations available to be assimilated.</td></tr>
    <tr><td valign="top"><b>Nused</b></td>
@@ -182,14 +168,21 @@
        <td>the number of observations that were above or below the highest 
            or lowest model level, respectively;</td></tr> 
    <tr><td valign="top"><b>rmse</b></td>
-       <td>the rmse of the ensemble;</td></tr> 
-   <tr><td valign="top"><b>bias</b></td> 
-       <td>the bias of the ensemble (forecast-observation);</td></tr> 
+       <td>The root-mean-squared error (the horizontal wind components 
+           are also used to calculate the vector wind velocity 
+           and its RMS error).</td></tr>
+   <tr><td valign="top"><b>bias</b></td>
+       <td>The simple sum of forecast - observation. The bias of the 
+           horizontal wind speed (not velocity) is also computed.</td></tr>
    <tr><td valign="top"><b>spread</b></td>
-       <td>the spread of the ensemble; and the</td></tr> 
-   <tr><td valign="top"><b>totalspread</b></td>
-       <td>pooled spread of the observation (knowing its observational error) 
-           and the ensemble.</td></tr> 
+       <td>The standard deviation of the univariate obs.
+           DART does not exploit the bivariate nature of U,V winds
+           and so the spread of the horizontal wind is defined as
+           the sum of the spreads of the U and V components.</td></tr>
+   <tr><td valign="top"><b>totalspread&nbsp;&nbsp;&nbsp;</b></td>
+       <td>The total standard deviation of the estimate. We pool the
+           ensemble variance of the observation plus the observation error 
+           variance and take the square root.</td></tr>
    <tr><td valign="top"><b>NbadDARTQC&nbsp;&nbsp;&nbsp;</b></td>
        <td>the number of observations that had a DART QC value 
            (&gt; 1 for a prior, &gt; 3 for a posterior)</td></tr> 
@@ -265,8 +258,6 @@
    namelist variables (the counts of these are available in the 
    netCDF file as <b>NbadQC</b> and <b>NbadIZ</b>, respectively).</LI>
 
-   <LI>removed the <em class="removed">obs_select</em> namelist variable,</LI>
-
 </OL>
 
 <P>
@@ -318,19 +309,20 @@
 <em class=call>namelist / obs_diag_nml / </em> 
       obs_sequence_name, obs_sequence_list, first_bin_center, last_bin_center, 
       bin_separation, bin_width, time_to_skip, max_num_bins, 
-      plevel, hlevel, mlevel, <em class=removed>obs_select</em>, 
+      plevel, hlevel, mlevel, 
       rat_cri, input_qc_threshold, 
       Nregions, lonlim1, lonlim2, latlim1, latlim2, 
-      reg_names, print_mismatched_locs, print_obs_locations, verbose
+      reg_names, print_mismatched_locs, print_obs_locations, verbose,
+      outliers_in_histogram
 </pre>
 </div>
 
 <div class=indent1>
 <!-- Description -->
 
-<P>This namelist is read in a file called <em class=file>input.nml</em>. <BR>
+<P>This namelist is read in a file called <em class=file>input.nml</em>. <br />
    The date-time integer arrays in this namelist have the form 
-   (YYYY,&nbsp;MO,&nbsp;DA,&nbsp;HR,&nbsp;MIN,&nbsp;SEC). <BR>
+   (YYYY,&nbsp;MO,&nbsp;DA,&nbsp;HR,&nbsp;MIN,&nbsp;SEC). <br />
    The allowable ranges for the region boundaries are: latitude [-90.,90], 
    longitude [0.,Inf.]
 <br />
@@ -349,7 +341,7 @@
 
 <TR><!--contents--><TD valign=top>   obs_sequence_name </TD>
     <!--  type  --><TD valign=top>   character         </TD>
-    <!--descript--><TD>Name of the observation sequence file(s). <BR>
+    <!--descript--><TD>Name of the observation sequence file(s). <br />
             This may be a relative or absolute filename. If the filename
             contains a '/', the filename is considered to be comprised of
 	    everything to the right, and a directory structure to the left.
@@ -359,21 +351,22 @@
 	    is to look for additional files to include until the files 
 	    are exhausted or an <em class=file>obs_seq.final</em> file 
 	    is found that contains observations beyond the timeframe of 
-	    interest.<BR>e.g. 'obsdir_001/obs_seq.final' will cause
+	    interest.<br />e.g. 'obsdir_001/obs_seq.final' will cause
 	    <em class="program">obs_diag</em> to look for
-	    'obsdir_002/obs_seq.final', and so on.<BR>
-            If this is set, 'obs_sequence_list' must be set to ''.
+	    'obsdir_002/obs_seq.final', and so on.<br />
+            If this is set, <em class=code>obs_sequence_list</em> must be set to ''.
+            <br />
 	    Default 'obs_seq.final'</TD></TR>  
 
 <TR><!--contents--><TD valign=top>   obs_sequence_list </TD>
     <!--  type  --><TD valign=top>   character         </TD>
     <!--descript--><TD>Name of an ascii text file which contains a list
             of one or more observation sequence files, one per line.
-            If this is specified, 'obs_sequence_name' must be set to ''. 
+            If this is specified, <em class=code>obs_sequence_name</em> must be set to ''. 
             Can be created by any method, including sending the output
             of the 'ls' command to a file, a text editor, or another program.
-            If this is set, 'obs_sequence_name' must be set to ''.
-	    <BR>
+            If this is set, <em class=code>obs_sequence_name</em> must be set to ''.
+	    <br />
 	    Default '' - an empty string.</TD></TR>  
 
 <TR><!--contents--><TD valign=top>   first_bin_center </TD>
@@ -385,7 +378,7 @@
 	    first and last observations in every observation sequence file. 
 	    Look for the string 'First observation date' in the logfile. 
 	    If the <em class="code">verbose</em> is 'true', it is also 
-	    written to the screen.<BR>
+	    written to the screen.<br />
             Default: 2003, 1, 1, 0, 0, 0</TD></TR>  
 
 <TR><!--contents--><TD valign=top>   last_bin_center  </TD>
@@ -398,13 +391,13 @@
 	    <em class=code>bin_separation</em> are used to populate the time 
 	    array and stop on or before the time defined by 
 	    <em class=code>last_bin_center</em>. 
-	    See also <em class=code>max_num_bins</em>.<BR>
+	    See also <em class=code>max_num_bins</em>.<br />
             Default: 2003, 1, 2, 0, 0, 0</TD></TR>  
 
 <TR><!--contents--><TD valign=top>   bin_separation       </TD>
     <!--  type  --><TD valign=top>  integer, dimension(6) </TD>
     <!--descript--><TD>Time between bin centers.
-	    The year and month values <em>must</em> be zero.<BR>
+	    The year and month values <em>must</em> be zero.<br />
             Default: 0, 0, 0, 6, 0, 0</TD></TR>  
 
 <TR><!--contents--><TD valign=top>   bin_width        </TD>
@@ -412,7 +405,7 @@
     <!--descript--><TD>Time span around bin centers in which obs will be 
 	    compared.  The year and month values <em>must</em> be zero. 
 	    Frequently, but not required to be, the same as the values 
-	    for bin_separation.<BR>
+	    for bin_separation.<br />
             Default: 0, 0, 0, 6, 0, 0</TD></TR>  
 
 <TR><!--contents--><TD valign=top>   time_to_skip     </TD>
@@ -422,7 +415,7 @@
 	    values <em>must</em> be zero. Useful because it takes some 
 	    time for the assimilation to settle down from the climatological 
 	    spread at the start. <em class=code>time_to_skip</em> is an
-            amount of time AFTER the first edge of the first bin.<BR>
+            amount of time AFTER the first edge of the first bin.<br />
             Default: 0, 0, 1, 0, 0, 0</TD></TR>  
 
 <TR><!--contents--><TD valign=top>   max_num_bins  </TD>
@@ -432,7 +425,7 @@
 	    If <em class=code>max_num_bins</em> is set to '10', only
 	    10 timesteps will be output - provided 
 	    <em class=code>last_bin_center</em> is set to some later date. 
-	    <BR>
+	    <br />
             Default: 1000000</TD></TR>  
 
 <TR><!--contents--><TD valign=top>   plevel     </TD>
@@ -442,36 +435,29 @@
 	    continuum is used. If a single midpoint value is entered, the 
 	    bin edges are +/- 10% of the midpoint value. If you'd like to 
 	    change that ... see the routine 
-	    <em class="routine">Rmidpoints2edges()</em>.<BR>
+	    <em class="routine">Rmidpoints2edges()</em>.<br />
             Default: (1000, 925, 850, 700, 500, 400, 300, 
 	    250, 200, 150, 100)</TD></TR>  
 
 <TR><!--contents--><TD valign=top>   hlevel     </TD>
     <!--  type  --><TD valign=top> integer, <em class="new">dimension(50)</em></TD>
     <!--descript--><TD>Same, but for observations that have height(m) as the 
-	    vertical coordinate.<BR>
+	    vertical coordinate.<br />
             Default: (1000, 2000, 3000, 4000, 5000, 6000,
 	    7000, 8000, 9000, 10000, 11000)</TD></TR>  
 
 <TR><!--contents--><TD valign=top>   mlevel     </TD>
     <!--  type  --><TD valign=top> integer, <em class="new">dimension(50)</em></TD>
     <!--descript--><TD>Same, but for observations on model level 
-            (good luck getting them).<BR>
+            (good luck getting them).<br />
             Default: (1, 2, 3, 4, 5, 6, 7, 8, 9, 10)</TD></TR>  
 
-<TR><!--contents--><TD valign=top><em class="removed">obs_select </em></TD>
-    <!--  type  --><TD valign=top>   integer    </TD>
-    <!--descript--><TD><em class="removed">Which obs subset to use; all(1), 
-            radiosonde only(2), everything but radiosondes(3).</em> 
-	    Removed because all observations are independently processed.<BR>
-            </TD></TR>  
-
 <TR><!--contents--><TD valign=top>   rat_cri    </TD>
     <!--  type  --><TD valign=top>   real       </TD>
     <!--descript--><TD>Ratio of 'distance' of the ensemble mean to the observation 
 	    in terms of the standard deviation of the ensemble and
-	    observational error.<br>
-	    <em class="code">abs(prior_mean-obsval) / sqrt(prior_spread+obs_error_variance)</em><br>
+	    observational error.<br />
+	    <em class="code">abs(prior_mean-obsval) / sqrt(prior_spread+obs_error_variance)</em><br />
 	    Since we generally don't have all the ensemble members (in the
 	    <em class="file">obs_seq.final</em> file) we cannot
 	    calculate the covariance term that should be in the denominator, 
@@ -481,25 +467,25 @@
 	    simply be counted in <b>NbadIZ</b>. 
 	    A histogram of the ratios is saved in <em
 	    class="file">LargeInnov.txt</em>, as are individual observations
-	    that have a large ratio (i.e. the outliers)<BR>
+	    that have a large ratio (i.e. the outliers)<br />
             Default: 5000.0</TD></TR>  
 
 <TR><!--contents--><TD valign=top>   input_qc_threshold </TD>
     <!--  type  --><TD valign=top>   real         </TD>
     <!--descript--><TD>Observations with quality control values greater than this
 	    <em class="removed">will be excluded from the diagnostics</em>
-	    will be counted in <b>NbadQC</b>.<BR>
+	    will be counted in <b>NbadQC</b>.<br />
             Default: 4.0</TD></TR>  
 
 <TR><!--contents--><TD valign=top>   Nregions   </TD>
     <!--  type  --><TD valign=top>   integer    </TD>
     <!--descript--><TD>Number of regions of the globe for which you'd like
-	    obs space diagnostics. Must be between [1,50].<BR>
+	    obs space diagnostics. Must be between [1,50].<br />
 	    Default: 4</TD></TR>  
 
 <TR><!--contents--><TD valign=top>   lonlim1    </TD>
     <!--  type  --><TD valign=top>   real, <em class="new">dimension(50)</em></TD>
-    <!--descript--><TD>Westernmost longitudes of each of the regions.<BR>
+    <!--descript--><TD>Westernmost longitudes of each of the regions.<br />
 	    Default: 0.0, 0.0, 0.0, 235.0</TD></TR>  
 
 <TR><!--contents--><TD valign=top>   lonlim2    </TD>
@@ -509,29 +495,29 @@
 	    the westernmost values, it defines a region that spans the 
 	    prime meridian. No problem.</em> It is perfectly acceptable 
 	    to specify lonlim1 = 330 , lonlim2 = 50 to identify a region
-	    like "Africa".<BR>
+	    like "Africa".<br />
 	    Default: 360.0, 360.0, 360.0, 295.0</TD></TR>  
 
 <TR><!--contents--><TD valign=top>   latlim1    </TD>
     <!--  type  --><TD valign=top>   real, <em class="new">dimension(50)</em></TD>
-    <!--descript--><TD>Southernmost latitudes of the regions.<BR>
+    <!--descript--><TD>Southernmost latitudes of the regions.<br />
 	    Default: 20.0, -80.0, -20.0, 25.0</TD></TR>  
 
 <TR><!--contents--><TD valign=top>   latlim2    </TD>
     <!--  type  --><TD valign=top>   real, <em class="new">dimension(50)</em></TD>
-    <!--descript--><TD>Northernmost latitudes of the regions.<BR>
+    <!--descript--><TD>Northernmost latitudes of the regions.<br />
 	    Default: 80.0, -20.0, 20.0, 55.0</TD></TR>  
 
 <TR><!--contents--><TD valign=top>   reg_names  </TD>
     <!--  type  --><TD valign=top>   character(len=129), <em class="new">dimension(50)</em>  </TD>
-    <!--descript--><TD>Array of names for the regions to be analyzed.<BR>
+    <!--descript--><TD>Array of names for the regions to be analyzed.<br />
 	    Default: 'Northern&nbsp;Hemisphere','Southern&nbsp;Hemisphere',
 	    'Tropics','North&nbsp;America'</TD></TR>  
 
 <TR><!--contents--><TD valign=top>   print_mismatched_locs </TD>
     <!--  type  --><TD valign=top>   logical               </TD>
     <!--descript--><TD>Print instances where U and V "pairs" don't have the 
-	    same location.<BR>
+	    same location.<br />
 	    Default: .false. </TD></TR>  
 
 <TR><!--contents--><TD valign=top>   print_obs_locations </TD>
@@ -545,18 +531,25 @@
 	    There is one output file per temporal bin.
 	    For pathological reasons, these files cannot exist prior to execution.
 	    If they do, an error will result, and execution halts..
-	    <BR>
+	    <br />
 	    Default: .false. </TD></TR>  
 
 <TR><!--contents--><TD valign=top>   verbose   </TD>
     <!--  type  --><TD valign=top>   logical   </TD>
-    <!--descript--><TD>Print extra info about the obs_diag run.<BR>
+    <!--descript--><TD>Print extra info about the obs_diag run.<br />
             Default: .false. </TD></TR>  
 
+<TR><!--contents--><TD valign=top>outliers_in_histogram</TD>
+    <!--  type  --><TD valign=top>   logical   </TD>
+    <!--descript--><TD>Determines whether or not to use observations that 
+            have been rejected by the outlier threshhold mechanism in the 
+            calculation of the rank histogram.<br />
+            Default: .false. </TD></TR>  
+
 </TABLE>
 
 </div>
-<br>
+<br />
 
 <!--==================================================================-->
 
@@ -610,7 +603,7 @@
 It all depends on what's needed by the <em class="file">model_mod.f90</em> 
 </P>
 
-<H3 class=indent1>Discussion of obs_diag_output.nc variables</H3>
+<H3 class=indent1>Discussion of obs_diag_output.nc</H3>
 
 <P>
 Every observation type encountered in the observation sequence file is 
@@ -667,14 +660,17 @@
           RADIOSONDE_V_WIND_COMPONENT_analy:missing_value = -888888.f ;
   ...
 </pre>
+
 <P>
-There are several things to note:<br>
-1) the 'WIND_VELOCITY' component is nowhere 'near' the corresponding U,V components.<br>
-2) all of the 'guess' variables come before the matching 'analy' variables.<br>
-3) surface variables (i.e. <tt>MARINE_SFC_ALTIMETER</tt> have a
+There are several things to note:
+</P>
+<ol>
+<li>the 'WIND_VELOCITY' component is nowhere 'near' the corresponding U,V components.</li>
+<li>all of the 'guess' variables come before the matching 'analy' variables.</li>
+<li>surface variables (i.e. <tt>MARINE_SFC_ALTIMETER</tt> have a
 coordinate called 'surface' as opposed to 'plevel' for the others in
-this example).
-</P>
+this example).</li>
+</ol>
 
 <H4 class=indent1>vertical profiles</H4>
 <P>
@@ -689,6 +685,7 @@
   float SAT_WIND_VELOCITY_VPguess(copy, plevel, region) ;
           SAT_WIND_VELOCITY_VPguess:_FillValue = -888888.f ;
           SAT_WIND_VELOCITY_VPguess:missing_value = -888888.f ;
+  ...
   float RADIOSONDE_U_WIND_COMPONENT_VPanaly(copy, plevel, region) ;
           RADIOSONDE_U_WIND_COMPONENT_VPanaly:_FillValue = -888888.f ;
           RADIOSONDE_U_WIND_COMPONENT_VPanaly:missing_value = -888888.f ;
@@ -701,6 +698,41 @@
 participate. There's a difference!
 </P>
 
+<H4 class=indent1>rank histograms</H4>
+<P>
+If it is possible to calculate a rank histogram, there will also be :</P>
+<pre>
+   ...
+   int RADIOSONDE_U_WIND_COMPONENT_guess_RankHi(time, rank_bins, plevel, region) ;
+   ...
+   int RADIOSONDE_V_WIND_COMPONENT_guess_RankHi(time, rank_bins, plevel, region) ;
+   ...
+   int MARINE_SFC_ALTIMETER_guess_RankHist(time, rank_bins, surface, region) ;
+   ...
+</pre>
+<P>
+as well as some global attributes. The attributes reflect the namelist settings
+and can be used by plotting routines to provide additional annotation for 
+the histogram. 
+</P>
+<pre>
+		:DART_QCs_in_histogram = 0, 1, 2, 3, 7 ;
+		:outliers_in_histogram = "TRUE" ;
+</pre>
+<P>
+Please note:
+</P>
+<ol>
+<li>netCDF restricts variable names to 40 characters, so '_Rank_Hist' may
+    be truncated.</li>
+<li>It is sufficiently vague to try to calculate a rank histogram for a 
+    velocity derived from the assimilation of U,V components such that 
+    NO rank histogram is created for velocity.  A run-time log message will 
+    inform as to which variables are NOT having a rank histogram variable 
+    preserved in the <em class=file>obs_diag_output.nc</em> file -
+    IFF it is possible to calculate a rank histogram in the first place.</li>
+</ol>
+
 <!--==================================================================-->
 <!-- Discuss  typical usage of obs_diag.                              -->
 <!--==================================================================-->
@@ -710,19 +742,19 @@
 <H2>USAGE</H2>
 
 <P>
-Obs_diag is built in .../DART/models/your_model/work, in the same way
+Obs_diag is built in .../DART/models/<em class=italic>your_model</em>/work, in the same way
 as the other DART components.
 </P>
 
 <H3 class=indent1>multiple observation sequence files</H3>
 <P>
-Perhaps the most user-friendly enhancement is the strategy for processing
-a series of observation sequence files. Everything keys off the 
-<em class=code>obs_sequence_name</em> namelist variable. Several examples
-might be the most effective description.
+There are two ways to specify input files for <em class=program>obs_diag</em>.
+You can either specify the name of a file containing a list of files (in
+<em class=code>obs_sequence_list</em>), or you may specify 
+<em class=code>obs_sequence_name</em>, which is descibed by the following illustration&nbsp;...
 </P>
 <table width="100%">
-<TR><TH width=50%>value</TH><TH>effect</TH></TR>
+<TR><TH width=50% align=left>value</TH><TH>effect</TH></TR>
 
 <TR><TD valign="top">"obs_seq.final"</TD>
     <TD>A single file is processed.
@@ -751,43 +783,113 @@
 	   of files to ingest, or if it encounters an observation sequence
 	   file whose first time is beyond the timeframe of interest.</TD></TR>
 </table>
+
+
+
+<H3 class=indent1>Example: creating a rank histogram from observation sequence files 
+spanning&nbsp;30&nbsp;days.</H3>
+<table width=90%>
+<tr><td>In this example, we are generating an <em class=file>obs_diag_output.nc</em> 
+   file that will have exactly ONE timestep in it - for the entire globe - and we 
+   will NOT be using outlier observations in the rank histogram. Lets presume 
+   that all your <em class=file>obs_seq.final</em> files are in alphabetically-nice 
+   directories:</td>
+   <td>&nbsp;&nbsp;<img src="../../doc/html/RankHistogram_matlab.png" width="200"></a></td>
+</tr>
+</table>
+<pre>
+/Exp1/Dir01/obs_seq.final
+/Exp1/Dir02/obs_seq.final
+/Exp1/Dir03/obs_seq.final
+...
+/Exp1/Dir99/obs_seq.final
+</pre>
+<P>The first step is to create a file containing the list of observation sequence
+files you want to use. This can be done with the unix command 'ls' with 
+the -1 option (that's a number one) to put one file per line.
+</P>
+<div class=unix>
+ls -1 /Exp1/Dir*/obs_seq.final &gt; obs_file_list.txt
+</div>
+<P>It is necessary to turn on the verbose option to 
+check the first/last times that will be used for the histogram.
+Then, the namelist settings for 2008 07 31 12Z through 2008 08 30 12Z are:
+</P>
+<div class="routine">
+<pre>
+&amp;obs_diag_nml
+   obs_sequence_name  = '',
+   obs_sequence_list  = <em class="input">'obs_file_list.txt',</em>
+   first_bin_center   = <em class="input"> 2008, 8,15,12, 0, 0 ,</em>
+   last_bin_center    = <em class="input"> 2008, 8,15,12, 0, 0 ,</em>
+   bin_separation     = <em class="input">    0, 0,30, 0, 0, 0 ,</em>
+   bin_width          = <em class="input">    0, 0,30, 0, 0, 0 ,</em>
+   time_to_skip       = <em class="input">    0, 0, 0, 0, 0, 0 ,</em>
+   max_num_bins       = 1000,
+   rat_cri            = 5000.0,
+   input_qc_threshold = 3.0,
+   Nregions   = 1,
+   lonlim1    = <em class="input">  0.0,</em>
+   lonlim2    = <em class="input">360.0,</em>
+   latlim1    = <em class="input">-90.0,</em>
+   latlim2    = <em class="input"> 90.0,</em>
+   reg_names  = <em class="input">'Entire Domain',</em>
+   print_mismatched_locs = .false.,
+   print_obs_locations   = .false.,
+   verbose               = <em class="input">.true.,</em>
+   outliers_in_histogram = <em class="input">.false.</em>,
+   /
+</pre>
+</div>
+<P>
+then, simply run <em class=program>obs_diag</em> in the usual manner - you may want
+to save the run-time output to a file. Here is a portion of the run-time output:
+</P>
+<div class="unix">
+<pre>
+...
+Region  1 Entire Domain                    (WESN):     0.0000   360.0000   -90.0000    90.0000
+ Requesting            1  assimilation periods.
+
+epoch      1  start day=148865, sec=43201
+epoch      1 center day=148880, sec=43200
+epoch      1    end day=148895, sec=43200
+epoch      1  start 2008 Jul 31 12:00:01
+epoch      1 center 2008 Aug 15 12:00:00
+epoch      1    end 2008 Aug 30 12:00:00
+...
+MARINE_SFC_HORIZONTAL_WIND_guess_RankHis has            0 "rank"able observations. 
+SAT_HORIZONTAL_WIND_guess_RankHist       has            0 "rank"able observations.
+...
+</pre>
+</div>
+<P>
+Discussion: It should be pretty clear that there is exactly 1 assimilation period,
+it may surprise you that the start is 1 second past 12Z. This is deliberate and reflects
+the DART convention of starting intervals 1 second after the end of the previous interval. 
+The times in the netCDF variables reflect the defined start/stop of the period, regardless
+of the time of the first/last observation.
+<br />
+<br />
+Please note that none of the 'horizontal_wind' variables will have a rank histogram,
+so they are not written to the netCDF file. ANY variable that does not have a rank 
+histogram with some observations will NOT have a rank histogram variable 
+in the netCDF file.
+<br />
+<br />
+Now that you have the <em class=file>obs_diag_output.nc</em>, you can explore it with
+<a href="http://meteora.ucsd.edu/~pierce/ncview_home_page.html">ncview</a> or
+<a href="https://proxy.subversion.ucar.edu/DAReS/DART/trunk/diagnostics/matlab/plot_rank_histogram.m">plot_rank_histogram.m</a>.
+</P>
+
   
-<H3 class=indent1>plotting locations</H3>
+<H3 class=indent1><em class=removed>Example: plotting locations</em> deprecated</H3>
 <P>
-<em class="program">obs_diag</em> does not plot your locations.
-Since <em class="program">obs_diag</em> already reads observation
-sequences, it is a tiny extension to write out files that contain 
-rudimentary location information. 
-If you set the value of the namelist parameter <em class="code">
-print_obs_locations</em> to <tt>.true.</tt>, you can create
-output files (of the form
-<em class="file">observation_locations.xxx.dat</em>) which may then be plotted 
-with just about anything ... the first line of the file (they are 
-plain text, despite the .dat suffix) explains it all. One of the 
-columns in the output files is the observation type - which can 
-be decoded by using the key in <em class="file">obs_diag_output.nc</em> -
-specifically, try the unix command:
+<strong>please convert your observation sequence file to netCDF format with <a href="https://proxy.subversion.ucar.edu/DAReS/DART/trunk/obs_sequence/obs_seq_to_netcdf.html">obs_seq_to_netcdf</a></strong> and use 
+<a href="https://proxy.subversion.ucar.edu/DAReS/DART/trunk/diagnostics/matlab/plot_obs_netcdf.m">plot_obs_netcdf.m</a></strong>.
 </P>
 
-<div class="unix">ncdump -v ObservationTypes obs_diag_output.nc</div>
-<BR>
 
-<table width="100%">
-<TR><TD>
-<img src="../../doc/html/obs_diag_location_example.png" width="500"
-     alt="map of the world with altimeter observation locations
-	plotted"></TD>
-<TD>The Matlab&#174; script <em class="file">plot_observation_locations.m</em>
-automatically decodes the observation types and creates a legend. 
-It also allows plotting of observations with select QC values. 
-So, for example, you can plot the location of all the observations 
-that WERE NOT assimilated successfully.<BR>
-Sometimes it is useful to plot observation sequence files that have
-not yet been assimilated, so they have no ensemble estimates. The
-rest of <em class="program">obs_diag</em> is pretty useless, but 
-you can still output the observation locations.</TD></TR>
-</table>
-
 <!--==================================================================-->
 <!-- Cite references, if need be.                                     -->
 <!--==================================================================-->
@@ -884,9 +986,9 @@
 <H2>Terms of Use</H2>
 
 <P>
-DART software - Copyright &#169; 2004 - 2010 UCAR.<br>
-This open source software is provided by UCAR, "as is",<br>
-without charge, subject to all terms of use at<br>
+DART software - Copyright &#169; 2004 - 2010 UCAR.<br />
+This open source software is provided by UCAR, "as is",<br />
+without charge, subject to all terms of use at<br />
 <a href="http://www.image.ucar.edu/DAReS/DART/DART_download">
 http://www.image.ucar.edu/DAReS/DART/DART_download</a>
 </P>


More information about the Dart-dev mailing list