[MITgcm-devel] netcdf units in diagnostics pkg

[Martin Losch]

Dear diagnostics-authors,

Sorry for bothering you with this, but I would like to understand the  
logic of the routine body getdiag (diagnostics_utils.F, see below):

The first case (mate=0) happens (according to my observations, not my  
understanding of the code, which I have to admit I don't understand)  
whenever diagnostics_fract_fill has been used to fill the diagnostics  
(e.g. SI_Thick), the second case happens when diagnostics_fill has  
been used to fill the diagnostics (e.g., SI_Fract). In the first  
case, any value above "undef" (hardwired to getcon('UNDEF') = 1D15)  
is reset to "undef", in the second case "undef" is only used if the  
diagnostics has been filled less than once (=not at all).

My questions:
1. In the above example (SI_Fract and SI_Thick), SI_Fract contains  
zeros where there is no ice (iceMask=0), good!. SI_Thick, which  
should be iceHeight*iceMask (and should be zero where iceMask=0),  
there are values of undef=getcon('UNDEF')=1D15 where there is no ice.  
Now, whether this is useful or not, I just don't understand, where  
these numbers come from (in getdiags, this value is only assigned if  
the diagnostics exceeds "undef", or has not been filled).
2. is it a good idea to reset a diagnostic to "undef", when it is  
larger than "undef", while "undef" is hardwired? As a counter  
example: I can image averaging some pressure variable (e.g. bottom  
pressure anomaly) over many timesteps (or maybe even the square of  
bottom pressure) and would expect to get really large numbers that  
are then capped by this logic.

Martin

>         IF ( mate.EQ.0 ) THEN
> C-      No counter diagnostics => average = Sum / ndiag :
>
>           ipnt = ip + lev - 1
>           factor = FLOAT(ndiag(ip,bi,bj))
>           IF (ndiag(ip,bi,bj).NE.0) factor = 1. _d 0 / factor
>
>           DO j = 1,sNy+1
>             DO i = 1,sNx+1
>               IF ( qdiag(i,j,ipnt,bi,bj) .LE. undef ) THEN
>                 qtmp(i,j) = qdiag(i,j,ipnt,bi,bj)*factor
>               ELSE
>                 qtmp(i,j) = undef
>               ENDIF
>             ENDDO
>           ENDDO
>
>         ELSE
> C-      With counter diagnostics => average = Sum / counter:
>
>           ipnt = ip + lev - 1
>           levCt= MIN(lev,kdiag(mate))
>           ipCt = im + levCt - 1
>           DO j = 1,sNy+1
>             DO i = 1,sNx+1
>               IF ( qdiag(i,j,ipCt,bi,bj) .NE. 0. ) THEN
>                 qtmp(i,j) = qdiag(i,j,ipnt,bi,bj)
>      &                    / qdiag(i,j,ipCt,bi,bj)
>               ELSE
>                 qtmp(i,j) = undef
>               ENDIF
>             ENDDO
>           ENDDO
>
>         ENDIF
>        ENDIF
>