Skip to content

Commit

Permalink
subset.py: support of --lat --lon in radar coordinates
Browse files Browse the repository at this point in the history
  • Loading branch information
Falk Amelung committed Jul 25, 2024
1 parent d2e551a commit febf7c5
Show file tree
Hide file tree
Showing 3 changed files with 51 additions and 19 deletions.
2 changes: 2 additions & 0 deletions src/mintpy/cli/subset.py
Original file line number Diff line number Diff line change
Expand Up @@ -20,6 +20,7 @@
EXAMPLE = """example:
subset.py inputs/ifgramStack.h5 -y 400 1500 -x 200 600
subset.py geo_velocity.h5 -l 30.5 30.8 -L 130.3 130.9
subset.py velocity.h5 --lat 25.894 25.896 --lon -80.124 -80.122
subset.py 030405_090801.unw -t SinabungT495F50AlosA.template
subset.py demLat*.dem.wgs84 --lat 32.5 33.0 --lon 130.2 130.6 -o srtm1.h5
Expand All @@ -28,6 +29,7 @@
# multiple files input
subset.py *velocity*.h5 timeseries*.h5 -y 400 1500 -x 200 600
subset.py geometryRadar.h5 slcStack.h5 --lat 25.894 25.896 --lon -80.124 -80.122
# crop to larger area with custom fill value
subset.py geo_velocity.h5 -l 32.2 33.5 --outfill-nan
Expand Down
8 changes: 6 additions & 2 deletions src/mintpy/objects/coord.py
Original file line number Diff line number Diff line change
Expand Up @@ -490,16 +490,20 @@ def bbox_radar2geo(self, pix_box, print_msg=False):
return geo_box


def bbox_geo2radar(self, geo_box, print_msg=False):
def bbox_geo2radar(self, geo_box, buf=None, print_msg=False):
"""Calculate bounding box in x/y for file in radar coord, based on input geo box.
Parameters: geo_box - tuple of 4 float, indicating the UL/LR lon/lat
buf - buffer [number of pixel]
Returns: pix_box - tuple of 4 int, indicating the UL/LR x/y of the bounding box in radar coord
for the corresponding lat/lon coverage.
"""
lat = np.array([geo_box[3], geo_box[3], geo_box[1], geo_box[1]])
lon = np.array([geo_box[0], geo_box[2], geo_box[0], geo_box[2]])
y, x, y_res, x_res = self.geo2radar(lat, lon, print_msg=print_msg)
buf = 2 * np.max(np.abs([x_res, y_res]))

if buf is None:
buf = 2 * np.max(np.abs([x_res, y_res]))

pix_box = (np.min(x) - buf, np.min(y) - buf,
np.max(x) + buf, np.max(y) + buf)
return pix_box
Expand Down
60 changes: 43 additions & 17 deletions src/mintpy/subset.py
Original file line number Diff line number Diff line change
Expand Up @@ -183,27 +183,53 @@ def subset_input_dict2box(subset_dict, meta_dict):

# Use subset_lat/lon input if existed, priority: lat/lon > y/x > len/wid
coord = ut.coordinate(meta_dict)
if subset_dict.get('subset_lat', None):
sub_y = coord.lalo2yx(subset_dict['subset_lat'], None)[0]
elif subset_dict['subset_y']:
sub_y = subset_dict['subset_y']
else:
sub_y = [0, length]

if subset_dict.get('subset_lon', None):
sub_x = coord.lalo2yx(None, subset_dict['subset_lon'])[1]
elif subset_dict['subset_x']:
sub_x = subset_dict['subset_x']
if 'Y_FIRST' in meta_dict.keys():
geocoded_flag =True
else:
sub_x = [0, width]
geocoded_flag =False

# Get subset box in y/x
sub_x = sorted(sub_x)
sub_y = sorted(sub_y)
pix_box = (sub_x[0], sub_y[0], sub_x[1], sub_y[1])
if geocoded_flag:
if subset_dict.get('subset_lat', None):
sub_y = coord.lalo2yx(subset_dict['subset_lat'], None)[0]
else:
sub_y = [0, length]

# Get subset box in lat/lon from subset box in y/x
geo_box = coord.box_pixel2geo(pix_box)
if subset_dict.get('subset_lon', None):
sub_x = coord.lalo2yx(None, subset_dict['subset_lon'])[1]
else:
sub_x = [0, width]

# Get subset box in y/x
sub_x = sorted(sub_x)
sub_y = sorted(sub_y)
pix_box = (sub_x[0], sub_y[0], sub_x[1], sub_y[1])

# Get subset box in lat/lon from subset box in y/x
geo_box = coord.box_pixel2geo(pix_box)

# for data in radar coordinates
else:
if subset_dict.get('subset_lat', None):
geo_box = (subset_dict['subset_lon'][0], subset_dict['subset_lat'][1], subset_dict['subset_lon'][1], subset_dict['subset_lat'][0])
pix_box = coord.bbox_geo2radar(geo_box, buf=0)
else:
if subset_dict['subset_y']:
sub_y = subset_dict['subset_y']
else:
sub_y = [0, length]

if subset_dict['subset_x']:
sub_x = subset_dict['subset_x']
else:
sub_x = [0, width]
# Get subset box in y/x
sub_x = sorted(sub_x)
sub_y = sorted(sub_y)
pix_box = (sub_x[0], sub_y[0], sub_x[1], sub_y[1])

# Get subset box in lat/lon from subset box in y/x
geo_box = coord.box_pixel2geo(pix_box)

return pix_box, geo_box

Expand Down

0 comments on commit febf7c5

Please sign in to comment.