Merge branch 'release/0.14.0'

* release/0.14.0:
  Update travis
  Update changelog
  Require eckit 0.20.0
  Version 0.14.0
  ATLAS-134: TransCacheMemoryEntry
  eckit::geometry::Sphere update (Willem, please confirm)

.travis.yml

@@ -196,7 +196,7 @@ install:
     ### Install eckit
     install-dep.sh --repo eckit --branch ${DEPS_BRANCH} --prefix ${DEPS_DIR}/eckit --cmake "-DENABLE_TESTS=OFF -DCMAKE_BUILD_TYPE=DEBUG ${ECKIT_CMAKE_OPTIONS}"
   - export ECKIT_PATH=${DEPS_DIR}/eckit
-  - ${DEPS_DIR}/eckit/bin/eckit_version
+  - ${DEPS_DIR}/eckit/bin/eckit-version
 
   #################################################################
   # Install fckit

CHANGELOG.md

@@ -7,6 +7,14 @@ This project adheres to [Semantic Versioning](http://semver.org/spec/v2.0.0.html
 
 ## [Unreleased]
 
+## [0.14.0] - 2018-03-22
+### Added
+- Spherical Harmonics transforms can receive a cache memory handle
+
+### Changed
+- Earth interface (C++)
+- Requires eckit 0.20.0, fckit 0.5.0
+
 ## [0.13.2] - 2018-03-20
 ### Fixed
 - C++ compilation using PGI 17, 18 and GCC 7
@@ -22,5 +30,6 @@ This project adheres to [Semantic Versioning](http://semver.org/spec/v2.0.0.html
 ## 0.13.0 - 2018-02-16
 
 [Unreleased]: https://github.com/ecmwf/atlas/compare/master...develop
+[0.14.0]: https://github.com/ecmwf/atlas/compare/0.13.2...0.14.0
 [0.13.2]: https://github.com/ecmwf/atlas/compare/0.13.1...0.13.2
 [0.13.1]: https://github.com/ecmwf/atlas/compare/0.13.0...0.13.1

CMakeLists.txt

@@ -31,7 +31,7 @@ ecbuild_declare_project()
 
 ### eckit
 
-ecbuild_use_package( PROJECT eckit VERSION 0.18.1 REQUIRED )
+ecbuild_use_package( PROJECT eckit VERSION 0.20.0 REQUIRED )
 ecbuild_debug( "   ECKIT_FEATURES : [${ECKIT_FEATURES}]" )
 
 # options & dependencies

VERSION.cmake

@@ -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.13.2" )
+set  ( ${PROJECT_NAME}_VERSION_STR  "0.14.0" )
 

src/apps/atlas-benchmark.cc

@@ -279,7 +279,7 @@ void AtlasBenchmark::execute( const Args& args ) {
 //----------------------------------------------------------------------------------------------------------------------
 
 void AtlasBenchmark::initial_condition( const Field& field, const double& beta ) {
-    const double radius  = util::Earth::radiusInMeters();
+    const double radius  = util::Earth::radius();
     const double USCAL   = 20.;
     const double pvel    = USCAL / radius;
     const double deg2rad = M_PI / 180.;

src/atlas/CMakeLists.txt

@@ -452,7 +452,6 @@ runtime/ErrorHandling.h
 util/Config.cc
 util/Config.h
 util/Constants.h
-util/Earth.cc
 util/Earth.h
 util/GaussianLatitudes.cc
 util/GaussianLatitudes.h
@@ -466,6 +465,7 @@ util/Rotation.cc
 util/Rotation.h
 util/SphericalPolygon.cc
 util/SphericalPolygon.h
+util/UnitSphere.h
 util/detail/BlackMagic.h
 util/detail/Cache.h
 util/detail/Debug.h

src/atlas/grid/detail/grid/Structured.cc

@@ -22,8 +22,8 @@
 #include "atlas/grid/detail/spacing/LinearSpacing.h"
 #include "atlas/runtime/ErrorHandling.h"
 #include "atlas/runtime/Log.h"
-#include "atlas/util/Earth.h"
 #include "atlas/util/Point.h"
+#include "atlas/util/UnitSphere.h"
 
 namespace atlas {
 namespace grid {
@@ -438,8 +438,11 @@ void Structured::computeTruePeriodicity() {
             const PointLonLat Pllmin = projection().lonlat( PointXY( xmin_[j], y_[j] ) );
             const PointLonLat Pllmax = projection().lonlat( PointXY( xmax_[j], y_[j] ) );
 
-            Point3 Pxmin = util::Earth::convertGeodeticToGeocentric( Pllmin, 1. );
-            Point3 Pxmax = util::Earth::convertGeodeticToGeocentric( Pllmax, 1. );
+            Point3 Pxmin;
+            util::UnitSphere::convertSphericalToCartesian(Pllmin, Pxmin );
+
+            Point3 Pxmax;
+            util::UnitSphere::convertSphericalToCartesian(Pllmax, Pxmax );
 
             periodic_x_ = points_equal( Pxmin, Pxmax );
         }

src/atlas/grid/detail/grid/Structured.h

@@ -140,8 +140,8 @@ class Structured : public Grid {
         }
 
     private:
-        Grid::IteratorXY::Predicate p_;
         const Structured& grid_;
+        Grid::IteratorXY::Predicate p_;
         size_t i_;
         size_t j_;
         size_t n_;

src/atlas/grid/detail/grid/Unstructured.h

@@ -112,8 +112,8 @@ class Unstructured : public Grid {
     private:
         const Unstructured& grid_;
         Grid::IteratorXY::Predicate p_;
-        size_t n_;
         size_t size_;
+        size_t n_;
     };
 
     class IteratorLonLat : public Grid::IteratorLonLat {

src/atlas/interpolation/method/FiniteElement.cc

@@ -65,7 +65,7 @@ void FiniteElement::setup( const FunctionSpace& source, const FunctionSpace& tar
         PointXYZ p2;
         for ( size_t n = 0; n < N; ++n ) {
             const PointLonLat p1( lonlat( n, 0 ), lonlat( n, 1 ) );
-            p2          = util::Earth::convertGeodeticToGeocentric( p1 );
+            util::Earth::convertSphericalToCartesian( p1, p2 );
             xyz( n, 0 ) = p2.x();
             xyz( n, 1 ) = p2.y();
             xyz( n, 2 ) = p2.z();
@@ -148,7 +148,8 @@ void FiniteElement::setup( const FunctionSpace& source ) {
                 failures.push_back( ip );
                 Log::debug() << "------------------------------------------------------"
                                 "---------------------\n";
-                const PointLonLat pll = util::Earth::convertGeocentricToGeodetic( p );
+                PointLonLat pll;
+                util::Earth::convertCartesianToSpherical( p, pll );
                 Log::debug() << "Failed to project point (lon,lat)=" << pll << '\n';
                 Log::debug() << failures_log.str();
             }
@@ -163,7 +164,8 @@ void FiniteElement::setup( const FunctionSpace& source ) {
         msg << "Rank " << eckit::mpi::comm().rank() << " failed to project points:\n";
         for ( std::vector<size_t>::const_iterator i = failures.begin(); i != failures.end(); ++i ) {
             const PointXYZ p{( *ocoords_ )( *i, 0 ), ( *ocoords_ )( *i, 1 ), ( *ocoords_ )( *i, 2 )};  // lookup point
-            const PointLonLat pll = util::Earth::convertGeocentricToGeodetic( p );
+            PointLonLat pll;
+            util::Earth::convertCartesianToSpherical(p, pll);
             msg << "\t(lon,lat) = " << pll << "\n";
         }
 

src/atlas/mesh/actions/BuildStatistics.cc

@@ -29,9 +29,11 @@
 #include "atlas/mesh/actions/BuildDualMesh.h"
 #include "atlas/parallel/Checksum.h"
 #include "atlas/runtime/ErrorHandling.h"
+#include "atlas/util/Constants.h"
 #include "atlas/util/CoordinateEnums.h"
 #include "atlas/util/Earth.h"
 #include "atlas/util/Point.h"
+#include "atlas/util/UnitSphere.h"
 
 namespace atlas {
 namespace mesh {
@@ -62,9 +64,11 @@ void quad_quality( double& eta, double& rho, const PointLonLat& p1, const PointL
                    const PointLonLat& p4 ) {
     // see http://geuz.org/gmsh/doc/preprints/gmsh_quad_preprint.pdf
 
-    const PointXYZ xyz[]{
-        util::Earth::convertGeodeticToGeocentric( p1, 1. ), util::Earth::convertGeodeticToGeocentric( p2, 1. ),
-        util::Earth::convertGeodeticToGeocentric( p3, 1. ), util::Earth::convertGeodeticToGeocentric( p4, 1. )};
+    PointXYZ xyz[4];
+    util::UnitSphere::convertSphericalToCartesian(p1, xyz[0]);
+    util::UnitSphere::convertSphericalToCartesian(p2, xyz[1]);
+    util::UnitSphere::convertSphericalToCartesian(p3, xyz[2]);
+    util::UnitSphere::convertSphericalToCartesian(p4, xyz[3]);
 
     PointXYZ l2m1( PointXYZ::sub( xyz[1], xyz[0] ) );
     PointXYZ l3m2( PointXYZ::sub( xyz[2], xyz[1] ) );
@@ -123,7 +127,7 @@ void build_statistics( Mesh& mesh ) {
             int ip2 = edge_nodes( jedge, 1 );
             PointLonLat p1( lonlat( ip1, LON ), lonlat( ip1, LAT ) );
             PointLonLat p2( lonlat( ip2, LON ), lonlat( ip2, LAT ) );
-            dist( jedge ) = util::Earth::distanceInMeters( p1, p2 ) * 1e-3;
+            dist( jedge ) = util::Earth::distance( p1, p2 ) * 1e-3;
         }
     }
 
@@ -209,8 +213,8 @@ void build_statistics( Mesh& mesh ) {
 
         for ( size_t jnode = 0; jnode < nodes.size(); ++jnode ) {
             const double lat        = util::Constants::degreesToRadians() * lonlat( jnode, LAT );
-            const double hx         = util::Constants::degreesToRadians() * util::Earth::radiusInKm() * std::cos( lat );
-            const double hy         = util::Constants::degreesToRadians() * util::Earth::radiusInKm();
+            const double hx         = util::Constants::degreesToRadians() * util::Earth::radius() * std::cos( lat );
+            const double hy         = util::Constants::degreesToRadians() * util::Earth::radius();
             dual_delta_sph( jnode ) = std::sqrt( dual_volumes( jnode ) * hx * hy );
         }
 

src/atlas/mesh/actions/BuildXYZField.cc

@@ -44,7 +44,7 @@ Field& BuildXYZField::operator()( mesh::Nodes& nodes ) const {
         PointXYZ p2;
         for ( size_t n = 0; n < nodes.size(); ++n ) {
             const PointLonLat p1( lonlat( n, 0 ), lonlat( n, 1 ) );
-            p2          = util::Earth::convertGeodeticToGeocentric( p1 );
+            util::Earth::convertSphericalToCartesian( p1, p2 );
             xyz( n, 0 ) = p2.x();
             xyz( n, 1 ) = p2.y();
             xyz( n, 2 ) = p2.z();

src/atlas/mesh/actions/ExtendNodesGlobal.cc

@@ -8,13 +8,13 @@
  * nor does it submit to any jurisdiction.
  */
 
-#include "eckit/exception/Exceptions.h"
+#include "atlas/mesh/actions/ExtendNodesGlobal.h"
 
+#include "eckit/exception/Exceptions.h"
 #include "atlas/field/Field.h"
 #include "atlas/grid/Grid.h"
 #include "atlas/mesh/Mesh.h"
 #include "atlas/mesh/Nodes.h"
-#include "atlas/mesh/actions/ExtendNodesGlobal.h"
 #include "atlas/util/CoordinateEnums.h"
 #include "atlas/util/Earth.h"
 
@@ -64,7 +64,10 @@ void ExtendNodesGlobal::operator()( const Grid& grid, Mesh& mesh ) const {
     for ( size_t i = 0; i < nb_extension_pts; ++i ) {
         const size_t n        = nb_real_pts + i;
         const PointLonLat pLL = grid.projection().lonlat( extended_pts[i] );
-        const PointXYZ pXYZ   = util::Earth::convertGeodeticToGeocentric( pLL );
+
+        PointXYZ pXYZ;
+        util::Earth::convertSphericalToCartesian( pLL, pXYZ );
+
         xyz( n, XX )          = pXYZ.x();
         xyz( n, YY )          = pXYZ.y();
         xyz( n, ZZ )          = pXYZ.z();

src/atlas/numerics/fvm/Method.cc

@@ -51,7 +51,7 @@ size_t get_levels( const eckit::Parametrisation& params ) {
 }
 
 double get_radius( const eckit::Parametrisation& params ) {
-    double radius = util::Earth::radiusInMeters();
+    double radius = util::Earth::radius();
     params.get( "radius", radius );
     return radius;
 }

src/atlas/option/Options.cc

@@ -58,7 +58,7 @@ radius::radius( double _radius ) {
 }
 
 radius::radius( const std::string& key ) {
-    if ( key == "Earth" ) { set( "radius", util::Earth::radiusInMeters() ); }
+    if ( key == "Earth" ) { set( "radius", util::Earth::radius() ); }
     else {
         NOTIMP;
     }

src/atlas/projection/detail/LambertProjection.cc

@@ -38,7 +38,7 @@ namespace detail {
 // constructors
 LambertProjection::LambertProjection( const eckit::Parametrisation& params ) {
     // check presence of radius
-    if ( !params.get( "radius", radius_ ) ) radius_ = util::Earth::radiusInMeters();
+    if ( !params.get( "radius", radius_ ) ) radius_ = util::Earth::radius();
     // check presence of lat1 and lat2
     if ( !params.get( "latitude1", lat1_ ) ) throw eckit::BadParameter( "latitude1 missing in Params", Here() );
     if ( !params.get( "latitude2", lat2_ ) ) lat2_ = lat1_;
@@ -95,7 +95,7 @@ LambertProjection::Spec LambertProjection::spec() const {
     proj_spec.set( "latitude1", lat1_ );
     proj_spec.set( "latitude2", lat2_ );
     proj_spec.set( "longitude0", lon0_ );
-    if ( radius_ != util::Earth::radiusInMeters() ) proj_spec.set( "radius", radius_ );
+    if ( radius_ != util::Earth::radius() ) proj_spec.set( "radius", radius_ );
     return proj_spec;
 }
 

src/atlas/projection/detail/MercatorProjection.cc

@@ -41,7 +41,7 @@ MercatorProjectionT<Rotation>::MercatorProjectionT( const eckit::Parametrisation
     ProjectionImpl(),
     rotation_( params ) {
     // check presence of radius
-    if ( !params.get( "radius", radius_ ) ) radius_ = util::Earth::radiusInMeters();
+    if ( !params.get( "radius", radius_ ) ) radius_ = util::Earth::radius();
     // check presence of lon0
     if ( !params.get( "longitude0", lon0_ ) ) lon0_ = 0.0;
 
@@ -74,7 +74,7 @@ typename MercatorProjectionT<Rotation>::Spec MercatorProjectionT<Rotation>::spec
     Spec proj_spec;
     proj_spec.set( "type", static_type() );
     proj_spec.set( "longitude0", lon0_ );
-    if ( radius_ != util::Earth::radiusInMeters() ) proj_spec.set( "radius", radius_ );
+    if ( radius_ != util::Earth::radius() ) proj_spec.set( "radius", radius_ );
     rotation_.spec( proj_spec );
     return proj_spec;
 }

src/atlas/trans/Trans.h

@@ -44,13 +44,13 @@ class TransCacheEntry {
     virtual const void* data() const = 0;
 };
 
-class EmptyCacheEntry : public TransCacheEntry {
+class EmptyCacheEntry final : public TransCacheEntry {
 public:
     virtual size_t size() const override { return 0; }
     virtual const void* data() const override { return nullptr; }
 };
 
-class TransCacheFileEntry : public TransCacheEntry {
+class TransCacheFileEntry final : public TransCacheEntry {
     eckit::Buffer buffer_;
 
 public:
@@ -64,6 +64,21 @@ class TransCacheFileEntry : public TransCacheEntry {
     virtual const void* data() const override { return buffer_.data(); }
 };
 
+class TransCacheMemoryEntry final : public TransCacheEntry {
+public:
+    TransCacheMemoryEntry(const void* data, size_t size) : data_(data), size_(size) {
+        ASSERT(data_);
+        ASSERT(size_);
+    }
+    virtual const void* data() const override { return data_; }
+    virtual size_t size() const override { return size_; }
+private:
+    const void* data_;
+    const size_t size_;
+};
+
+//-----------------------------------------------------------------------------
+
 class Cache {
 public:
     Cache() : legendre_( new EmptyCacheEntry() ), fft_( new EmptyCacheEntry() ) {}
@@ -86,6 +101,9 @@ class Cache {
 
 class LegendreCache : public Cache {
 public:
+    LegendreCache(const void* address, size_t size) :
+        Cache(std::make_shared<atlas::trans::TransCacheMemoryEntry>(address, size)) {
+    }
     LegendreCache( const eckit::PathName& path ) :
         Cache( std::shared_ptr<TransCacheEntry>( new TransCacheFileEntry( path ) ) ) {}
 };

src/atlas/trans/local/TransLocal.cc

@@ -18,6 +18,7 @@
 #include "atlas/trans/local/FourierTransforms.h"
 #include "atlas/trans/local/LegendrePolynomials.h"
 #include "atlas/trans/local/LegendreTransforms.h"
+#include "atlas/util/Constants.h"
 
 namespace atlas {
 namespace trans {

src/atlas/trans/local/TransLocal.h

@@ -108,8 +108,8 @@ class TransLocal : public trans::TransImpl {
                       const eckit::Configuration& = util::NoConfig() ) const;
 
 private:
-    int truncation_;
     Grid grid_;
+    int truncation_;
     bool precompute_;
     std::vector<double> legendre_;
     std::vector<size_t> legendre_begin_;

src/atlas/trans/local/VorDivToUVLocal.cc

@@ -8,6 +8,7 @@
  * nor does it submit to any jurisdiction.
  */
 
+#include <cmath>  // for std::sqrt
 #include "atlas/trans/local/VorDivToUVLocal.h"
 #include "atlas/functionspace/Spectral.h"
 #include "atlas/runtime/Log.h"
@@ -78,7 +79,7 @@ void vd2uv( const int truncation,               // truncation
     repsnm[0] = 0.;
 
     // rlapin: constant factor from eq.(2.2) and (2.3) in [Temperton 1991]
-    double ra = util::Earth::radiusInMeters();
+    double ra = util::Earth::radius();
     std::vector<double> rlapin( truncation + 3 );
     for ( int jn = 1; jn <= truncation + 2; ++jn ) {
         rlapin[jn] = -ra * ra / ( jn * ( jn + 1. ) );
@@ -145,7 +146,7 @@ void vd2uv( const int truncation,               // truncation
     int ilcm = truncation - km;
     int ioff = ( 2 * truncation - km + 3 ) * km;
     // ioff: start index of zonal wavenumber km in spectral data
-    double za_r = 1. / util::Earth::radiusInMeters();
+    double za_r = 1. / util::Earth::radius();
     for ( int j = 0; j <= ilcm; ++j ) {
         // ilcm-j = total wavenumber
         int inm = ioff + ( ilcm - j ) * 2;