Skip to content

Commit

Permalink
#2218 Check and adjust semi_major_axis in meters
Browse files Browse the repository at this point in the history
  • Loading branch information
Howard Soh committed Jan 31, 2023
1 parent 256ee5c commit 34fe360
Showing 1 changed file with 22 additions and 22 deletions.
44 changes: 22 additions & 22 deletions src/libcode/vx_data2d_nccf/nccf_file.cc
View file
Original file line number Diff line number Diff line change
Expand Up @@ -2198,24 +2198,14 @@ void NcCfFile::get_grid_mapping_polar_stereographic(const NcVar *grid_mapping_va

if (is_eq(semi_major_axis, bad_data_double))
semi_major_axis = EARTH_MAJOR_AXIS_km * m_per_km;
else semi_major_axis *= x_coord_to_m_cf; // meters
else if (semi_major_axis < 10000.0) semi_major_axis *= m_per_km; // meters

// Calculate the pin indices. The pin will be located at the grid's reference
// location since that's the only lat/lon location we know about.

double x_pin = -(x_values[0] / dx_m);
double y_pin = -(y_values[0] / dy_m);

if(!has_scale_factor && has_standard_parallel) {
if(is_eq(eccentricity, 0.0)) {
x_pin = y_pin = 0.;
}
else {
x_pin = x_values[0] / dx_m_a;
y_pin = y_values[0] / dy_m_a;
}
}

// Fill in the data structure. Remember to negate the longitude values.

StereographicData data;
Expand All @@ -2227,9 +2217,9 @@ void NcCfFile::get_grid_mapping_polar_stereographic(const NcVar *grid_mapping_va
data.y_pin = y_pin;
data.scale_lat = proj_origin_lat;
data.lon_orient = -1.0 * proj_vertical_lon;
data.d_km = dx_m / 1000.0;
data.dy_km = dy_m / 1000.0;
data.r_km = semi_major_axis / 1000.0;
data.d_km = dx_m / m_per_km;
data.dy_km = dy_m / m_per_km;
data.r_km = semi_major_axis / m_per_km;
data.nx = _xDim->getSize();
data.ny = _yDim->getSize();

Expand Down Expand Up @@ -2258,14 +2248,24 @@ void NcCfFile::get_grid_mapping_polar_stereographic(const NcVar *grid_mapping_va
data.lat_pin = lat;
data.lon_pin = -lon;

if(!has_scale_factor && has_standard_parallel) {
if(is_eq(eccentricity, 0.0)) {
data.x_pin = data.y_pin = 0.;
}
else {
data.x_pin = x / dx_m_a;
data.y_pin = y / dy_m_a;
}
}

st_latlon_to_xy_func(lat, lon, x2, y2, scale_factor, proj_vertical_lon,
semi_major_axis, false_east, false_north,
eccentricity, is_north_hemisphere);
mlog << Debug(6) << method_name
<< " x: between " << x_values[0] <<" and " << x_values[x_counts-1]
<< ", y: between " << y_values[0] <<" and " << y_values[y_counts-1] << "\n";
mlog << Debug(6) << method_name
<< "pin: (x,y)m -> (lat,lon) -> (x2,y2)m: (" << x << "," << y << ") -> ("
<< "pin: (x,y) -> (lat,lon) -> (x2,y2): (" << x << "," << y << ") -> ("
<< lat << "," << lon << ") -> (" << x2 << "," << y2
<< ") Diff (m): x=" << (x-x2) <<", y=" << (y-y2) << "\n";

Expand Down Expand Up @@ -2297,7 +2297,7 @@ void NcCfFile::get_grid_mapping_polar_stereographic(const NcVar *grid_mapping_va
semi_major_axis, false_east, false_north,
eccentricity, is_north_hemisphere);
mlog << Debug(10) << method_name
<< "bottom left: (x,y)m -> (lat,lon) -> (x2,y2)m: ("
<< " bottom left: (x,y) -> (lat,lon) -> (x2,y2): ("
<< x1 << "," << y1 << ") -> ("
<< lat1 << "," << lon1 << ") -> (" << x2 << "," << y2
<< ") Diff (m): x=" << (x1-x2) <<", y=" << (y1-y2)
Expand All @@ -2312,7 +2312,7 @@ void NcCfFile::get_grid_mapping_polar_stereographic(const NcVar *grid_mapping_va
semi_major_axis, false_east, false_north,
eccentricity, is_north_hemisphere);
mlog << Debug(10) << method_name
<< " top left: (x,y)m -> (lat,lon) -> (x2,y2)m: ("
<< " top left: (x,y) -> (lat,lon) -> (x2,y2): ("
<< x1 << "," << y1 << ") -> ("
<< lat1 << "," << lon1 << ") -> (" << x2 << "," << y2
<< ") Diff (m): x=" << (x1-x2) <<", y=" << (y1-y2)
Expand All @@ -2327,7 +2327,7 @@ void NcCfFile::get_grid_mapping_polar_stereographic(const NcVar *grid_mapping_va
semi_major_axis, false_east, false_north,
eccentricity, is_north_hemisphere);
mlog << Debug(10) << method_name
<< " center: (x,y)m -> (lat,lon) -> (x2,y2)m: (" << x1
<< " center: (x,y) -> (lat,lon) -> (x2,y2): (" << x1
<< "," << y1 << ") -> ("
<< lat1 << "," << lon1 << ") -> (" << x2 << "," << y2
<< ") Diff (m): x=" << (x1-x2) <<", y=" << (y1-y2)
Expand All @@ -2342,7 +2342,7 @@ void NcCfFile::get_grid_mapping_polar_stereographic(const NcVar *grid_mapping_va
semi_major_axis, false_east, false_north,
eccentricity, is_north_hemisphere);
mlog << Debug(10) << method_name
<< " top right: (x,y)m -> (lat,lon) -> (x2,y2)m: (" << x1
<< " top right: (x,y) -> (lat,lon) -> (x2,y2): (" << x1
<< "," << y1 << ") -> ("
<< lat1 << "," << lon1 << ") -> (" << x2 << "," << y2
<< ") Diff (m): x=" << (x1-x2) <<", y=" << (y1-y2)
Expand All @@ -2357,7 +2357,7 @@ void NcCfFile::get_grid_mapping_polar_stereographic(const NcVar *grid_mapping_va
semi_major_axis, false_east, false_north,
eccentricity, is_north_hemisphere);
mlog << Debug(10) << method_name
<< "bottom right: (x,y)m -> (lat,lon) -> (x2,y2)m: (" << x1
<< "bottom right: (x,y) -> (lat,lon) -> (x2,y2): (" << x1
<< "," << y1 << ") -> ("
<< lat1 << "," << lon1 << ") -> (" << x2 << "," << y2
<< ") Diff (m): x=" << (x1-x2) <<", y=" << (y1-y2)
Expand All @@ -2384,8 +2384,8 @@ void NcCfFile::get_grid_mapping_polar_stereographic(const NcVar *grid_mapping_va
data.eccentricity, is_north_hemisphere);
x_diff = (x1-x2);
y_diff = (y1-y2);
mlog << Debug(15) << method_name
<< "index=("<< ix << "," << iy << "): (x,y)m -> (lat,lon) -> (x2,y2)m: (" << x1
mlog << Debug(15)
<< " index=("<< ix << "," << iy << "): (x,y) -> (lat,lon) -> (x2,y2): (" << x1
<< "," << y1 << ") -> ("
<< lat1 << "," << lon1 << ") -> (" << x2 << "," << y2
<< ") Diff (m): x=" << x_diff <<", y=" << y_diff
Expand Down

0 comments on commit 34fe360

Please sign in to comment.