Skip to content

Commit

Permalink
Merge pull request #69 from gshiroma/changes_after_calval_point_relea…
Browse files Browse the repository at this point in the history
…se_5

Changes for final delivery
  • Loading branch information
gshiroma authored Sep 5, 2023
2 parents e1bd926 + d057983 commit 81b51fe
Show file tree
Hide file tree
Showing 6 changed files with 63 additions and 32 deletions.
2 changes: 1 addition & 1 deletion Docker/Dockerfile
Original file line number Diff line number Diff line change
Expand Up @@ -2,7 +2,7 @@ FROM oraclelinux:8

LABEL author="OPERA ADT" \
description="RTC cal/val release R4" \
version="0.4.2-calval"
version="1.0.0-final"

RUN yum -y update &&\
yum -y install curl &&\
Expand Down
2 changes: 1 addition & 1 deletion Docker/environment.yml
Original file line number Diff line number Diff line change
@@ -1,4 +1,4 @@
name: rtc_sas_final
name: rtc_s1_sas_final
channels:
- conda-forge
dependencies:
Expand Down
2 changes: 1 addition & 1 deletion build_docker_image.sh
Original file line number Diff line number Diff line change
Expand Up @@ -2,7 +2,7 @@

REPO=opera
IMAGE=rtc
TAG=final_0.5.0
TAG=final_1.0.0

echo "IMAGE is $REPO/$IMAGE:$TAG"

