[metview] 04/13: New upstream release 4.6.1
Alastair McKinstry
mckinstry at moszumanska.debian.org
Thu Feb 11 19:24:28 UTC 2016
This is an automated email from the git hooks/post-receive script.
mckinstry pushed a commit to tag debian/4.6.1-1
in repository metview.
commit b8dbd30421b7ad14fb545f833c5be6bed2295530
Author: Alastair McKinstry <mckinstry at debian.org>
Date: Mon Jan 4 03:48:35 2016 +0000
New upstream release 4.6.1
---
VERSION.cmake | 2 +-
share/metview/etc/CMakeLists.txt | 1 +
share/metview/etc/GeoToGribDef | 26 +-
share/metview/etc/GeoToGribRules | 7 +
share/metview/etc/ObjectList | 3 +-
share/metview/etc/ecmwf.def | 5 +-
src/Macro/grib.cc | 77 ++++-
src/Macro/misc.cc | 1 +
src/MvApp/geo_to_grib.cc | 541 +++++++++++++++++++------------
src/XSection/Average.cc | 2 +-
src/XSection/CrossS.cc | 2 +-
src/libMars/CMakeLists.txt | 14 +-
src/libMarsClient/version.c | 9 +-
src/libMetview/MvGrid.h | 4 +
src/libUtil/MvVersionInfo.cc | 4 +
src/libUtil/MvVersionInfo.h | 28 +-
src/uPlot/OutputFormatAction.cc | 24 +-
src/uPlot/PlotAction.cc | 6 +-
test/macros/CMakeLists.txt | 6 +
test/macros/gpt_to_grib_with_template.mv | 98 ++++++
test/macros/grib_to_gpt.mv | 49 +++
21 files changed, 659 insertions(+), 250 deletions(-)
diff --git a/VERSION.cmake b/VERSION.cmake
index 14a8e54..05233b2 100644
--- a/VERSION.cmake
+++ b/VERSION.cmake
@@ -1,5 +1,5 @@
set(${PROJECT_NAME}_MAJOR_VERSION_STR "4")
set(${PROJECT_NAME}_MINOR_VERSION_STR "6")
-set(${PROJECT_NAME}_REVISION_VERSION_STR "0")
+set(${PROJECT_NAME}_REVISION_VERSION_STR "1")
set(${PROJECT_NAME}_VERSION_STR "${${PROJECT_NAME}_MAJOR_VERSION_STR}.${${PROJECT_NAME}_MINOR_VERSION_STR}.${${PROJECT_NAME}_REVISION_VERSION_STR}")
diff --git a/share/metview/etc/CMakeLists.txt b/share/metview/etc/CMakeLists.txt
index aaf4cb9..3c634a5 100644
--- a/share/metview/etc/CMakeLists.txt
+++ b/share/metview/etc/CMakeLists.txt
@@ -30,6 +30,7 @@ set(files
GenAppDef
GenAppRules
GeoToGribDef
+ GeoToGribRules
GeoViewDef GeoViewRules
GraphDef
GraphRules
diff --git a/share/metview/etc/GeoToGribDef b/share/metview/etc/GeoToGribDef
index 13f6543..16cd08c 100644
--- a/share/metview/etc/GeoToGribDef
+++ b/share/metview/etc/GeoToGribDef
@@ -14,15 +14,22 @@ GEO_TO_GRIB; Transform Geopoints to regular grid; metview
}
+ GRID_DEFINITION_MODE
+ {
+ USER
+ GRIB
+ } = USER
+
+
AREA
[
- help = help_input,
- input_type = area,
- input_window = '/Metview/Defaults/Input Window'
+ help = help_input,
+ input_type = area,
+ input_window = '/Metview/Defaults/Input Window'
]
{
- *
- /
+ *
+ /
} = 90/-180/-90/180
@@ -33,6 +40,15 @@ GEO_TO_GRIB; Transform Geopoints to regular grid; metview
} = 1.5/1.5
+ TEMPLATE_GRIB
+ [
+ interface = icon, class = GRIB, exclusive = true
+ ]
+ {
+ @
+ }
+
+
TOLERANCE
{
*
diff --git a/share/metview/etc/GeoToGribRules b/share/metview/etc/GeoToGribRules
new file mode 100644
index 0000000..ab8fff2
--- /dev/null
+++ b/share/metview/etc/GeoToGribRules
@@ -0,0 +1,7 @@
+
+%if GRID_DEFINITION_MODE <> USER %then
+ %unset AREA
+ %unset GRID
+
+%if GRID_DEFINITION_MODE <> GRIB %then
+ %unset TEMPLATE_GRIB
diff --git a/share/metview/etc/ObjectList b/share/metview/etc/ObjectList
index 0736570..a32be7b 100644
--- a/share/metview/etc/ObjectList
+++ b/share/metview/etc/ObjectList
@@ -1102,10 +1102,11 @@ object,
can_be_created = True,
check = True,
type = Data,
- macro = geo_to_grib,
+ macro = geo_to_grib,
pixmap = '$METVIEW_DIR_SHARE/icons/GEO_TO_GRIB.icon',
editor_type = SimpleEditor,
definition_file = '$METVIEW_DIR_SHARE/etc/GeoToGribDef',
+ rules_file = '$METVIEW_DIR_SHARE/etc/GeoToGribRules',
default_object = False,
help_page = Geopoints_To_GRIB,
default_name = Geopoints to GRIB
diff --git a/share/metview/etc/ecmwf.def b/share/metview/etc/ecmwf.def
index 20a7c26..f9f7614 100755
--- a/share/metview/etc/ecmwf.def
+++ b/share/metview/etc/ecmwf.def
@@ -1070,7 +1070,7 @@ GRIB IMAGES #= (0/1/10/20)
# Extractions
-BUFR = (1/3/9/11/13/19/21/22/23/28/61/62/63/82/83/84/85/86/87/91/92/95/96/101/102/103/106/110/109/111/112/113/140/142/144/145/146/147/164/170/172/176/178/180)
+BUFR = (1/3/9/11/13/19/21/22/23/28/61/62/63/82/83/84/85/86/87/91/92/95/96/101/102/103/106/110/109/111/112/113/140/142/144/145/146/147/149/164/170/172/176/178/180)
WAVESCAT = (121)
SCAT = (122)
RALT = (123)
@@ -1211,7 +1211,7 @@ SATELLITE = (VSS/SSBT/SSMI/SSMIS/QSCAT/REO3/WSAT/VASS/ATMS)
# Type SLNS (NSD)
-SINGLE LEVEL UPPER AIR NOT SATELLITE ; SLNS = (141/142/143/144/145/146)
+SINGLE LEVEL UPPER AIR NOT SATELLITE ; SLNS = (141/142/143/144/145/146/149)
CODAR ; COD ; 141
AIREP ; AIR ; 142
@@ -1219,6 +1219,7 @@ SINGLE LEVEL UPPER AIR NOT SATELLITE ; SLNS = (141/142/143/144/145/146)
AMDAR ; AM ; 144
ACARS ; AC ; 145
E-AMDAR ; EAM ; 146
+ ACARS with mixing ratio ; ACMR ; 149
# Type SLS (NSD)
diff --git a/src/Macro/grib.cc b/src/Macro/grib.cc
index 1033a92..4ade188 100644
--- a/src/Macro/grib.cc
+++ b/src/Macro/grib.cc
@@ -2137,9 +2137,9 @@ Value SetGridValsFunction::Execute(int /*arity*/,Value *arg)
release_mem(h->values);
h->values = (double*)reserve_mem(sizeof(double)*h->value_count);
}
- else if (v->Count() != h->value_count) // otherwise, issue an error
+ else if (v->Count() != (signed)h->value_count) // otherwise, issue an error
{
- return Error( "set_gridvals: input vector has %d points, field has %d - they should be the same.", v->Count(), h->value_count);
+ return Error( "set_gridvals: input vector has %d points, field has %d - they should be the same.", v->Count(), h->value_count);
}
@@ -3369,8 +3369,8 @@ Value GFindFunction::Execute(int arity,Value *arg)
grd->lat_y(),
grd->lon_x(),
0.0, // level
- 0, // date
- 0, // time
+ 0L, // date
+ 0L, // time
grd->value() );
}
}
@@ -5658,6 +5658,73 @@ Value FrequenciesFunction::Execute(int arity,Value *arg)
}
//=============================================================================
+// Fill missing values along the horizontal line.
+// Implemented only for regular latlong grid format.
+// For each latitude line, analyses each point from left to right. If a point is
+// a missing value then replaces it by either the previous point, if the next
+// point is a missing value, or by the next point.
+
+class FillMVEWFunction : public Function {
+public:
+ FillMVEWFunction(const char *n) : Function(n,1,tgrib)
+ { info = "Fill missing values along the horizontal line"; }
+ virtual Value Execute(int arity,Value *arg);
+};
+
+Value FillMVEWFunction::Execute(int,Value *arg)
+{
+ double v1;
+ fieldset *v;
+
+ arg[0].GetValue(v);
+ fieldset *z = copy_fieldset(v,v->count,false);
+
+ for( int i = 0; i < v->count; i++)
+ {
+ auto_ptr<MvGridBase> grd( MvGridFactory( v->fields[i] ) );
+ if(! grd->isRegularLatLongGrid() )
+ return Error( "fill_missing_values_ew: implemented only for regular latlong format" );
+
+ auto_ptr<MvGridBase> newGrd( MvGridFactory( z->fields[i] ) );
+
+ // Copy first point
+ v1 = grd->value();
+ newGrd->value(v1);
+ grd->advance();
+ newGrd->advance();
+
+ // Loop from the second to the one before the last point
+ for( int j = 1; j < grd->length()-1; j++ )
+ {
+ if ( grd->hasValue() )
+ {
+ v1 = grd->value();
+ newGrd->value(v1);
+ grd->advance();
+ }
+ else
+ {
+ grd->advance();
+ if ( grd->hasValue() )
+ newGrd->value(grd->value());
+ else
+ newGrd->value(v1);
+ }
+
+ newGrd->advance();
+ }
+
+ // Compute the last point
+ if ( grd->hasValue() )
+ newGrd->value(grd->value());
+ else
+ newGrd->value(v1);
+ }
+
+ return Value(z);
+}
+
+//=============================================================================
static void install(Context *c)
{
@@ -5752,7 +5819,7 @@ static void install(Context *c)
c->AddFunction(new GribHeaderFunctionW("grib_set_string", GRIB_STRING));
c->AddFunction(new GribHeaderFunctionWGeneric("grib_set"));
-
+ c->AddFunction(new FillMVEWFunction("fill_missing_values_ew"));
}
static Linkage linkage(install);
diff --git a/src/Macro/misc.cc b/src/Macro/misc.cc
index 5d27e17..421fc7a 100644
--- a/src/Macro/misc.cc
+++ b/src/Macro/misc.cc
@@ -775,6 +775,7 @@ Value MetviewVersionFunction::Execute(int arity,Value *arg)
set_value(r, "metview_version", "%d", mvInfo.version());
set_value(r, "metview_major", "%d", mvInfo.majorVersion());
set_value(r, "metview_minor", "%d", mvInfo.minorVersion());
+ set_value(r, "metview_revision","%d", mvInfo.revision());
}
else
diff --git a/src/MvApp/geo_to_grib.cc b/src/MvApp/geo_to_grib.cc
index 64e7fad..84ededa 100644
--- a/src/MvApp/geo_to_grib.cc
+++ b/src/MvApp/geo_to_grib.cc
@@ -36,182 +36,217 @@ double makeDateNumber( long date, long time );
class ToMatrix
{
- MvGeoPoints GPoints;
-
- double North;
- double South;
- double West;
- double East;
- double GridLat;
- double GridLon;
- double Tolerance;
- int Weight_;
- int NbLat;
- int NbLon;
- double Date;
- int Parameter;
- int Table2;
- double *Matrix;
- vector< list<MvGeoP1> > LatListVec;
+ MvGeoPoints GPoints;
+
+ double North;
+ double South;
+ double West;
+ double East;
+ double GridLat;
+ double GridLon;
+ double Tolerance;
+ int Weight_;
+ int NbLat;
+ int NbLon;
+ int numValsInTemplateGrib;
+ double Date;
+ int Parameter;
+ int Table2;
+ double *Matrix;
+ vector< list<MvGeoP1> > LatListVec;
+ MvField *templateField; // template field - defines the output grid
public:
- ToMatrix(MvRequest&);
- int load(const char*);
- int save(const char*);
- void estimate();
- double value(float lat, float lon);
- void sortPoints();
+ enum GribDefinitionMode {User, Grib};
+
+ ToMatrix(MvRequest&);
+ int load(const char*);
+ int save(string &);
+ bool estimate();
+ bool estimateToUserGrid();
+ bool estimateToTemplateGrid();
+ double value(float lat, float lon);
+ void sortPoints();
+ string &errorMessage() {return errorMessage_;};
+ GribDefinitionMode mode() {return gribDefMode;};
protected:
- float long180( float x ){ return x>180 ? x-360 : x; }
- float long360( float x ){ return x<0 ? x+360 : x; }
- float minDistance( float x1, float x2 );
- int latBand( float lat );
+ GribDefinitionMode gribDefMode;
+
+ float long180( float x ){ return x>180 ? x-360 : x; }
+ float long360( float x ){ return x<0 ? x+360 : x; }
+ float minDistance( float x1, float x2 );
+ int latBand( float lat );
+ string errorMessage_;
};
//_____________________________________________________________________
class GeoToGRIB : public MvService {
public:
- GeoToGRIB() : MvService("GEO_TO_GRIB") {};
- void serve(MvRequest&,MvRequest&);
+ GeoToGRIB() : MvService("GEO_TO_GRIB") {};
+ void serve(MvRequest&,MvRequest&);
};
//_____________________________________________________________________
double makeDateNumber( long date, long time )
{
- double myDate = date;
- double myTime = (double)time;
+ double myDate = date;
+ double myTime = (double)time;
- if( myDate < 470620 ) //-- must be relative => no time
- return myDate;
+ if( myDate < 470620 ) //-- must be relative => no time
+ return myDate;
- if( myTime < 0 ) //-- illegal or rubbish
- return myDate;
+ if( myTime < 0 ) //-- illegal or rubbish
+ return myDate;
- if( myTime > 2400 ) //-- illegal or rubbish
- return myDate;
+ if( myTime > 2400 ) //-- illegal or rubbish
+ return myDate;
- if( myTime >= 100 ) //-- 01:00...24:00
- return myDate + myTime / 2400.;
+ if( myTime >= 100 ) //-- 01:00...24:00
+ return myDate + myTime / 2400.;
- if( myTime > 59 ) //-- illegal or rubbish
- return myDate;
+ if( myTime > 59 ) //-- illegal or rubbish
+ return myDate;
- if( myTime > 24 ) //-- 00:25...00:59
- return myDate + myTime/60. / 24.;
+ if( myTime > 24 ) //-- 00:25...00:59
+ return myDate + myTime/60. / 24.;
- return myDate + myTime / 24.;
+ return myDate + myTime / 24.;
}
//_____________________________________________________________________
ToMatrix::ToMatrix(MvRequest& def)
{
- North = def("AREA", 0);
- West = def("AREA", 1);
- South = def("AREA", 2);
- East = def("AREA", 3);
- GridLon = def("GRID");
- GridLat = def("GRID", 1);
- Tolerance = def("TOLERANCE");
- Parameter = def("PARAMETER");
- Table2 = def("GRIB_TABLE2_VERSION");
- string st = (const char*)def("INTERPOLATION_METHOD");
-
- if ( st == "RECIPROCAL" )
- Weight_ = G2G_RECIPROCAL;
- else if ( st == "EXPONENTIAL_MEAN" )
- Weight_ = G2G_EXPONENTIAL_MEAN;
- else
- Weight_ = G2G_EXPONENTIAL_SUM;
-
- if(North < South)
- {
- double tmp = North;
- North = South;
- South = tmp;
- }
+ string mode = (const char*)def("GRID_DEFINITION_MODE");
+
+ if (mode == "USER")
+ {
+ gribDefMode = User;
+ North = def("AREA", 0);
+ West = def("AREA", 1);
+ South = def("AREA", 2);
+ East = def("AREA", 3);
+ GridLon = def("GRID");
+ GridLat = def("GRID", 1);
+
+ if(North < South)
+ {
+ double tmp = North;
+ North = South;
+ South = tmp;
+ }
+
+ if( !GridLon )
+ GridLon = 1.5;
+ if( !GridLat )
+ GridLat = GridLon;
+ }
+ else
+ {
+ gribDefMode = Grib; // user supplied a template GRIB file to define the output grid
+
+ //-- data request => MvField
+ MvRequest grb;
+ def.getValue(grb, "TEMPLATE_GRIB");
+ fieldset* fs = request_to_fieldset((request*)grb);
+ field* f = get_field(fs,0,expand_mem);
+ templateField = new MvField(f);
+
+ if (fs->count > 1) //-- warning if several fields
+ {
+ marslog( LOG_INFO, "Info: only first field in the template GRIB will be considered" );
+ }
+ }
+
- if( !GridLon )
- GridLon = 1.5;
- if( !GridLat )
- GridLat = GridLon;
+ Tolerance = def("TOLERANCE");
+ Parameter = def("PARAMETER");
+ Table2 = def("GRIB_TABLE2_VERSION");
+ string st = (const char*)def("INTERPOLATION_METHOD");
+ if ( st == "RECIPROCAL" )
+ Weight_ = G2G_RECIPROCAL;
+ else if ( st == "EXPONENTIAL_MEAN" )
+ Weight_ = G2G_EXPONENTIAL_MEAN;
+ else
+ Weight_ = G2G_EXPONENTIAL_SUM;
- MvRequest data;
- def.getValue(data,"GEOPOINTS");
- const char* path = data("PATH");
- GPoints.load( path );
- Date = makeDateNumber( GPoints[0].date(), GPoints[0].time() );
-
- list<MvGeoP1> emptyList;
- LatListVec.assign( cLatBandCount+1, emptyList );
+ MvRequest data;
+ def.getValue(data,"GEOPOINTS");
+ const char* path = data("PATH");
+
+ GPoints.load( path );
+ Date = makeDateNumber( GPoints[0].date(), GPoints[0].time() );
+
+ list<MvGeoP1> emptyList;
+ LatListVec.assign( cLatBandCount+1, emptyList );
}
//_____________________________________________________________________
double ToMatrix::value(float lat, float lon)
{
- double val = MISSING_DATA;
- double dist, coef;
- double sigma = 0;
-
- int band1 = latBand( min( (int)((lat+Tolerance)+0.5), 90) );
- int band2 = latBand( max( (int)((lat-Tolerance)-0.5),-90) );
-
- for( int b=band1; b>=band2; --b)
- {
- for( list<MvGeoP1>::iterator gp_iter = LatListVec[b].begin();
- gp_iter != LatListVec[b].end();
- ++gp_iter )
- {
- float pi_lat = gp_iter->lat_y();
-
- if( ( fabs(lat - pi_lat) ) > Tolerance )
- continue;
-
- float pi_lon = gp_iter->lon_x();
- if( minDistance( pi_lon, lon ) > Tolerance )
- continue;
-
- //-- Here we have found a point inside the Interval;
- if( ! gp_iter->value_missing() )
- {
- double x = minDistance( lon, pi_lon );
- double y = lat - pi_lat;
- dist = (x*x) + (y*y);
-
- // Compute weight
- if ( Weight_ == G2G_RECIPROCAL )
- {
- if( dist == 0 )
- return gp_iter->value(); //-- Here the point is on the Grid
- dist = 1/dist;
- }
- else // exponential
- {
- if ( Tolerance != 0 )
- dist = exp(-(dist/(pow(Tolerance,2))));
- else
- dist = dist ? 0 : 1;
- }
-
- sigma +=dist;
- coef = dist * gp_iter->value();
-
- if( val == MISSING_DATA )
- val = coef;
- else
- val += coef;
- }
- }
- }
+ double val = MISSING_DATA;
+ double dist, coef;
+ double sigma = 0;
+
+ int band1 = latBand( min( (int)((lat+Tolerance)+0.5), 90) );
+ int band2 = latBand( max( (int)((lat-Tolerance)-0.5),-90) );
- if (sigma && Weight_ != G2G_EXPONENTIAL_SUM )
- val = val / sigma;
+ for( int b=band1; b>=band2; --b)
+ {
+ for( list<MvGeoP1>::iterator gp_iter = LatListVec[b].begin();
+ gp_iter != LatListVec[b].end();
+ ++gp_iter )
+ {
+ float pi_lat = gp_iter->lat_y();
+
+ if( ( fabs(lat - pi_lat) ) > Tolerance )
+ continue;
+
+ float pi_lon = gp_iter->lon_x();
+ if( minDistance( pi_lon, lon ) > Tolerance )
+ continue;
+
+ //-- Here we have found a point inside the Interval;
+ if( ! gp_iter->value_missing() )
+ {
+ double x = minDistance( lon, pi_lon );
+ double y = lat - pi_lat;
+ dist = (x*x) + (y*y);
+
+ // Compute weight
+ if ( Weight_ == G2G_RECIPROCAL )
+ {
+ if( dist == 0 )
+ return gp_iter->value(); //-- Here the point is on the Grid
+ dist = 1/dist;
+ }
+ else // exponential
+ {
+ if ( Tolerance != 0 )
+ dist = exp(-(dist/(pow(Tolerance,2))));
+ else
+ dist = dist ? 0 : 1;
+ }
+
+ sigma +=dist;
+ coef = dist * gp_iter->value();
+
+ if( val == MISSING_DATA )
+ val = coef;
+ else
+ val += coef;
+ }
+ }
+ }
- return val;
+ if (sigma && Weight_ != G2G_EXPONENTIAL_SUM )
+ val = val / sigma;
+
+ return val;
}
//_____________________________________________________________________
@@ -219,18 +254,18 @@ double ToMatrix::value(float lat, float lon)
float
ToMatrix::minDistance( float x1, float x2 )
{
- // calculate two ways, in case given values are
- // on the different sides of discontinuity line!
+ // calculate two ways, in case given values are
+ // on the different sides of discontinuity line!
- float min1 = long180( x1 ) - long180( x2 );
- if( min1 < 0.0 )
- min1 = -min1;
+ float min1 = long180( x1 ) - long180( x2 );
+ if( min1 < 0.0 )
+ min1 = -min1;
- float min2 = long360( x1 ) - long360( x2 );
- if( min2 < 0.0 )
- min2 = -min2;
+ float min2 = long360( x1 ) - long360( x2 );
+ if( min2 < 0.0 )
+ min2 = -min2;
- return min1 < min2 ? min1 : min2;
+ return min1 < min2 ? min1 : min2;
}
//_____________________________________________________________________
@@ -243,21 +278,30 @@ ToMatrix::latBand( float lat )
void ToMatrix::sortPoints()
{
- cout << "ToMatrix::sortPoints(): Lat Band Size & Count: "
- << cLatBandSize
- << ", "
- << cLatBandCount
- << endl;
-
- for (int s = 0; s < GPoints.count(); s++)
- {
- int band = latBand( GPoints[s].lat_y() );
- LatListVec[band].push_back(GPoints[s]);
- }
+ cout << "ToMatrix::sortPoints(): Lat Band Size & Count: "
+ << cLatBandSize
+ << ", "
+ << cLatBandCount
+ << endl;
+
+ for (int s = 0; s < GPoints.count(); s++)
+ {
+ int band = latBand( GPoints[s].lat_y() );
+ LatListVec[band].push_back(GPoints[s]);
+ }
}
//_____________________________________________________________________
-void ToMatrix::estimate()
+bool ToMatrix::estimate()
+{
+ if (gribDefMode == User)
+ return estimateToUserGrid();
+ else
+ return estimateToTemplateGrid();
+}
+//_____________________________________________________________________
+
+bool ToMatrix::estimateToUserGrid()
{
int i, j;
float lat,lon;
@@ -302,46 +346,111 @@ void ToMatrix::estimate()
{
lon = West;
for (i = 0; i < NbLon; i++)
- {
- Matrix[i + ( j*NbLon) ] = value(lat, lon);
- lon += GridLon;
- }
+ {
+ Matrix[i + ( j*NbLon) ] = value(lat, lon);
+ lon += GridLon;
+ }
lat -= GridLat;
}
+
+ return true;
}
//_____________________________________________________________________
-int ToMatrix::save(const char * path)
+bool ToMatrix::estimateToTemplateGrid()
{
- FILE *f = fopen(path,"w");
- if(f == NULL) return 1;
-
- fprintf(f, "#LLMATRIX\n");
- fprintf(f, "GRID = %g/%g\n", GridLon, GridLat);
- fprintf(f, "NORTH = %g\n", North);
- fprintf(f, "SOUTH = %g\n", South);
- fprintf(f, "EAST = %g\n", East);
- fprintf(f, "WEST = %g\n", West);
- fprintf(f, "NLAT = %d\n", NbLat);
- fprintf(f, "NLON = %d\n", NbLon);
- fprintf(f, "DATE = %fl\n", Date);
- fprintf(f, "MISSING = '%g'\n", MISSING_DATA);
- fprintf(f, "PARAM = %d\n", Parameter );
- fprintf(f, "TABLE2 = %d\n", Table2 );
- fprintf(f, "#DATA\n");
-
- int last;
-
- for (int j = 0; j < NbLat; j++) {
- for (int i = 0; i < NbLon; i++) {
- last = i + ( j*NbLon);
- fprintf(f, "%g\n", Matrix[i + ( j*NbLon)]);
- }
+ MvGridBase *grd = templateField->mvGrid(); // be careful with this pointer because MvField is in charge!
+
+ if (!grd->hasLocationInfo()) //-- check that it is ok
+ {
+ errorMessage_ = "Template GRIB field grid locations cannot be read by Metview.";
+ return false;
+ }
+
+
+ // for each gridpoint, get the location and compute the geopoint value there
+
+ Matrix = new double[templateField->countValues()];
+
+ MvGridPoint gp;
+ int i = 0;
+
+ do
+ {
+ gp = grd->gridPoint();
+
+ Matrix[i] = value(gp.loc_.latitude(), gp.loc_.longitude());
+
+ if (Matrix[i] == GEOPOINTS_MISSING_VALUE) // convert missing value
+ Matrix[i] = grd->missingValue();
+
+ i++;
}
- fclose(f);
- printf("indice = %d \n", last);
+ while (grd->advance());
+
+ numValsInTemplateGrib = i;
+
+ return true;
+}
+//_____________________________________________________________________
- return 0;
+int ToMatrix::save(string &path) // returns path to caller
+{
+
+ if (gribDefMode == User)
+ {
+ // create an LLMatrix text file to be later converted to GRIB
+
+ path = marstmp();
+ FILE *f = fopen(path.c_str(),"w");
+ if(f == NULL) return 1;
+
+ fprintf(f, "#LLMATRIX\n");
+ fprintf(f, "GRID = %g/%g\n", GridLon, GridLat);
+ fprintf(f, "NORTH = %g\n", North);
+ fprintf(f, "SOUTH = %g\n", South);
+ fprintf(f, "EAST = %g\n", East);
+ fprintf(f, "WEST = %g\n", West);
+ fprintf(f, "NLAT = %d\n", NbLat);
+ fprintf(f, "NLON = %d\n", NbLon);
+ fprintf(f, "DATE = %fl\n", Date);
+ fprintf(f, "MISSING = '%g'\n", MISSING_DATA);
+ fprintf(f, "PARAM = %d\n", Parameter );
+ fprintf(f, "TABLE2 = %d\n", Table2 );
+ fprintf(f, "#DATA\n");
+
+ int last;
+
+ for (int j = 0; j < NbLat; j++) {
+ for (int i = 0; i < NbLon; i++) {
+ last = i + ( j*NbLon);
+ fprintf(f, "%g\n", Matrix[i + ( j*NbLon)]);
+ }
+ }
+ fclose(f);
+ printf("indice = %d \n", last);
+ }
+ else
+ {
+ // clone the template GRIB field and copy the computed values to it
+
+ field *z = copy_field(templateField->libmars_field(), false);
+ memcpy(z->values, Matrix, numValsInTemplateGrib*sizeof(double));
+ MvField resultField(z);
+
+ if (Parameter != 255)
+ {
+ resultField.setGribKeyValueLong("table2Version", Table2);
+ resultField.setGribKeyValueLong("paramId", Parameter); // paramId to ensure GRIB2 ok
+ }
+
+ MvFieldSet fs;
+ fs += resultField;
+ MvRequest req (fs.getRequest()); // this forces save-to-file
+ path = (const char*)req("PATH"); // return path to caller
+ }
+
+ return 0;
}
//_____________________________________________________________________
@@ -355,29 +464,51 @@ void GeoToGRIB::serve( MvRequest& in, MvRequest& out )
matrix.sortPoints();
timer.lapTime( "sortPoints" );
- matrix.estimate();
+ if (!matrix.estimate())
+ {
+ setError(13, matrix.errorMessage().c_str());
+ return;
+ }
timer.lapTime( "estimate" );
- string path = marstmp();
- if( matrix.save( path.c_str() ) )
- {
- setError( 1, "geo_to_grib: can not write matrix into file %s", path.c_str() );
- return;
- }
- string tmp = marstmp();
+ string tmpGribFile;
- if( LLMatrixToGRIB( path.c_str(), tmp.c_str() ) )
- {
- setError( 1, "geo_to_grib: convertion from geo/LLMatrix to GRIB failed" );
- return;
- }
+ if (matrix.mode() == ToMatrix::User)
+ {
+ // create an intermediate LLMatrix file, then convert to GRIB
- unlink(path.c_str()); // remove the temporary LLMATRIX file
+ string path;
+ if (matrix.save(path))
+ {
+ setError( 1, "geo_to_grib: can not write matrix into file %s", path.c_str() );
+ return;
+ }
+
+ tmpGribFile = marstmp();
+
+ if( LLMatrixToGRIB( path.c_str(), tmpGribFile.c_str() ) )
+ {
+ setError( 1, "geo_to_grib: convertion from geo/LLMatrix to GRIB failed" );
+ return;
+ }
+
+ unlink(path.c_str()); // remove the temporary LLMATRIX file
+ }
+ else
+ {
+ // write the GRIB directly rather than creating an intermediate LLMatrix file
+
+ if (matrix.save(tmpGribFile))
+ {
+ setError( 1, "geo_to_grib: can not write GRIB into file %s", tmpGribFile.c_str() );
+ return;
+ }
+ }
MvRequest grib = "GRIB";
- grib("PATH") = tmp.c_str();
+ grib("PATH") = tmpGribFile.c_str();
grib("TEMPORARY") = 1;
out = grib;
diff --git a/src/XSection/Average.cc b/src/XSection/Average.cc
index 7c27c1c..17d6e1f 100644
--- a/src/XSection/Average.cc
+++ b/src/XSection/Average.cc
@@ -188,6 +188,7 @@ MvRequest Average::createOutputRequest( ApplicationInfo& appInfo, ParamInfo* par
out1("NETCDF_VALUE_VARIABLE") = varname.c_str();
out1("NETCDF_Y_VARIABLE") = getNetcdfLevelVarname(varname).c_str();
out1("NETCDF_X_VARIABLE") = (inputMode_ == AVERAGE_NS ) ? "lon" : "lat";
+ out1("NETCDF_DATA") = xs;
// Create customised Cartesian View request.
// If an icon was dropped into a view, uPlot will ignore it.
@@ -214,7 +215,6 @@ MvRequest Average::createOutputRequest( ApplicationInfo& appInfo, ParamInfo* par
}
// Final output request
- out1("NETCDF_DATA") = xs;
MvRequest out = viewReq + out1;
return out;
}
diff --git a/src/XSection/CrossS.cc b/src/XSection/CrossS.cc
index 4607435..5c9379e 100644
--- a/src/XSection/CrossS.cc
+++ b/src/XSection/CrossS.cc
@@ -788,6 +788,7 @@ MvRequest CrossS::createOutputRequest( ApplicationInfo& appInfo, ParamInfo* parI
// Add flag to indicate that a MWIND visdef should be applied to this data
out1("_VISDEF") = "MWIND";
}
+ out1("NETCDF_DATA") = xs;
// Add the orography curve visualiser
if ( appInfo.haveLNSP() && (appInfo.levelType() == XS_ML_LNSP || appInfo.levelType() == XS_PL) )
@@ -845,7 +846,6 @@ MvRequest CrossS::createOutputRequest( ApplicationInfo& appInfo, ParamInfo* parI
}
// Final output request
- out1("NETCDF_DATA") = xs;
MvRequest out = viewReq + out1;
return out;
}
diff --git a/src/libMars/CMakeLists.txt b/src/libMars/CMakeLists.txt
index d6fb4f7..e5a8407 100644
--- a/src/libMars/CMakeLists.txt
+++ b/src/libMars/CMakeLists.txt
@@ -47,6 +47,17 @@ foreach( file ${_libMars_srcs} )
list( APPEND libMars_srcs ${MARS_SOURCE}/${file} )
endforeach()
+
+# get the MARS client version number (just invent it for now since we don't use it)
+
+set(MARS_CLIENT_MAJOR_VERSION 0)
+set(MARS_CLIENT_MINOR_VERSION 0)
+set(MARS_CLIENT_PATCH_VERSION 0)
+configure_file( ${MARS_SOURCE}/mars_client_version.h.in mars_client_version.h )
+configure_file( ${MARS_SOURCE}/mars_client_version.c.in mars_client_version.c )
+
+
+
# copy all the MARS source files to a subdirectory of the build directory; this is
# so that we have a clean version of the MARS client source code in its own directory
# which can be packaged into the export tarball.
@@ -120,7 +131,8 @@ include_directories(BEFORE ${CMAKE_CURRENT_BINARY_DIR} ${METVIEW_STANDARD_INCL
ecbuild_add_library( TARGET MvMars
TYPE STATIC
- SOURCES ${libMars_srcs}
+ GENERATED mars_client_version.c
+ SOURCES ${libMars_srcs} mars_client_version.c
TEMPLATES ${common_templates}
DEFINITIONS ${METVIEW_EXTRA_DEFINITIONS}
LIBS ${METVIEW_EXTRA_LIBRARIES} ${METVIEW_ODB_API_LIBRARIES} ${CURL_LIBRARIES}
diff --git a/src/libMarsClient/version.c b/src/libMarsClient/version.c
index 0c18bcc..f860d51 100755
--- a/src/libMarsClient/version.c
+++ b/src/libMarsClient/version.c
@@ -8,9 +8,16 @@
* does it submit to any jurisdiction.
*/
+#include "mars_client_version.h"
+#include "mars.h"
long marsversion() {
- static long version = 20150414;
+ static long version = -1;
+
+ if(version == -1)
+ {
+ version = (long)(atol(mars_client_buildstamp())/1000000);
+ }
return version;
}
diff --git a/src/libMetview/MvGrid.h b/src/libMetview/MvGrid.h
index ade7f2b..9f1c8b5 100644
--- a/src/libMetview/MvGrid.h
+++ b/src/libMetview/MvGrid.h
@@ -200,6 +200,10 @@ class MvGridBase
//! Returns the grid type as a string
string gridType() const { return gridType_; }
+ //! Check grid type (implement for the other grids, if it is needed)
+ bool isRegularLatLongGrid()
+ { return gridType_ == cLatLonGrid; }
+
//! Returns the value of the current grid point [iterator]
double value() const
{ return valueAt(currentIndex_); }
diff --git a/src/libUtil/MvVersionInfo.cc b/src/libUtil/MvVersionInfo.cc
index b22ce23..4a7ea10 100644
--- a/src/libUtil/MvVersionInfo.cc
+++ b/src/libUtil/MvVersionInfo.cc
@@ -52,6 +52,10 @@ MvVersionInfo::MvVersionInfo()
copyright_period_ = mv_period;
release_date_ = mv_release_date;
+ char charNameAndVersion[64];
+ sprintf(charNameAndVersion, "%s %d.%d.%d", name_.c_str(), major_, minor_, revision_);
+ nameAndVersion_ = charNameAndVersion;
+
fclose (fp);
info_found_ = true;
}
diff --git a/src/libUtil/MvVersionInfo.h b/src/libUtil/MvVersionInfo.h
index b19b439..13a4a05 100644
--- a/src/libUtil/MvVersionInfo.h
+++ b/src/libUtil/MvVersionInfo.h
@@ -23,19 +23,20 @@ class MvAbstractVersionInfo
public:
virtual ~MvAbstractVersionInfo(){};
- bool infoFound() {return info_found_;};
- int fileVersion() {return file_version_;};
- int majorVersion() {return major_;};
- int minorVersion() {return minor_;};
- int revision() {return revision_;};
- int version() {return version_;};
- int year() {return year_;};
- const string &name() {return name_;};
- const string &label() {return label_;};
- const string &period() {return copyright_period_;};
- const string &releaseDate() {return release_date_;};
- const string &errorMessage() {return error_message_;};
- const string &installDir() {return install_dir_;};
+ bool infoFound() {return info_found_;};
+ int fileVersion() {return file_version_;};
+ int majorVersion() {return major_;};
+ int minorVersion() {return minor_;};
+ int revision() {return revision_;};
+ int version() {return version_;};
+ int year() {return year_;};
+ const string &name() {return name_;};
+ const string &nameAndVersion() {return nameAndVersion_;};
+ const string &label() {return label_;};
+ const string &period() {return copyright_period_;};
+ const string &releaseDate() {return release_date_;};
+ const string &errorMessage() {return error_message_;};
+ const string &installDir() {return install_dir_;};
protected:
MvAbstractVersionInfo(){};
@@ -48,6 +49,7 @@ protected:
int version_;
int year_;
string name_;
+ string nameAndVersion_;
string label_;
string copyright_period_;
string release_date_;
diff --git a/src/uPlot/OutputFormatAction.cc b/src/uPlot/OutputFormatAction.cc
index 5e75960..e10242e 100644
--- a/src/uPlot/OutputFormatAction.cc
+++ b/src/uPlot/OutputFormatAction.cc
@@ -14,6 +14,8 @@
#include "OutputFormatAction.h"
#include "PlotMod.h"
+#include "MvVersionInfo.h"
+
void
OutputFormatAction::Execute (PmContext& context)
{
@@ -47,10 +49,11 @@ OutputFormatAction::Execute (PmContext& context)
}
#endif
- // Magics requires the output filename to have an absolute path
- // Hidden parameter _CWD contains the absolute path. If it is not
- // given then nothing can be done
- bool flag = false;
+ // Update OUTPUT_DEVICES request:
+ // a) Magics requires the output filename to have an absolute path.
+ // Hidden parameter _CWD contains the absolute path. If it is not
+ // given then nothing can be done
+ // b) Add some meta-data about who created the plot
while ( reqdev )
{
if ( (const char*)reqdev("OUTPUT_FULLNAME") && (const char*)req("_CWD") )
@@ -60,7 +63,6 @@ OutputFormatAction::Execute (PmContext& context)
{
name = (string)req("_CWD") + "/" + name; // absolute path
reqdev("OUTPUT_FULLNAME") = name.c_str();
- flag = true;
}
}
@@ -72,19 +74,19 @@ OutputFormatAction::Execute (PmContext& context)
{
name = (string)req("_CWD") + "/" + name; // absolute path
reqdev("OUTPUT_NAME") = name.c_str();
- flag = true;
}
}
+ // Add meta-data OUTPUT_CREATOR
+ MvVersionInfo mvInfo; // the constructor populates the class with information
+ reqdev("OUTPUT_CREATOR") = mvInfo.nameAndVersion().c_str();
+
reqdev.advance();
}
// Update output definition
- if ( flag )
- {
- reqdev.rewind();
- req("OUTPUT_DEVICES") = reqdev;
- }
+ reqdev.rewind();
+ req("OUTPUT_DEVICES") = reqdev;
// Save output definition
PlotMod::Instance().OutputFormat(req);
diff --git a/src/uPlot/PlotAction.cc b/src/uPlot/PlotAction.cc
index 7b69d14..693853c 100644
--- a/src/uPlot/PlotAction.cc
+++ b/src/uPlot/PlotAction.cc
@@ -1,6 +1,6 @@
/***************************** LICENSE START ***********************************
- Copyright 2012 ECMWF and INPE. This software is distributed under the terms
+ Copyright 2015 ECMWF and INPE. This software is distributed under the terms
of the Apache License version 2.0. In applying this license, ECMWF does not
waive the privileges and immunities granted to it by virtue of its status as
an Intergovernmental Organization or submit itself to any jurisdiction.
@@ -24,8 +24,8 @@ PlotAction::Execute (PmContext& context)
MvRequest fullReq = PlotMod::Instance().OutputFormat()("OUTPUT_DEVICES");
fullReq = fullReq + context.InRequest();
- // Call MAGICS to process the request
- MagPlusService::Instance().SetRequest ( fullReq );
+ // Call MAGICS to process the request
+ MagPlusService::Instance().SetRequest ( fullReq );
// Check where to send the output file(s)
// If destination is file, Magics has already created it
diff --git a/test/macros/CMakeLists.txt b/test/macros/CMakeLists.txt
index 674f7c8..98cd93e 100644
--- a/test/macros/CMakeLists.txt
+++ b/test/macros/CMakeLists.txt
@@ -48,6 +48,12 @@ metview_macro_test(MACRO interpolation_sh_to_rgg.mv
metview_macro_test(MACRO interpolation_rgg_to_ll.mv
RESOURCES ../data/t1000_SH_to_RGG_reference.grb)
+metview_macro_test(MACRO grib_to_gpt.mv
+ RESOURCES ../data/t1000_LL_2x2.grb)
+
+metview_macro_test(MACRO gpt_to_grib_with_template.mv
+ RESOURCES ../data/t1000_LL_2x2.grb ../data/10U_GG.grb)
+
#metview_macro_test(MACRO interpolation_gg_to_subarea.mv
# RESOURCES ../data/10U_GG.grb)
diff --git a/test/macros/gpt_to_grib_with_template.mv b/test/macros/gpt_to_grib_with_template.mv
new file mode 100644
index 0000000..0d56430
--- /dev/null
+++ b/test/macros/gpt_to_grib_with_template.mv
@@ -0,0 +1,98 @@
+# Metview Macro
+
+# **************************** LICENSE START ***********************************
+#
+# Copyright 2015 ECMWF. This software is distributed under the terms
+# of the Apache License version 2.0. In applying this license, ECMWF does not
+# waive the privileges and immunities granted to it by virtue of its status as
+# an Intergovernmental Organization or submit itself to any jurisdiction.
+#
+# ***************************** LICENSE END ************************************
+
+
+# Overview: take an existing GRIB file, convert to geopoints, then
+# convert to GRIB using another GRIB field as a template for the
+# output grid.
+
+
+# read the original data file
+grib = read('t1000_LL_2x2.grb')
+
+
+# convert to geopoints
+gpt = grib_to_geo
+(
+ data : grib
+)
+
+# ----------------------------------------------------------------------
+# TEST 1 - convert to a regular Gaussian grid and check the co-ordinates
+# of the output field
+# ----------------------------------------------------------------------
+
+# read the GRIB whose grid we will use as a template for the result
+template_grib = read('10U_GG.grb')
+
+
+# convert to GRIB using the template
+new_grib = geo_to_grib(
+ grid_definition_mode : "grib",
+ tolerance : 1,
+ geopoints : gpt,
+ template_grib : template_grib
+ )
+
+
+# write the result so that we can manually compare it if needed
+write ('gpt_to_grib_result_1.grb', new_grib)
+
+
+# now check whether the geopoints co-ordinates are (almost) the same as in the GRIB
+check_column('latitudes 1', latitudes(gpt), latitudes(grib))
+check_column('longitudes 1', longitudes(gpt), longitudes(grib))
+
+
+# ----------------------------------------------------------------------
+# TEST 2 - convert to the same as the original input GRIB so that we
+# can compare more thoroughly
+# ----------------------------------------------------------------------
+
+# read the GRIB whose grid we will use as a template for the result
+template_grib = read('t1000_LL_2x2.grb')
+
+
+# convert to GRIB using the template
+new_grib = geo_to_grib(
+ grid_definition_mode : "grib",
+ tolerance : 1,
+ geopoints : gpt,
+ template_grib : template_grib
+ )
+
+
+# write the result so that we can manually compare it if needed
+write ('gpt_to_grib_result_2.grb', new_grib)
+
+
+# now check whether the geopoints co-ordinates are (almost) the same as in the GRIB
+check_column('latitudes 2', latitudes(gpt), latitudes(grib))
+check_column('longitudes 2', longitudes(gpt), longitudes(grib))
+check_column('values 2', values(gpt), values(grib))
+
+
+
+# compare a column (e.g. 'values') and fail if the difference is
+# a large enough proportion of the max value of the original GRIB
+
+function check_column(colname: string, gpt_data:vector, grb_data:vector)
+
+ maxdiff = maxvalue(abs(gpt_data-grb_data))
+ maxdiff_as_proportion_of_max = maxdiff / maxvalue(abs(grb_data))
+ print(colname, ' maxdiff: ', maxdiff, ' prop: ', maxdiff_as_proportion_of_max)
+
+ if maxdiff > 0.001 then # more then 0.1 percent difference
+ msg = colname & " has proportional max difference of " & maxdiff_as_proportion_of_max
+ fail(msg)
+ end if
+
+end check_column
diff --git a/test/macros/grib_to_gpt.mv b/test/macros/grib_to_gpt.mv
new file mode 100644
index 0000000..afdddde
--- /dev/null
+++ b/test/macros/grib_to_gpt.mv
@@ -0,0 +1,49 @@
+# Metview Macro
+
+# **************************** LICENSE START ***********************************
+#
+# Copyright 2015 ECMWF. This software is distributed under the terms
+# of the Apache License version 2.0. In applying this license, ECMWF does not
+# waive the privileges and immunities granted to it by virtue of its status as
+# an Intergovernmental Organization or submit itself to any jurisdiction.
+#
+# ***************************** LICENSE END ************************************
+
+
+# read the original data file
+grib = read('t1000_LL_2x2.grb')
+
+
+# convert to geopoints
+gpt = grib_to_geo
+(
+ data : grib
+)
+
+
+# write the result so that we can manually compare it if needed
+write ('grib_to_gpt_result.gpt', gpt)
+
+
+# now check whether the geopoints points are (almost) the same as the GRIB values
+check_column('latitudes ', latitudes(gpt), latitudes(grib))
+check_column('longitudes', longitudes(gpt), longitudes(grib))
+check_column('values ', values(gpt), values(grib))
+
+
+
+# compare a column (e.g. 'values') and fail if the difference is
+# a large enough proportion of the max value of the original GRIB
+
+function check_column(colname: string, gpt_data:vector, grb_data:vector)
+
+ maxdiff = maxvalue(abs(gpt_data-grb_data))
+ maxdiff_as_proportion_of_max = maxdiff / maxvalue(abs(grb_data))
+ print(colname, ' maxdiff: ', maxdiff, ' prop: ', maxdiff_as_proportion_of_max)
+
+ if maxdiff > 0.001 then # more then 0.1 percent difference
+ msg = colname & " has proportional max difference of " & maxdiff_as_proportion_of_max
+ fail(msg)
+ end if
+
+end check_column
--
Alioth's /usr/local/bin/git-commit-notice on /srv/git.debian.org/git/debian-science/packages/metview.git
More information about the debian-science-commits
mailing list