Skip to content

Commit

Permalink
Incorporating geojson2txt utilities from TectonicUtils
Browse files Browse the repository at this point in the history
  • Loading branch information
kmaterna committed Aug 7, 2023
1 parent 4b8caaa commit 929bcda
Show file tree
Hide file tree
Showing 5 changed files with 100 additions and 4 deletions.
1 change: 0 additions & 1 deletion src/Downsample/__init__.py
Original file line number Diff line number Diff line change
@@ -1 +0,0 @@
from . import quadtree_downsample_kite
72 changes: 72 additions & 0 deletions src/Downsample/geojson2txt.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,72 @@
"""
This script contains utility functions to convert between two different InSAR data representations:
1. geojson, produced by Kite after downsampling
2. text representation, used for inversions
"""

import numpy as np
import json
from .pixel_object import Downsampled_pixel


def read_geojson(infile):
"""
Read a geojson as created by Kite downsampling into downsampled pixel objects.
:param infile: name of geojson file
:type infile: string
:return: list of Downsampled_pixel objects
:rtype: list
"""
# "features" is a list with n elements
# Each element is a downsampled pixel and is stored as a dictionary
# containing 'type','id','geometry','properties'
# These properties will be unpacked into a named tuple that contains helpful information.

with open(infile) as f:
data = json.load(f);
features = data["features"]; # beginning to unpack the geojson
pixel_list = [];
for feature in features:
# Unpacking what's in each pixel
bl = feature["geometry"]["coordinates"][0][0];
tr = feature["geometry"]["coordinates"][0][2];
mean = feature["properties"]["mean"] * 0.001
median = feature["properties"]["median"] * 0.001
std = feature["properties"]["std"] * 0.001
unitE = feature["properties"]["unitE"]
unitN = feature["properties"]["unitN"]
unitU = feature["properties"]["unitU"]
onePixel = Downsampled_pixel(mean=mean, median=median, std=std, BL_corner=bl, TR_corner=tr,
unitE=unitE, unitN=unitN, unitU=unitU);
pixel_list.append(onePixel);
return pixel_list;


def pixels_to_txt(pixel_list, text_file, bbox=(-180, 180, -90, 90), std_min=0.001):
"""
Write InSAR pixels in basic text format: Lon Lat Value unitE unitN unitU.
:param pixel_list: list of Downsampled_pixel objects
:type pixel_list: list
:param text_file: name of output file
:type text_file: str
:param bbox: tuple of (W, E, S, N) bounding box, defaults to (-180, 180, -90, 90)
:type bbox: array_like, optional
:param std_min: minimum value of uncertainty (m), defaults to 0.001
:type std_min: float, optional
"""
ofile = open(text_file, 'w');
ofile.write("# Header: lon, lat, disp(m), sig(m), unitE, unitN, unitU from ground to satellite\n")
for pixel in pixel_list:
lon = np.mean([pixel.BL_corner[0], pixel.TR_corner[0]]);
lat = np.mean([pixel.BL_corner[1], pixel.TR_corner[1]]);
std = np.max([pixel.std, std_min]); # don't let uncertainty get unreasonably small
if bbox[0] <= lon <= bbox[1]:
if bbox[2] <= lat <= bbox[3]:
ofile.write("%.5f %.5f %.5f %.5f %.5f %.5f %.5f\n" % (
lon, lat, pixel.median, std, pixel.unitE, pixel.unitN, pixel.unitU));
ofile.close();
print("Writing %s with %d pixels " % (text_file, len(pixel_list)));
return;
25 changes: 25 additions & 0 deletions src/Downsample/pixel_object.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,25 @@
# the downsampled pixel object

class Downsampled_pixel:
"""
An object containing a quadtree-downsampled pixel, including
downsampled pixel footprint, downsampled pixel look vector, and downsampled deformation values.
"""
def __init__(self, mean, median, std, BL_corner, TR_corner, unitE, unitN, unitU):
self.mean = mean; # mean of LOS values within the pixel (meters)
self.median = median; # median of LOS values within the pixel (meters)
self.std = std; # standard deviation of LOS values within the pixel (meters)
self.BL_corner = BL_corner; # Coordinates of Bottom Left corner (longitude, latitude)
self.TR_corner = TR_corner; # Coordinates of Top Right corner (longitude, latitude)
self.unitE = unitE; # east component of unit vector from ground to satellite
self.unitN = unitN; # north component of unit vector from ground to satellite
self.unitU = unitU; # up component of unit vector from ground to satellite

def set_mean(self, new_mean):
self.mean = new_mean;

def set_median(self, new_median):
self.median = new_median;

def set_std(self, new_std):
self.std = new_std;
3 changes: 1 addition & 2 deletions src/Downsample/quadtree_downsample_kite.py
Original file line number Diff line number Diff line change
Expand Up @@ -12,8 +12,7 @@
Not used: fig = qt.plot() # this will open a plot.show(), won't save it.
"""

from Tectonic_Utils.geodesy import geojson2txt
from . import plot_geojson
from . import plot_geojson, geojson2txt


def kite_downsample_isce_unw(datafile, outname,
Expand Down
3 changes: 2 additions & 1 deletion src/bin/Projects/Mendo_inv/mendo_driver.py
Original file line number Diff line number Diff line change
Expand Up @@ -103,7 +103,8 @@ def run_main():
fault_traces_from_file=outdir+"/temp-outfile.txt",
scale_arrow=(1.0, 0.001, "1 mm"),
region=(-124.7, -122, 39.05, 41.8), vert_mult=1000,
vert_disp_units='mm', vmin=-2, vmax=6);
vert_disp_units='mm', vmin=-2, vmax=6,
slip_cbar_opts=(-0.1, 0.1, 0.001));
fso.plot_fault_slip.map_source_slip_distribution([], outdir+'/obs_disps.png', disp_points=obs_data_points,
scale_arrow=(1.0, 0.001, "1 mm"),
region=(-124.7, -122, 39.05, 41.8), vert_mult=1000,
Expand Down

0 comments on commit 929bcda

Please sign in to comment.