Expand Down
85 changes: 58 additions & 27 deletions src/rtc/core.py
Original file line number Diff line number Diff line change
Expand Up @@ -248,9 +248,12 @@ def get_tile_srs_bbox(tile_min_y_projected, tile_max_y_projected,
tile_max_x = np.max(tile_x_array)

# handles antimeridian: tile_max_x around +180 and tile_min_x around -180
# add 360 to tile_min_x, so it becomes a little greater than +180
if tile_max_x > tile_min_x + 340:
tile_min_x, tile_max_x = tile_max_x, tile_min_x + 360
if polygon_srs.IsGeographic() and tile_max_x - tile_min_x > 180:
# unwrap negative longitude values
# move longitude range from [-180, 180] to [0, 360]
new_tile_x_array = [x + (x < 0) * 360 for x in tile_x_array]
tile_min_x = np.min(new_tile_x_array)
tile_max_x = np.max(new_tile_x_array)

tile_ring = ogr.Geometry(ogr.wkbLinearRing)
tile_ring.AddPoint(tile_min_x, tile_max_y)
Expand All @@ -265,7 +268,7 @@ def get_tile_srs_bbox(tile_min_y_projected, tile_max_y_projected,


def _antimeridian_crossing_requires_special_handling(
file_srs, file_min_x, tile_min_x, tile_max_x):
file_srs, ancillary_min_x, ancillary_max_x, tile_min_x, tile_max_x):
'''
Check if ancillary input requires special handling due to
the antimeridian crossing
Expand All @@ -274,8 +277,10 @@ def _antimeridian_crossing_requires_special_handling(
----------
file_srs: osr.SpatialReference
Ancillary file spatial reference system (SRS)
file_min_x: float
ancillary_min_x: float
Ancillary file min longitude value in degrees
ancillary_max_x: float
Ancillary file max longitude value in degrees
tile_min_x: float
Tile min longitude value in degrees
tile_max_x: float
Expand All @@ -287,24 +292,36 @@ def _antimeridian_crossing_requires_special_handling(
Flag that indicate if the ancillary input requires special handling
'''

# Flag to indicate if the if the tile crosses the antimeridian.
flag_tile_crosses_antimeridian = tile_min_x < 180 and tile_max_x >= 180

# Flag to test if the ancillary input file is in geographic
# coordinates and if its longitude domain is represented
# within the [-180, +180] range, rather than, for example, inside
# the [0, +360] interval.
# This is verified by the test `min_x < -165`. There's no specific reason
# why -165 is used. It could be -160, or even 0. However, testing for
# -165 is more general than -160 or 0, but still not too close to -180.
flag_input_geographic_and_longitude_lt_m165 = \
file_srs.IsGeographic() and file_min_x < -165

# If both are true, tile requires special handling due to the
# antimeridian crossing
flag_requires_special_handling = (
flag_tile_crosses_antimeridian and
flag_input_geographic_and_longitude_lt_m165)
# The coordinate wrapping around the ancillary file
# discontinuities such as the antimeridian only happens
# if the DEM is provided in geographic coordinates
if not file_srs.IsGeographic():
flag_requires_special_handling = False
return flag_requires_special_handling

# Check whether the ancillary file covers the entire longitude range
# Mark flag as True if the longitude coordinates span completes
# or is close to complete the circle (360 degrees).
# We use 359 instead of 360 degrees to have some buffer
# for incomplete ancillary files
flag_ancillary_covers_entire_longitude_range = \
(ancillary_max_x - ancillary_min_x) > 359

if not flag_ancillary_covers_entire_longitude_range:
flag_requires_special_handling = False
return flag_requires_special_handling

# Flag to indicate if the tile crosses the ancillary file
# discontinuity (ancillary_max_x).
# Notice that`ancillary_max_x` refers to the eastern edge
# of the ancillary file, which, for the antimeridian case,
# may not lay exactly on +/- 180. For example, it may be
# 179.9998611111111
flag_tile_crosses_discontinuity = \
tile_min_x < ancillary_max_x < tile_max_x

flag_requires_special_handling = \
flag_tile_crosses_discontinuity

return flag_requires_special_handling

Expand Down Expand Up @@ -430,7 +447,8 @@ def check_ancillary_inputs(check_ancillary_inputs_coverage,

# If needed, test for antimeridian ("dateline") crossing
if _antimeridian_crossing_requires_special_handling(
ancillary_srs, ancillary_x0, geogrid_x0, geogrid_xf):
ancillary_srs, ancillary_x0, ancillary_xf,
geogrid_x0, geogrid_xf):

logger.info(f'The input RTC-S1 product crosses the antimeridian'
' (dateline). Verifying the'
Expand All @@ -446,15 +464,28 @@ def check_ancillary_inputs(check_ancillary_inputs_coverage,
logger.info(f' left side (-180 -> +180): {check_1_str}')

# Right side of the antimeridian crossing: +180 -> +360
ancillary_polygon_2 = _get_ogr_polygon(
#
# Get the intersection between the geogrid and the right side
# of the antimeridian (with a litter buffer represented by
# ANTIMERIDIAN_CROSSING_RIGHT_SIDE_TEST_BUFFER)
antimeridian_right_side_polygon = _get_ogr_polygon(
ancillary_xf + ANTIMERIDIAN_CROSSING_RIGHT_SIDE_TEST_BUFFER,
90, ancillary_xf + 360, -90, ancillary_srs)
intersection_2 = geogrid_polygon.Intersection(ancillary_polygon_2)

intersection_2 = geogrid_polygon.Intersection(
antimeridian_right_side_polygon)

# Create a polygon representing the ancillary dataset
# at the right side of the antimeridian
ancillary_polygon_2 = _get_ogr_polygon(
ancillary_x0 + 360, ancillary_y0,
ancillary_xf + 360, ancillary_yf,
ancillary_srs)
flag_2_ok = intersection_2.Within(ancillary_polygon_2)

# Check if the geogrid-intersected area (if valid) is within
# the ancillary polygon
flag_2_ok = (intersection_2.IsEmpty() or
intersection_2.Within(ancillary_polygon_2))
check_2_str = 'ok' if flag_2_ok else 'fail'
logger.info(f' right side (+180 -> +360): {check_2_str}')

Expand Down
2 changes: 1 addition & 1 deletion src/rtc/h5_prep.py
Original file line number Diff line number Diff line change
Expand Up @@ -23,7 +23,7 @@
DATE_TIME_FILENAME_FORMAT = '%Y%m%dT%H%M%SZ'

DATA_BASE_GROUP = '/data'
PRODUCT_SPECIFICATION_VERSION = 0.1
PRODUCT_SPECIFICATION_VERSION = 1.0

logger = logging.getLogger('rtc_s1')

Expand Down
2 changes: 1 addition & 1 deletion src/rtc/version.py
Original file line number Diff line number Diff line change
@@ -1 +1 @@
VERSION='0.4.2'
VERSION='1.0.0'

0 comments on commit 81b51fe

Please sign in to comment.