NOAA-PMEL/Ferret

Ensuring that a given axis is always used in a consistent grid position

Closed this issue · 5 comments

Reported by @AnsleyManke on 15 Dec 2010 23:44 UTC
Some datasets have dimensions or coordinate variables that don't explicitly tell Ferret their direction. In such cases the flag line_direction(line) gets the value of 'NA' and then when the grids are created, the axis might take on different directions for variables in the dataset. This can cause troubles in f-tds.

Here's a short cdl file, where the axes bbb and ccc have no defined direction. The behavior is similar if the dimensions do not have corresponding coordinate variables - though we need to work out the behavior in both cases.

netcdf ambiguous {
dimensions:
	bbb = 3 ;
	ccc = 5 ;
variables:
	float E(ccc, bbb) ;
		E:missing_value = -1.e+34f ;
		E:_FillValue = -1.e+34f ;
		E:long_name = "2D variable (ccc,bbb)" ;
	float F(bbb) ;
		F:missing_value = -1.e+34f ;
		F:_FillValue = -1.e+34f ;
		F:long_name = "1D in bbb dimension" ;
	float G(ccc) ;
		G:missing_value = -1.e+34f ;
		G:_FillValue = -1.e+34f ;
		G:long_name = "1D in ccc dimension" ;
	double bbb(bbb);
		bbb:point_spacing = "even" ;
		bbb:point_spacing = "even" ;
	double ccc(ccc);
		ccc:point_spacing = "even" ;


// global attributes:
		:history = "FERRET V6.65   15-Dec-10" ;
		:Conventions = "CF-1.0" ;
data:

 bbb =
  1, 2, 3 ;

 ccc =
  1, 2, 3, 4, 5 ;

 E =
  2, 3, 4,
  3, 4, 5,
  4, 5, 6,
  5, 6, 7,
  6, 7, 8 ;

 F = 2, 4, 6 ;

 G = 1, 2, 3, 4, 5 ;

}

Opening the corresponding nc file, the SHOW DATA has the ccc dimension being used as the X axis in some variables, and the Y axis in another.

yes? use ambiguous.nc 
yes? sh dat
     currently SET data sets:
    1> ./ambiguous.nc  (default)
 name     title                             I         J         K         L
 E        2D variable (ccc,bbb)            1:3       1:5       ...       ...
 F        1D in bbb dimension              1:3       ...       ...       ...
 G     
 



