Skip to content

Commit

Permalink
Merge branch 'hotfix/0.19.1'
Browse files Browse the repository at this point in the history 
* hotfix/0.19.1:
  Version 0.19.1
  ATLAS-256, PGEN-274: fix different double array hashes and west=0 (different cache files, requires Legendre cache version increment)
  ATLAS-255: regression test fix
  ATLAS-256, PGEN-274: unit tests for both 'local' and 'ifs' spectral transforms cache file naming (more tests)
  ATLAS-256, PGEN-274: unit tests for both 'local' and 'ifs' spectral transforms cache file naming
  ATLAS-255: regression test fix
  ATLAS-255: regression test
  • Loading branch information
@wdeconinck
wdeconinck committed Dec 19, 2019
2 parents e058524 + 69a5c1f commit aae9187
Show file tree
Hide file tree
Showing 8 changed files with 343 additions and 118 deletions.
5 changes: 5 additions & 0 deletions CHANGELOG.md
View file
Original file line number Diff line number Diff line change
Expand Up @@ -8,6 +8,11 @@ This project adheres to [Semantic Versioning](http://semver.org/spec/v2.0.0.html
## [Unreleased]


## [0.19.1] - 2019-12-19
### Fixed
- Keep Gaussian identity of a Gaussian grid if a given domain does not crop any latitudes
- Fix naming for LegendreCache, to be more specific, and platform independent

## [0.19.0] - 2019-10-01
### Fixed
- Lambert ( conformal conic ) projection xy coordinates are now corrected
Expand Down
2 changes: 1 addition & 1 deletion VERSION.cmake
View file
Original file line number Diff line number Diff line change
Expand Up @@ -6,5 +6,5 @@
# granted to it by virtue of its status as an intergovernmental organisation nor
# does it submit to any jurisdiction.

set ( ${PROJECT_NAME}_VERSION_STR "0.19.0" )
set ( ${PROJECT_NAME}_VERSION_STR "0.19.1" )

191 changes: 95 additions & 96 deletions src/atlas/grid/detail/grid/Structured.cc
View file
Original file line number Diff line number Diff line change
Expand Up @@ -357,115 +357,114 @@ class Normalise {
} // namespace

void Structured::crop( const Domain& dom ) {
if ( dom ) {
ATLAS_ASSERT( dom.units() == projection().units() );
if ( !dom ) {
domain_ =
RectangularDomain( {xspace_.min(), xspace_.max()}, {yspace_.min(), yspace_.max()}, projection_.units() );
return;
}

ATLAS_ASSERT( dom.units() == projection().units() );
auto rect_domain = RectangularDomain( dom );

auto rect_domain = RectangularDomain( dom );
if ( !rect_domain ) {
std::stringstream errmsg;
errmsg << "Cannot crop the grid with domain " << dom;
throw_Exception( errmsg.str(), Here() );
}

if ( rect_domain ) {
const double cropped_ymin = rect_domain.ymin();
const double cropped_ymax = rect_domain.ymax();
const double cropped_ymin = rect_domain.ymin();
const double cropped_ymax = rect_domain.ymax();

// Cropping in Y
idx_t jmin = ny();
idx_t jmax = 0;
for ( idx_t j = 0; j < ny(); ++j ) {
if ( rect_domain.contains_y( y( j ) ) ) {
jmin = std::min( j, jmin );
jmax = std::max( j, jmax );
}
}
ATLAS_ASSERT( jmax >= jmin );

idx_t cropped_ny = jmax - jmin + 1;
std::vector<double> cropped_y( y_.begin() + jmin, y_.begin() + jmin + cropped_ny );
std::vector<double> cropped_dx( dx_.begin() + jmin, dx_.begin() + jmin + cropped_ny );

std::vector<double> cropped_xmin( cropped_ny, std::numeric_limits<double>::max() );
std::vector<double> cropped_xmax( cropped_ny, -std::numeric_limits<double>::max() );
std::vector<idx_t> cropped_nx( cropped_ny );

// Cropping in X
Normalise normalise( rect_domain );
for ( idx_t j = jmin, jcropped = 0; j <= jmax; ++j, ++jcropped ) {
idx_t n = 0;
for ( idx_t i = 0; i < nx( j ); ++i ) {
const double _x = normalise( x( i, j ) );
if ( rect_domain.contains_x( _x ) ) {
cropped_xmin[jcropped] = std::min( cropped_xmin[jcropped], _x );
cropped_xmax[jcropped] = std::max( cropped_xmax[jcropped], _x );
++n;
}
}
cropped_nx[jcropped] = n;
// Cropping in Y
idx_t jmin = ny();
idx_t jmax = 0;
for ( idx_t j = 0; j < ny(); ++j ) {
if ( rect_domain.contains_y( y( j ) ) ) {
jmin = std::min( j, jmin );
jmax = std::max( j, jmax );
}
}
ATLAS_ASSERT( jmax >= jmin );

idx_t cropped_ny = jmax - jmin + 1;
ATLAS_ASSERT( cropped_ny <= ny() );

std::vector<double> cropped_y( y_.begin() + jmin, y_.begin() + jmin + cropped_ny );
std::vector<double> cropped_dx( dx_.begin() + jmin, dx_.begin() + jmin + cropped_ny );

std::vector<double> cropped_xmin( cropped_ny, std::numeric_limits<double>::max() );
std::vector<double> cropped_xmax( cropped_ny, -std::numeric_limits<double>::max() );
std::vector<idx_t> cropped_nx( cropped_ny );

// Cropping in X
Normalise normalise( rect_domain );
for ( idx_t j = jmin, jcropped = 0; j <= jmax; ++j, ++jcropped ) {
idx_t n = 0;
for ( idx_t i = 0; i < nx( j ); ++i ) {
const double _x = normalise( x( i, j ) );
if ( rect_domain.contains_x( _x ) ) {
cropped_xmin[jcropped] = std::min( cropped_xmin[jcropped], _x );
cropped_xmax[jcropped] = std::max( cropped_xmax[jcropped], _x );
++n;
}
bool endpoint = true;
if ( ZonalBandDomain( rect_domain ) ) {
for ( idx_t j = 0; j < cropped_ny; ++j ) {
if ( eckit::types::is_approximately_equal( cropped_xmax[j] + cropped_dx[j], cropped_xmin[j] + 360.,
1.e-10 ) ) {
cropped_xmax[j] = cropped_xmin[j] + 360.;
endpoint = false;
}
}
}
cropped_nx[jcropped] = n;
}
bool endpoint = true;
if ( ZonalBandDomain( rect_domain ) ) {
for ( idx_t j = 0; j < cropped_ny; ++j ) {
if ( eckit::types::is_approximately_equal( cropped_xmax[j] + cropped_dx[j], cropped_xmin[j] + 360.,
1.e-10 ) ) {
cropped_xmax[j] = cropped_xmin[j] + 360.;
endpoint = false;
}
}
}

// Complete structures
// Complete structures

idx_t cropped_nxmin, cropped_nxmax;
cropped_nxmin = cropped_nxmax = cropped_nx.front();
idx_t cropped_nxmin, cropped_nxmax;
cropped_nxmin = cropped_nxmax = cropped_nx.front();

for ( idx_t j = 1; j < cropped_ny; ++j ) {
cropped_nxmin = std::min( cropped_nx[j], cropped_nxmin );
cropped_nxmax = std::max( cropped_nx[j], cropped_nxmax );
}
idx_t cropped_npts = std::accumulate( cropped_nx.begin(), cropped_nx.end(), idx_t{0} );

Spacing cropped_yspace(
new spacing::CustomSpacing( cropped_ny, cropped_y.data(), {cropped_ymin, cropped_ymax} ) );
for ( idx_t j = 1; j < cropped_ny; ++j ) {
cropped_nxmin = std::min( cropped_nx[j], cropped_nxmin );
cropped_nxmax = std::max( cropped_nx[j], cropped_nxmax );
}
idx_t cropped_npts = std::accumulate( cropped_nx.begin(), cropped_nx.end(), idx_t{0} );

std::vector<Spacing> cropped_xspacings( cropped_ny );
for ( idx_t j = 0; j < cropped_ny; ++j ) {
cropped_xspacings[j] =
new spacing::LinearSpacing( cropped_xmin[j], cropped_xmax[j], cropped_nx[j], endpoint );
}
XSpace cropped_xspace( cropped_xspacings );

for ( idx_t j = 0; j < cropped_ny; ++j ) {
ATLAS_ASSERT( eckit::types::is_approximately_equal( cropped_xspace.xmin()[j], cropped_xmin[j] ) );
if ( cropped_nx[j] > 1 ) {
ATLAS_ASSERT(
eckit::types::is_approximately_equal( cropped_xspace.dx()[j], cropped_dx[j], 1.e-10 ) );
}
ATLAS_ASSERT( cropped_xspace.nx()[j] == cropped_nx[j] );
ATLAS_ASSERT( cropped_xspace.nxmin() == cropped_nxmin );
ATLAS_ASSERT( cropped_xspace.nxmax() == cropped_nxmax );
}
std::vector<Spacing> cropped_xspacings( cropped_ny );
for ( idx_t j = 0; j < cropped_ny; ++j ) {
cropped_xspacings[j] = new spacing::LinearSpacing( cropped_xmin[j], cropped_xmax[j], cropped_nx[j], endpoint );
}
XSpace cropped_xspace( cropped_xspacings );

// Modify grid
{
yspace_ = cropped_yspace;
xspace_ = cropped_xspace;
xmin_ = cropped_xmin;
xmax_ = cropped_xmax;
dx_ = cropped_dx;
nx_ = cropped_nx;
nxmin_ = cropped_nxmin;
nxmax_ = cropped_nxmax;
npts_ = cropped_npts;
y_ = cropped_y;
}
}
else {
std::stringstream errmsg;
errmsg << "Cannot crop the grid with domain " << dom;
throw_Exception( errmsg.str(), Here() );
for ( idx_t j = 0; j < cropped_ny; ++j ) {
ATLAS_ASSERT( eckit::types::is_approximately_equal( cropped_xspace.xmin()[j], cropped_xmin[j] ) );
if ( cropped_nx[j] > 1 ) {
ATLAS_ASSERT( eckit::types::is_approximately_equal( cropped_xspace.dx()[j], cropped_dx[j], 1.e-10 ) );
}
ATLAS_ASSERT( cropped_xspace.nx()[j] == cropped_nx[j] );
ATLAS_ASSERT( cropped_xspace.nxmin() == cropped_nxmin );
ATLAS_ASSERT( cropped_xspace.nxmax() == cropped_nxmax );
}
domain_ =
dom ? dom
: RectangularDomain( {xspace_.min(), xspace_.max()}, {yspace_.min(), yspace_.max()}, projection_.units() );

// Modify grid
if (ny() != cropped_ny) {
// keep specialised spacing (Gaussian, Linear, ...) unless cropping happens
yspace_ = new spacing::CustomSpacing( cropped_ny, cropped_y.data(), {cropped_ymin, cropped_ymax} );
}

xspace_ = cropped_xspace;
xmin_ = cropped_xmin;
xmax_ = cropped_xmax;
dx_ = cropped_dx;
nx_ = cropped_nx;
nxmin_ = cropped_nxmin;
nxmax_ = cropped_nxmax;
npts_ = cropped_npts;
y_ = cropped_y;

domain_ = dom;
}

void Structured::computeTruePeriodicity() {
Expand Down
37 changes: 29 additions & 8 deletions src/atlas/trans/ifs/LegendreCacheCreatorIFS.cc
View file
Original file line number Diff line number Diff line change
Expand Up @@ -9,12 +9,17 @@
*/

#include "atlas/trans/ifs/LegendreCacheCreatorIFS.h"

#include <cmath>
#include <sstream>
#include <string>

#include "eckit/utils/MD5.h"

#include "atlas/grid.h"
#include "atlas/option.h"
#include "atlas/runtime/Exception.h"
#include "atlas/trans/Trans.h"
#include "eckit/utils/MD5.h"

namespace atlas {
namespace trans {
Expand All @@ -32,16 +37,32 @@ std::string truncate( const std::string& str ) {

std::string hash( const Grid& grid ) {
eckit::MD5 h;
if ( StructuredGrid( grid ) && not grid.projection() ) {
auto g = StructuredGrid( grid );
h.add( g.y().data(), g.y().size() * sizeof( double ) );

StructuredGrid structured( grid );
if ( structured && not grid.projection() ) {
for ( auto& y : structured.y() ) {
h.add( std::lround( y * 1.e8 ) );
}
}
else {
grid.hash( h );
}
return truncate( h.digest() );
}

std::string hash_pl( const Grid& grid ) {
eckit::MD5 h;

StructuredGrid structured( grid );
ATLAS_ASSERT( structured );

for ( auto& n : structured.nx() ) {
h.add( long( n ) );
}

return truncate( h.digest() );
}

std::string hash( const eckit::Configuration& config ) {
eckit::MD5 h;

Expand All @@ -57,13 +78,13 @@ std::string LegendreCacheCreatorIFS::uid() const {
if ( unique_identifier_.empty() ) {
std::ostringstream stream;
stream << "ifs-T" << truncation_ << "-";
if ( GaussianGrid( grid_ ) ) {
GaussianGrid gaussian( grid_ );
if ( gaussian ) {
if ( RegularGaussianGrid( grid_ ) ) {
stream << "RegularGaussianN" << GaussianGrid( grid_ ).N();
stream << "RegularGaussianN" << gaussian.N();
}
else {
stream << "ReducedGaussianN" << GaussianGrid( grid_ ).N() << "-PL";
stream << hash( grid_ );
stream << "ReducedGaussianN" << gaussian.N() << "-PL" << hash_pl( grid_ );
}
}
else if ( RegularLonLatGrid( grid_ ) ) {
Expand Down
17 changes: 12 additions & 5 deletions src/atlas/trans/local/LegendreCacheCreatorLocal.cc
View file
Original file line number Diff line number Diff line change
Expand Up @@ -10,6 +10,7 @@

#include "atlas/trans/local/LegendreCacheCreatorLocal.h"

#include <cmath>
#include <sstream>
#include <string>

Expand Down Expand Up @@ -38,9 +39,12 @@ std::string truncate( const std::string& str ) {

std::string hash( const Grid& grid ) {
eckit::MD5 h;
if ( StructuredGrid( grid ) && not grid.projection() ) {
auto g = StructuredGrid( grid );
h.add( g.y().data(), g.y().size() * sizeof( double ) );

StructuredGrid structured( grid );
if ( structured && not grid.projection() ) {
for ( auto& y : structured.y() ) {
h.add( std::lround( y * 1.e8 ) );
}
}
else {
grid.hash( h );
Expand Down Expand Up @@ -68,7 +72,10 @@ std::string LegendreCacheCreatorLocal::uid() const {
};
stream << "local-T" << truncation_ << "-";
StructuredGrid structured( grid_ );
if ( GaussianGrid( grid_ ) ) {
if ( grid_.projection() ) {
give_up();
}
else if ( GaussianGrid( grid_ ) ) {
// Same cache for any global Gaussian grid
stream << "GaussianN" << GaussianGrid( grid_ ).N();
}
Expand All @@ -94,7 +101,7 @@ std::string LegendreCacheCreatorLocal::uid() const {
give_up();
}
}
else if ( RegularGrid( grid_ ) && not grid_.projection() && structured.yspace().type() == "linear" ) {
else if ( RegularGrid( grid_ ) && structured.yspace().type() == "linear" ) {
RectangularDomain domain( grid_.domain() );
ATLAS_ASSERT( domain );
stream << "Regional";
Expand Down
11 changes: 11 additions & 0 deletions src/tests/grid/test_grids.cc
View file
Original file line number Diff line number Diff line change
Expand Up @@ -196,6 +196,17 @@ CASE( "Create unstructured from unstructured" ) {
EXPECT( ugrid.size() == StructuredGrid( agrid, domain ).size() );
}

CASE( "ATLAS-255: regular Gaussian grid with global domain" ) {
GlobalDomain globe;
Grid grid( "F80", globe );
EXPECT( GaussianGrid( grid ) );
}

CASE( "ATLAS-255: reduced Gaussian grid with global domain" ) {
GlobalDomain globe;
Grid grid( "O80", globe );
EXPECT( GaussianGrid( grid ) );
}

//-----------------------------------------------------------------------------

Expand Down
Loading

0 comments on commit aae9187

Please sign in to comment.