diff --git a/.travis.yml b/.travis.yml new file mode 100644 index 0000000000000000000000000000000000000000..3f9a603b7583e263655f616a43dccdd18fed55e8 --- /dev/null +++ b/.travis.yml @@ -0,0 +1,36 @@ +language: python +sudo: false +os: + - linux + - osx + + env: + - PYTHON_VERSION=3.5 + - PYTHON_VERSION=3.6 + - PYTHON_VERSION=3.7 + + before_install: + # We do this conditionally because it saves us some downloading if the + # version is the same. + - if [ "$TRAVIS_OS_NAME" == "linux" ]; then + wget https://repo.continuum.io/miniconda/Miniconda3-latest-Linux-x86_64.sh -O miniconda.sh; + else + curl -o miniconda.sh https://repo.continuum.io/miniconda/Miniconda3-latest-MacOSX-x86_64.sh; + fi + - bash miniconda.sh -b -p $HOME/miniconda + - export PATH="$HOME/miniconda/bin:$PATH" + - hash -r + - conda config --set always_yes yes --set changeps1 no + - conda update -q conda + # Useful for debugging any issues with conda + - conda info -a + # Create the env + - conda create -q -n test python=$PYTHON_VERSION + +script: + - pytest knoten + +after_success: + - coveralls + - conda install conda-build anaconda-client + - conda build --token $CONDA_UPLOAD_TOKEN --python $PYTHON_VERSION -c conda-forge -c menpo -c usgs-astrogeology conda diff --git a/environment.yml b/environment.yml new file mode 100644 index 0000000000000000000000000000000000000000..dd6ed96a4f8b294d492c4fa78134d40fafc970c8 --- /dev/null +++ b/environment.yml @@ -0,0 +1,13 @@ +name: knoten +channels: + - conda-forge + - usgs-astrogeology +dependencies: + - python >= 3 + - coveralls + - csmapi + - gdal + - jinja + - numpy + - requests + - scipy \ No newline at end of file diff --git a/spaetzle/__init__.py b/knoten/__init__.py similarity index 100% rename from spaetzle/__init__.py rename to knoten/__init__.py diff --git a/spaetzle/csm.py b/knoten/csm.py similarity index 66% rename from spaetzle/csm.py rename to knoten/csm.py index 706944b0127c1e87465738217fd1bf7e3e1525c5..6a23480bddbec19c749e5baa35e62dea38da2847 100644 --- a/spaetzle/csm.py +++ b/knoten/csm.py @@ -4,80 +4,102 @@ import os from csmapi import csmapi import jinja2 - +from gdal import ogr import numpy as np -import scipy.stats import pyproj -import gdal -from gdal import ogr -import pvl +import requests +import scipy.stats +class NumpyEncoder(json.JSONEncoder): + def default(self, obj): + if isinstance(obj, np.ndarray): + return obj.tolist() + elif isinstance(obj, datetime.date): + return obj.isoformat() + return json.JSONEncoder.default(self, obj) -def generate_boundary(isize, nnodes=5, n_points=10): - ''' - Generates a bounding box given a camera model +def create_camera(label, url='http://pfeffer.wr.usgs.gov/api/1.0/pds/'): + """ + Given an ALE supported label, create a CSM compliant ISD file. This func + defaults to supporting PDS labels. The URL kwarg can be used to point + to the appropriate pfeffernusse endpoint for a given input data type. Parameters ---------- - camera : object - csmapi generated camera model - - nnodes : int - Not sure + label : str + The image label + + url : str + The service endpoint what is being acessed to generate an ISD. - n_points : int + Returns + ------- + model : object + A CSM sensor model (or None if no associated model is available.) + """ + data = json.dumps(label, cls=NumpyEncoder) + response = requests.post(url, data=data) + fname = '' + with open(fname, 'w') as f: + json.dump(response.json(), f) + isd = csmapi.Isd(fname) + + plugin = csmapi.Plugin.findPlugin('UsgsAstroPluginCSM') + + model_name = response.json()['name_model'] + if plugin.canModelBeConstructedFromISD(isd, model_name): + model = plugin.constructModelFromISD(isd, model_name) + return model + +def generate_boundary(isize, npoints=10): + ''' + Generates a bounding box given a camera model + Parameters + ---------- + isize : list + image size in the form xsize, ysize + npoints : int Number of points to generate between the corners of the bounding box per side. - Returns ------- boundary : list List of full bounding box ''' - x = np.linspace(0, isize[0], n_points) - y = np.linspace(0, isize[1], n_points) + x = np.linspace(0, isize[0], npoints) + y = np.linspace(0, isize[1], npoints) boundary = [(i,0.) for i in x] + [(isize[0], i) for i in y[1:]] +\ [(i, isize[1]) for i in x[::-1][1:]] + [(0.,i) for i in y[::-1][1:]] return boundary -def generate_latlon_boundary(camera, nnodes=5, semi_major=3396190, semi_minor=3376200, n_points=10): +def generate_latlon_boundary(camera, boundary, semi_major=3396190, semi_minor=3376200): ''' Generates a latlon bounding box given a camera model - Parameters ---------- camera : object csmapi generated camera model - + boundary : list + of boundary coordinates nnodes : int Not sure - semi_major : int Semimajor axis of the target body - semi_minor : int Semiminor axis of the target body - n_points : int Number of points to generate between the corners of the bounding box per side. - Returns ------- lons : list List of longitude values - lats : list List of latitude values - alts : list List of altitude values ''' - isize = camera.getImageSize() - isize = (isize.line, isize.samp) - - boundary = generate_boundary(isize, nnodes=nnodes, n_points = n_points) ecef = pyproj.Proj(proj='geocent', a=semi_major, b=semi_minor) lla = pyproj.Proj(proj='latlon', a=semi_major, b=semi_minor) @@ -96,38 +118,33 @@ def generate_latlon_boundary(camera, nnodes=5, semi_major=3396190, semi_minor=33 lons, lats, alts = pyproj.transform(ecef, lla, gnds[:,0], gnds[:,1], gnds[:,2]) return lons, lats, alts -def generate_gcps(camera, nnodes=5, semi_major=3396190, semi_minor=3376200, n_points=10): +def generate_gcps(camera, boundary, semi_major=3396190, semi_minor=3376200): ''' Generates an area of ground control points formated as: <GCP Id="" Info="" Pixel="" Line="" X="" Y="" Z="" /> per record - Parameters ---------- camera : object csmapi generated camera model - + boundary : list + of image boundary coordinates nnodes : int Not sure - semi_major : int Semimajor axis of the target body - semi_minor : int Semiminor axis of the target body - n_points : int Number of points to generate between the corners of the bounding box per side. - Returns ------- gcps : list List of all gcp records generated ''' - lons, lats, alts = generate_latlon_boundary(camera, nnodes=nnodes, + lons, lats, alts = generate_latlon_boundary(camera, boundary, semi_major=semi_major, - semi_minor=semi_minor, - n_points=n_points) + semi_minor=semi_minor) lla = np.vstack((lons, lats, alts)).T @@ -140,40 +157,33 @@ def generate_gcps(camera, nnodes=5, semi_major=3396190, semi_minor=3376200, n_po return gcps -def generate_latlon_footprint(camera, nnodes=5, semi_major=3396190, semi_minor=3376200, n_points=10): +def generate_latlon_footprint(camera, boundary, semi_major=3396190, semi_minor=3376200): ''' Generates a latlon footprint from a csmapi generated camera model - Parameters ---------- camera : object csmapi generated camera model - nnodes : int Not sure - semi_major : int Semimajor axis of the target body - semi_minor : int Semiminor axis of the target body - n_points : int Number of points to generate between the corners of the bounding box per side. - Returns ------- : object ogr multipolygon containing between one and two polygons ''' - lons, lats, _ = generate_latlon_boundary(camera, nnodes=nnodes, - semi_major=semi_major, - semi_minor=semi_minor, - n_points=n_points) + lons, lats, _ = generate_latlon_boundary(camera, boundary, + semi_major=semi_major, + semi_minor=semi_minor) # Transform coords from -180, 180 to 0, 360 - # Makes crossing the maridian easier to identify + # Makes crossing the meridian easier to identify ll_coords = [*zip(((lons + 180) % 360), lats)] ring = ogr.Geometry(ogr.wkbLinearRing) @@ -185,7 +195,7 @@ def generate_latlon_footprint(camera, nnodes=5, semi_major=3396190, semi_minor=3 multipoly = ogr.Geometry(ogr.wkbMultiPolygon) current_ring = ring - switch_point = None + switch_point = (None, None) previous_point = ll_coords[0] for coord in ll_coords: @@ -225,28 +235,24 @@ def generate_latlon_footprint(camera, nnodes=5, semi_major=3396190, semi_minor=3 return multipoly -def generate_bodyfixed_footprint(camera, nnodes=5, semi_major=3396190, semi_minor=3376200, n_points=10): +def generate_bodyfixed_footprint(camera, boundary, semi_major=3396190, semi_minor=3376200): ''' Generates a bodyfixed footprint from a csmapi generated camera model - Parameters ---------- camera : object csmapi generated camera model - nnodes : int Not sure - n_points : int Number of points to generate between the corners of the bounding box per side. - Returns ------- : object ogr polygon ''' - latlon_fp = generate_latlon_footprint(camera, nnodes=5, semi_major = semi_major, semi_minor = semi_minor, n_points = n_points) + latlon_fp = generate_latlon_footprint(camera, boundary, semi_major = semi_major, semi_minor = semi_minor) ecef = pyproj.Proj(proj='geocent', a=semi_major, b=semi_minor) lla = pyproj.Proj(proj='latlon', a=semi_major, b=semi_minor) @@ -265,15 +271,45 @@ def generate_bodyfixed_footprint(camera, nnodes=5, semi_major=3396190, semi_mino return latlon_fp -def generate_vrt(camera, raster_size, fpath, outpath=None, no_data_value=0): - gcps = generate_gcps(camera) - xsize, ysize = raster_size +def generate_vrt(raster_size, gcps, fpath, + no_data_value=0, + proj='+proj=longlat +a=3396190 +b=3376200 +no_defs'): + """ + Create a GDAl VRT string from a list of ground control points and + a projection. + + Parameters + ---------- + raster_size : iterable + in the form xsize, ysize + + gcps : list + of GCPs (likely created by the generate_gcp function) + + fpath : str + path to the source file that the VRT points to - if outpath is None: - outpath = os.path.dirname(fpath) - outname = os.path.splitext(os.path.basename(fpath))[0] + '.vrt' - outname = os.path.join(outpath, outname) + no_data_value : numeric + the no data value for the VRT (default=0) + + proj : str + A proj4 projection string for the VRT. + Returns + ------- + template : str + The rendered VRT string + + Example + -------- + >>> camera = create_camera(label) # Get a camera + >>> boundary = generate_boundary((100,100), 10) # Compute the boundary points in image space + >>> gcps = generate_gcps(camera, boundary) # Generate GCPs using the boundary in lat/lon + >>> vrt = generate_vrt((100,100), gcps, 'my_original_image.img') # Create the vrt + >>> # Then optionally, warp the VRT to render to disk + >>> warp_options = gdal.WarpOptions(format='VRT', dstNodata=0) + >>> gdal.Warp(some_output.tif, vrt, options=warp_options) + """ xsize, ysize = raster_size vrt = r'''<VRTDataset rasterXSize="{{ xsize }}" rasterYSize="{{ ysize }}"> <Metadata/> @@ -303,6 +339,6 @@ def generate_vrt(camera, raster_size, fpath, outpath=None, no_data_value=0): 'fpath':fpath, 'no_data_value':no_data_value} template = jinja2.Template(vrt) - tmp = template.render(context) - warp_options = gdal.WarpOptions(format='VRT', dstNodata=0) - gdal.Warp(outname, tmp, options=warp_options) + return template.render(context) + + diff --git a/recipe/meta.yml b/recipe/meta.yml new file mode 100644 index 0000000000000000000000000000000000000000..2e2bbbae35d22c0468da4fd5fa63f729be9d6c0f --- /dev/null +++ b/recipe/meta.yml @@ -0,0 +1,38 @@ +{% set data = load_setup_py_data() %} + +package: + name: knoten + version: {{ data.get('version') }} + +source: + git_url: https://github.com/USGS-Astrogeology/knoten + +build: + number : {{ GIT_DESCRIBE_NUMBER }} + string: dev + script: "{{ PYTHON }} -m pip install . --no-deps -vv" + +extra: + channels: + - conda-forge + +requirements: + host: + - pip + - python + - setuptools + run: + - python + - csmapi + - gdal + - jinja + - numpy + - requests + - scipy + +test: + imports: + - knoten + + + diff --git a/setup.py b/setup.py index 02074f820b6a8a8f27b49e3a122067bf3c7cbc27..1030f2c33ea7a943ebd419ba917059857bc339c7 100644 --- a/setup.py +++ b/setup.py @@ -1,6 +1,6 @@ from setuptools import setup, find_packages setup( - name='spaetzle', + name='knoten', version='0.1.0', long_description='', packages=find_packages(),