The approach advocated below is based on the following:

   1. line (aka axis) orientations of "NA" are OK.  A non-"NA" line orientation should be assigned (in cd_get_1_axis) only when there is firm evidence to support the conclusion
   2. if a netCDF file is interpreted correctly then any given axis should be found in the same idim position in every grid
   3. we will create new code to resolve simple ambiguities in the placement of "NA" axes into grids.  We can correct only problems that can be fixed by "sliding" an axis to one side or another in the grid.
          * this is the case in Rich's file where we have two "NA" axes, "boundary" and "tracer".  The ambiguity of having (boundary, tracer) in grid1 and (tracer) alone in grid2 can be resolved by sliding "tracer" over to slot 2 in grid2 

   1. in cd_get_1_axis when determining the orientations of axes (setting of line_direction)
         1. test appropriate standard_name strings (e.g. for sigma levels)
         2. remove weak inferences of orientation
                * remove the block of code below:  "* If still undefined, 
                  try the starting letter."  (though shall we still interpret axes 
                  named with the single letters X, Y, Z, T as axes in those 
                  directions?)
                * If you find other weak inferences of orientation lets talk about
                  them.  

   2. in cd_init_dset just after cd_get_parent_grid.F has been called, call this
      new routine: CD_CONSISTENT_AXIS_ORIENTATION
         1. scan the grids just created for this dataset to generate a list of any 
            axes that occupy different idim positions in different grids.
         2. if multiply-oriented axes are found
                * (START)  if the axis is of known orientation issue a warning and 
                  do not try to fix
                * if the axis is of "NA" orientation then try this simple (tho not 
                  fool-proof) algorithm:
               1. choose the first grid where the axis is found as grid1
               2. let pos1 be the position of the axis in grid1
               3. set adjustment_success flag to TRUE
               4. for each other grid (always including the full list of this 
                  dataset's grids )
                      * if grid has this axis in a different position then
                            o let pos2 be the position of the axis in grid2
                            o in grid2 if pos1 and all positions between pos1 and 
                              pos2 are unoccupied
                              (** or are occupied by other "NA" axes)
                                  + then
                                        # shift the axis in question from pos2 to 
                                          pos1 
                                  + else
                                        # set adjustment_success flag to FALSE
               5. if adjustment_success is FALSE, find the next grid where the axis
                  is found and go back to #2.  If no next grid exists, give up on 
                  adjustments and issue a warning about ambiguous orientations of 
                  this axis.
               6. if adjustment_success is TRUE, move to the next multiply-oriented
                  axis and go back to START.  IF no more, then we're done.

** Note above where it says "** or are occupied by other "NA" axes" that it can be desirable to shift other "NA" axes.   For example if grid1 has axes (na1, na2, na3)  and grid2 has (na2, na3) it would be desirable to shift both na2 and na3 axes to the right by one idim position in grid2.  The simple algorithm above needs to be enhanced to allow this. Suppress the warnings, and repeat the entire algorithm from #1 above START any time you shift an axis other than the current multi-oriented axis.  And then you have to include a trap to make sure that you don't repeat the process more than (say) 4 times.  Else (I suspect) you could become caught in an infinite loop.   

The above neglects the possibility of a permuted grid. Probably implement non-permuted grids, first and then look at the possibility of permuted ones. The logic of "sliding an axis to one side or the other in the grid" needs to be considered in the context of the permutation.

The bad news is that the permutation logic is, itself, not clean.  (The routine tm_axis_order.F has the gory details.  The good news is that permutation is rarely anything other than the default 1-2-3-4).   idim slots must be neighbors in the sense of the permutation array to allow shifting.   This is a nasty problem.  Key to understanding it will be having some concrete examples to work from.  For example

example file #1:

    * var1(tracer, lon,lat)
    * var2(tracer)

example file #2:

    * var1(tracer,lat,lon)
    * var2(tracer)




Migrated-From: http://dunkel.pmel.noaa.gov/trac/ferret/ticket/1774

Comment by @AnsleyManke on 20 Jan 2011 21:27 UTC
I have this working. The grids for the dataset are checked, and axis positions moved if they are inconsistent. The issue of permutations needs more testing; what I have done is to make axis positions consistent and then apply the user-supplied ordering. This works for the cases I've tried, but as Steve notes it's not entirely clear how that should be done.

Comment by @AnsleyManke on 28 Jan 2011 01:03 UTC
Fixed in the ferret_v67beta executable installed on Dunkel and Stout. The issue of permutations needs further testing but I think it works consistently.

Comment by @AnsleyManke on 29 Feb 2012 00:26 UTC
Reopening. Reported by Josefine Ghattas, a colleague of Patrick Brockmann.

When an axis in the netCDF file is reversed ( here a Y-axis listed with coordinates as follows: rlatu = 1.55, 1.52, 1.48, ..., -1.48, -1.52, -1.55 ; )

and two datasets that use the same grids are opened, the reversed-character of the axis is not applied to the second dataset. The information is lost when doing the consistent-axis-ordering fixes here. The dataset has a dimension with ambiguous direction, SIGS, so the code to figure out a consistent axis ordering is called.

I created a subset of the datset, see short_data.nc attached to this ticket. Here is an example:

yes? use short_data.nc 
yes? sh dat
     currently SET data sets:
    1> ./short_data.nc  (default)
 name     title                             I         J         K         L
 UCOV     UCOV                             1:6       1:9       1:3       1:1
 CONTROLE CONTROLE                         1:6       ...       ...       ...
 NIVSIGS  NIVSIGS                          ...       ...       1:3       ...
 NIVSIG   NIVSIG                           1:6       ...       ...       ...
 AP       AP                               1:6       ...       ...       ...
 
yes? sp cp short_data.nc short_copy.nc
yes? use short_copy.ncyes? list/k=3/i=1 ucov[d=1]
             VARIABLE : UCOV
             FILENAME : short_data.nc
             SUBSET   : 9 points (Y)
             X        : 1
             Z        : 3
             TIME     : 08-JAN-1980 00:00
                1.014  
                 1
 -1.554 / 1: -642176.
 -1.521 / 2: -621550.
 -1.488 / 3: -700065.
 -1.455 / 4: -640073.
 0.033  / 5: -375011.
 1.455  / 6: -153460.
 1.488  / 7:   10663.
 1.521  / 8:   60693.
 1.554  / 9:       0.
yes? list/k=3/i=1 ucov[d=2]
             VARIABLE : UCOV
             FILENAME : short_copy.nc
             SUBSET   : 9 points (Y)
             X        : 1
             Z        : 3
             TIME     : 08-JAN-1980 00:00
                1.014  
                 1
 -1.554 / 1:       0.
 -1.521 / 2:   60693.
 -1.488 / 3:   10663.
 -1.455 / 4: -153460.
 0.033  / 5: -375011.
 1.455  / 6: -640073.
 1.488  / 7: -700065.
 1.521  / 8: -621550.
 1.554  / 9: -642176.

What happens is that cd_consistent_axis_orient.F has been trying to use the an array of logicals marking the axis "reversed". That's set when initializing the first datset from which the axis is read. When looking at the grid for the variable in the second dataset, Ferret has determined that the axes match, but that "reversed" flag is no longer set for an axis that came in from the first datset. The grid has been created, and the info we need is stored in another variable, ds_ordering, which cd_consistent_axis_orient.F can use instead.

Attachment from @AnsleyManke on 29 Feb 2012 00:25 UTC
the file for reversed-axis example
short_data.nc.zip

This last issue with reversed axes is the same as the report in #1941, closed.