# (C) Copyright 2005- ECMWF.

include "grib2/template.3.grid.def"

#  Di - i direction increment
unsigned[4] iDirectionIncrement : can_be_missing;
alias Di = iDirectionIncrement;

#  N - number of parallels between a pole and the equator
unsigned[4] N  : dump;
alias numberOfParallelsBetweenAPoleAndTheEquator=N ;
alias geography.N=N;

include "grib2/template.3.scanning_mode.def"

modify Ni : can_be_missing,dump;

meta g2grid g2grid(
    latitudeOfFirstGridPoint,
    longitudeOfFirstGridPoint,
    latitudeOfLastGridPoint,
    longitudeOfLastGridPoint,
    iDirectionIncrement,
    null,
    basicAngleOfTheInitialProductionDomain,
    subdivisionsOfBasicAngle
);

meta geography.latitudeOfFirstGridPointInDegrees   g2latlon(g2grid,0) : dump;
meta geography.longitudeOfFirstGridPointInDegrees  g2latlon(g2grid,1) : dump;
meta geography.latitudeOfLastGridPointInDegrees    g2latlon(g2grid,2) : dump;
meta geography.longitudeOfLastGridPointInDegrees   g2latlon(g2grid,3) : dump;
meta geography.iDirectionIncrementInDegrees        g2latlon(g2grid,4,iDirectionIncrementGiven) : can_be_missing,dump;

meta global global_gaussian(
    N,Ni,iDirectionIncrement,
    latitudeOfFirstGridPoint,
    longitudeOfFirstGridPoint,
    latitudeOfLastGridPoint,
    longitudeOfLastGridPoint,
    PLPresent, pl,
    basicAngleOfTheInitialProductionDomain,
    subdivisionsOfBasicAngle
) = 0 : dump;

alias xFirst=longitudeOfFirstGridPointInDegrees;
alias yFirst=latitudeOfFirstGridPointInDegrees;
alias xLast=longitudeOfLastGridPointInDegrees;
alias yLast=latitudeOfLastGridPointInDegrees;

alias latitudeFirstInDegrees  = latitudeOfFirstGridPointInDegrees;
alias longitudeFirstInDegrees = longitudeOfFirstGridPointInDegrees;
alias latitudeLastInDegrees   = latitudeOfLastGridPointInDegrees;
alias longitudeLastInDegrees  = longitudeOfLastGridPointInDegrees;
alias DiInDegrees             = iDirectionIncrementInDegrees;

if (missing(Ni) && PLPresent == 1) {
  iterator gaussian_reduced(numberOfPoints,missingValue,values,
        latitudeOfFirstGridPointInDegrees,longitudeOfFirstGridPointInDegrees,
        latitudeOfLastGridPointInDegrees,longitudeOfLastGridPointInDegrees,
        N,pl,Nj);
  nearest reduced(values,radius,Nj,pl);

  #meta sumPlArray sum(pl);
  #meta dataGlobal evaluate( sumPlArray == (numberOfValues+numberOfMissing) );
} else {
  iterator gaussian(numberOfPoints,missingValue,values,
        longitudeFirstInDegrees,DiInDegrees,
        Ni,Nj,iScansNegatively,
        latitudeFirstInDegrees, latitudeLastInDegrees,
        N,jScansPositively);
  nearest regular(values,radius,Ni,Nj);
}
meta latLonValues latlonvalues(values);
alias latitudeLongitudeValues=latLonValues;
meta latitudes latitudes(values,0);
meta longitudes longitudes(values,0);
meta distinctLatitudes latitudes(values,1);
meta distinctLongitudes longitudes(values,1);

meta isOctahedral octahedral_gaussian(N, Ni, PLPresent, pl) = 0 : no_copy,dump;

meta gaussianGridName gaussian_grid_name(N, Ni, isOctahedral);
alias gridName=gaussianGridName;


# For sub-areas
# Uses new algorithm for counting. No support for legacy mode
meta numberOfDataPointsExpected number_of_points_gaussian(Ni,Nj,PLPresent,pl,N,
  latitudeOfFirstGridPointInDegrees,longitudeOfFirstGridPointInDegrees,
  latitudeOfLastGridPointInDegrees,longitudeOfLastGridPointInDegrees,zero) : dump;

meta legacyGaussSubarea evaluate(numberOfDataPoints != numberOfDataPointsExpected);