initial burst metadata reader (#1)

* Create ISCE3-compatible Sentinel1 burst class given S1 SAFE, orbit file, subswath index, and polarization.
* Monotonically increasing bursts IDs.

Co-authored-by: Zhang Yunjun <yunjunzgeo@gmail.com>

.gitignore

@@ -0,0 +1,131 @@
+*.DS_Store
+
+# Byte-compiled / optimized / DLL files
+__pycache__/
+*.py[cod]
+*$py.class
+
+## C extensions
+*.so
+
+# Distribution / packaging
+.Python
+build/
+develop-eggs/
+dist/
+downloads/
+eggs/
+.eggs/
+lib/
+lib64/
+parts/
+sdist/
+var/
+wheels/
+pip-wheel-metadata/
+share/python-wheels/
+*.egg-info/
+.installed.cfg
+*.egg
+MANIFEST
+
+# PyInstaller
+#  Usually these files are written by a python script from a template
+#  before PyInstaller builds the exe, so as to inject date/other infos into it.
+*.manifest
+*.spec
+
+# Installer logs
+pip-log.txt
+pip-delete-this-directory.txt
+
+# Unit test / coverage reports
+htmlcov/
+.tox/
+.nox/
+.coverage
+.coverage.*
+.cache
+nosetests.xml
+coverage.xml
+*.cover
+*.py,cover
+.hypothesis/
+.pytest_cache/
+
+# Translations
+*.mo
+*.pot
+
+# Django stuff:
+*.log
+local_settings.py
+db.sqlite3
+db.sqlite3-journal
+
+# Flask stuff:
+instance/
+.webassets-cache
+
+# Scrapy stuff:
+.scrapy
+
+# Sphinx documentation
+docs/_build/
+
+# PyBuilder
+target/
+
+# Jupyter Notebook
+.ipynb_checkpoints
+
+# IPython
+profile_default/
+ipython_config.py
+
+# pyenv
+.python-version
+
+# pipenv
+#   According to pypa/pipenv#598, it is recommended to include Pipfile.lock in version control.
+#   However, in case of collaboration, if having platform-specific dependencies or dependencies
+#   having no cross-platform support, pipenv may install dependencies that don't work, or not
+#   install all needed dependencies.
+#Pipfile.lock
+
+# PEP 582; used by e.g. github.com/David-OConnor/pyflow
+__pypackages__/
+
+# Celery stuff
+celerybeat-schedule
+celerybeat.pid
+
+# SageMath parsed files
+*.sage.py
+
+# Environments
+.env
+.venv
+env/
+venv/
+ENV/
+env.bak/
+venv.bak/
+
+# Spyder project settings
+.spyderproject
+.spyproject
+
+# Rope project settings
+.ropeproject
+
+# mkdocs documentation
+/site
+
+# mypy
+.mypy_cache/
+.dmypy.json
+dmypy.json
+
+# Pyre type checker
+.pyre/

README.md

@@ -1,2 +1,45 @@
-# sentinel1-reader
-Sentinel-1 reader
+### Features
+
++ Create ISCE3-compatible Sentinel1 burst class given:
+
+  - S1 SAFE
+  - subswath index
+  - polarization
+  - path to orbit directory
+
++ Monotonically increasing bursts IDs.
+
+### Install
+
+1. Set up and activate virtual environment with ISCE3.
+2. Clone repository.
+
+```bash
+git clone https://github.com/LiangJYu/sentinel1-reader.git
+```
+
+3. Install into virtual environment with pip. From clone directory:
+
+```bash
+cd sentinel1-reader
+pip install .
+```
+
+### Usage
+
+The following sample code demonstrates how to process a single burst from a S1 SAFE zip:
+
+```python
+from sentinel1_reader import sentinel1_reader, sentinel1_orbit_reader
+
+zip_path = "S1A_IW_SLC__1SDV_20190909T134419_20190909T134446_028945_03483B_B9E1.zip"
+i_subswath = 2
+pol = "HH"
+
+# read orbits
+orbit_dir = '/home/user/data/sentinel1_orbits'
+orbit_path = sentinel1_orbit_reader.get_swath_orbit_file(zip_path, orbit_dir)
+
+# returns the list of the bursts
+bursts = sentinel1_reader.burst_from_zip(zip_path, orbit_path, i_subswath, pol)
+```

bin/s1_read.py

@@ -0,0 +1,29 @@
+import os
+import sys
+
+from sentinel1_reader import sentinel1_reader, sentinel1_orbit_reader
+
+if __name__ == "__main__":
+    # TODO replace with argparse
+    zip_path = sys.argv[1]
+    if not os.path.isfile(zip_path):
+        raise FileNotFoundError(f"{zip_path} does not exist")
+
+    i_subswath = int(sys.argv[2])
+    if i_subswath < 1  or i_subswath > 3:
+        raise ValueError("i_subswath not <1 or >3")
+
+    pol = sys.argv[3]
+    pols = ['vv', 'vh', 'hh', 'hv']
+    if pol not in pols:
+        raise ValueError("polarization not in {pols}")
+
+    orbit_dir = sys.argv[4]
+    if not os.path.isdir(orbit_dir):
+        raise NotADirectoryError(f"{orbit_dir} not found")
+    orbit_path = sentinel1_orbit_reader.get_swath_orbit_file(zip_path, orbit_dir)
+
+    bursts = sentinel1_reader.burst_from_zip(zip_path, orbit_path, i_subswath, pol)
+
+    for i, burst in enumerate(bursts):
+        print(burst.burst_id, burst.center)

pyproject.toml

@@ -0,0 +1,5 @@
+[build-system]
+requires = [
+    "setuptools>=42"
+]
+build-backend = "setuptools.build_meta"

setup.cfg

@@ -0,0 +1,22 @@
+[metadata]
+name = sentinel-reader
+version = 0.1.0
+description = A Sentinel1 metadata reader
+long_description = file: README.md
+long_description_content_type = text/markdown
+url = https://github.com/opera-adt/sentinel1-reader
+classifiers =
+    Development Status :: 2 - Pre-Alpha
+    Intended Audience :: Science/Research
+    Programming Language :: Python :: 3
+    License = file : LICENSE
+    Operating System :: OS Independent
+
+[options]
+package_dir =
+    = src
+packages = find:
+python_requires = >=3.7
+
+[options.packages.find]
+where = src

src/sentinel1_reader/sentinel1_burst_slc.py

@@ -0,0 +1,119 @@
+from dataclasses import dataclass
+import datetime
+
+from osgeo import gdal
+
+import isce3
+
+@dataclass(frozen=True)
+class Doppler:
+    poly1d: isce3.core.Poly1d
+    lut2d: isce3.core.LUT2d
+
+@dataclass(frozen=True)
+class Sentinel1BurstSlc:
+    '''Raw values extracted from SAFE XML.
+    '''
+    sensing_start: datetime.datetime# *
+    radar_center_frequency: float
+    wavelength: float
+    azimuth_steer_rate: float
+    azimuth_time_interval: float
+    slant_range_time: float
+    starting_range: float
+    range_sampling_rate: float
+    range_pixel_spacing: float
+    shape: tuple()
+    azimuth_fm_rate: isce3.core.Poly1d
+    doppler: Doppler
+    range_bandwidth: float
+    polarization: str # {VV, VH, HH, HV}
+    burst_id: str # t{track_number}_iw{1,2,3}_{burst_index}
+    platform_id: str # S1{A,B}
+    center: tuple # {center lon, center lat} in degrees
+    border: list # list of lon, lat coordinate tuples (in degrees) representing burst border
+    orbit: isce3.core.Orbit
+    # VRT params
+    tiff_path: str
+    i_burst: int
+    first_valid_sample: int
+    last_valid_sample: int
+    first_valid_line: int
+    last_valid_line: int
+
+    def as_isce3_radargrid(self):
+        '''Init and return isce3.product.RadarGridParameters.
+
+        Returns:
+        --------
+        _ : RadarGridParameters
+            RadarGridParameters constructed from class members.
+        '''
+
+        prf = 1 / self.azimuth_time_interval
+
+        length, width = self.shape
+
+        time_delta = datetime.timedelta(days=2)
+        ref_epoch = isce3.core.DateTime(self.sensing_start - time_delta)
+        # sensing start with respect to reference epoch
+        sensing_start = time_delta.total_seconds()
+
+        # init radar grid
+        return isce3.product.RadarGridParameters(sensing_start,
+                                                 self.wavelength,
+                                                 prf,
+                                                 self.starting_range,
+                                                 self.range_pixel_spacing,
+                                                 isce3.core.LookSide.Right,
+                                                 length,
+                                                 width,
+                                                 ref_epoch)
+
+    def to_vrt_file(self, out_path):
+        '''Write burst to VRT file.
+
+        Parameters:
+        -----------
+        out_path : string
+            Path of output VRT file.
+        '''
+        line_offset = self.i_burst * self.shape[0]
+
+        inwidth = self.last_valid_sample - self.first_valid_sample + 1
+        inlength = self.last_valid_line - self.first_valid_line + 1
+        outlength, outwidth = self.shape
+        yoffset = line_offset + self.first_valid_line
+        localyoffset = self.first_valid_line
+        xoffset = self.first_valid_sample
+        gdal_obj = gdal.Open(self.tiff_path, gdal.GA_ReadOnly)
+        fullwidth = gdal_obj.RasterXSize
+        fulllength = gdal_obj.RasterYSize
+
+        # TODO maybe cleaner to write with ElementTree
+        tmpl = f'''<VRTDataset rasterXSize="{outwidth}" rasterYSize="{outlength}">
+    <VRTRasterBand dataType="CFloat32" band="1">
+        <NoDataValue>0.0</NoDataValue>
+        <SimpleSource>
+            <SourceFilename relativeToVRT="1">{self.tiff_path}</SourceFilename>
+            <SourceBand>1</SourceBand>
+            <SourceProperties RasterXSize="{fullwidth}" RasterYSize="{fulllength}" DataType="CInt16"/>
+            <SrcRect xOff="{xoffset}" yOff="{yoffset}" xSize="{inwidth}" ySize="{inlength}"/>
+            <DstRect xOff="{xoffset}" yOff="{localyoffset}" xSize="{inwidth}" ySize="{inlength}"/>
+        </SimpleSource>
+    </VRTRasterBand>
+</VRTDataset>'''
+
+        with open(out_path, 'w') as fid:
+            fid.write(tmpl)
+
+    def get_sensing_mid(self):
+        '''Returns sensing mid as datetime object.
+
+        Returns:
+        --------
+        _ : datetime
+            Sensing mid as datetime object.
+        '''
+        d_seconds = 0.5 * (self.shape[0] - 1) * self.azimuth_time_interval
+        return self.sensing_start + datetime.timedelta(seconds=d_seconds)