src/atlas/grid/detail/tiles/FV3Tiles.cc
// Calculates the FV3 panel
// Input (crd) is Longitude and Latitude in Radians
// Output is tile number
Calculates the FV3 panel
Input (crd) is Longitude and Latitude in Radians
Output is tile number
Calculates the FV3 panel Input ( crd ) is Longitude and Latitude in Radians Output is tile number
idx_t FV3CubedSphereTiles::indexFromLonLat(const double crd[]) const {
    idx_t t(-1);  
    double xyz[3];
    sphericalToCartesian(crd, xyz);
    static const double cornerLat = std::asin(1. / std::sqrt(3.0)) * rad2deg;
    const double& lon = crd[LON];
    const double& lat = crd[LAT];
    double zPlusAbsX  = xyz[ZZ] + std::abs(xyz[XX]);
    double zPlusAbsY  = xyz[ZZ] + std::abs(xyz[YY]);
    double zMinusAbsX = xyz[ZZ] - std::abs(xyz[XX]);
    double zMinusAbsY = xyz[ZZ] - std::abs(xyz[YY]);
    if (is_tiny(zPlusAbsX)) {
        zPlusAbsX = 0.;
    }
    if (is_tiny(zPlusAbsY)) {
        zPlusAbsY = 0.;
    }
    if (is_tiny(zMinusAbsX)) {
        zMinusAbsX = 0.;
    }
    if (is_tiny(zMinusAbsY)) {
        zMinusAbsY = 0.;
    }
    if (is_geq(lon, 315.) || is_less(lon, 45.)) {
        if ((zPlusAbsX <= 0.) && (zPlusAbsY <= 0.)) {
            t = 2;
        }
        else if ((zMinusAbsX > 0.) && (zMinusAbsY > 0.)) {
            t = 5;
        }
        else {
            t = 0;
        }
        if (is_same(lon, -45.) && is_same(lat, cornerLat)) {
            t = 0;
        }
        if (is_same(lon, 315.) && is_same(lat, cornerLat)) {
            t = 0;
        }
    }
    if (is_geq(lon, 45.) && is_less(lon, 135.)) {
        if ((zPlusAbsX <= 0.) && (zPlusAbsY <= 0.)) {
            t = 2;
        }
        else if ((zMinusAbsX > 0.) && (zMinusAbsY > 0.)) {
            t = 5;
        }
        else {
            t = 1;
        }
    }
    if (is_geq(lon, 135.) && is_less(lon, 225.)) {
        if ((zPlusAbsX < 0.) && (zPlusAbsY < 0.)) {
            t = 2;
        }
        else if ((zMinusAbsX >= 0.) && (zMinusAbsY >= 0.)) {
            t = 5;
        }
        else {
            t = 3;
        }
        if (is_same(lon, 135.) && is_same(lat, -cornerLat)) {
            t = 1;
        }
    }
    if (is_geq(lon, 225.) && is_less(lon, 315.)) {
        if ((zPlusAbsX < 0.) && (zPlusAbsY < 0.)) {
            t = 2;
        }
        else if ((zMinusAbsX >= 0.) && (zMinusAbsY >= 0.)) {
            t = 5;
        }
        else {
            t = 4;
        }
    }
    if (debug) {
        Log::info() << "tileFromLonLat:: lonlat abs xyz t = "
                    << std::setprecision(std::numeric_limits<double>::digits10 + 1) << zPlusAbsX << " " << zPlusAbsY
                    << " " << zMinusAbsX << " " << zMinusAbsY << " " << crd[LON] << " " << crd[LAT] << " " << xyz[XX]
                    << " " << xyz[YY] << " " << xyz[ZZ] << " " << t << std::endl;
    }
    return t;
}
idx_t FV3CubedSphereTiles :: indexFromLonLat ( const double crd [ ] ) const { idx_t t ( -1 ) ; double xyz [ 3 ] ; sphericalToCartesian ( crd , xyz ) ; static const double cornerLat = std :: asin ( 1. / std :: sqrt ( 3.0 ) ) * rad2deg ; const double & lon = crd [ LON ] ; const double & lat = crd [ LAT ] ; double zPlusAbsX = xyz [ ZZ ] + std :: abs ( xyz [ XX ] ) ; double zPlusAbsY = xyz [ ZZ ] + std :: abs ( xyz [ YY ] ) ; double zMinusAbsX = xyz [ ZZ ] - std :: abs ( xyz [ XX ] ) ; double zMinusAbsY = xyz [ ZZ ] - std :: abs ( xyz [ YY ] ) ; if ( is_tiny ( zPlusAbsX ) ) { zPlusAbsX = 0. ; } if ( is_tiny ( zPlusAbsY ) ) { zPlusAbsY = 0. ; } if ( is_tiny ( zMinusAbsX ) ) { zMinusAbsX = 0. ; } if ( is_tiny ( zMinusAbsY ) ) { zMinusAbsY = 0. ; } if ( is_geq ( lon , 315. ) || is_less ( lon , 45. ) ) { if ( ( zPlusAbsX <= 0. ) && ( zPlusAbsY <= 0. ) ) { t = 2 ; } else if ( ( zMinusAbsX > 0. ) && ( zMinusAbsY > 0. ) ) { t = 5 ; } else { t = 0 ; } if ( is_same ( lon , -45. ) && is_same ( lat , cornerLat ) ) { t = 0 ; } if ( is_same ( lon , 315. ) && is_same ( lat , cornerLat ) ) { t = 0 ; } } if ( is_geq ( lon , 45. ) && is_less ( lon , 135. ) ) { if ( ( zPlusAbsX <= 0. ) && ( zPlusAbsY <= 0. ) ) { t = 2 ; } else if ( ( zMinusAbsX > 0. ) && ( zMinusAbsY > 0. ) ) { t = 5 ; } else { t = 1 ; } } if ( is_geq ( lon , 135. ) && is_less ( lon , 225. ) ) { if ( ( zPlusAbsX < 0. ) && ( zPlusAbsY < 0. ) ) { t = 2 ; } else if ( ( zMinusAbsX >= 0. ) && ( zMinusAbsY >= 0. ) ) { t = 5 ; } else { t = 3 ; } if ( is_same ( lon , 135. ) && is_same ( lat , - cornerLat ) ) { t = 1 ; } } if ( is_geq ( lon , 225. ) && is_less ( lon , 315. ) ) { if ( ( zPlusAbsX < 0. ) && ( zPlusAbsY < 0. ) ) { t = 2 ; } else if ( ( zMinusAbsX >= 0. ) && ( zMinusAbsY >= 0. ) ) { t = 5 ; } else { t = 4 ; } } if ( debug ) { Log :: info ( ) << " " << std :: setprecision ( std :: numeric_limits < double > :: digits10 + 1 ) << zPlusAbsX << " " << zPlusAbsY << " " << zMinusAbsX << " " << zMinusAbsY << " " << crd [ LON ] << " " << crd [ LAT ] << " " << xyz [ XX ] << " " << xyz [ YY ] << " " << xyz [ ZZ ] << " " << t << std :: endl ; } return t ; }
Calculates the FV3 panel
Input (crd) is Longitude and Latitude in Radians
Output is tile number
Calculates the FV3 panel Input ( crd ) is Longitude and Latitude in Radians Output is tile number
// tile index // the magnitude of the latitude at the corners of the cube (not including the sign) // in radians. // Note that this method can lead to roundoff errors that can // cause the tile selection to fail. // To this end we enforce that tiny values close (in roundoff terms) // to a boundary should end up exactly on the boundary. // extra point corner point (considering two different longitudes depending on // whether we are on longitude points [0,360] or [-45.,45.] // interior // interior // extra point corner point // interior