Merge pull request #15 from acpaquette/master

import pyproj

import requests

import scipy.stats

from singledispatch import singledispatch

from plio.io.io_gdal import GeoDataset

class NumpyEncoder(json.JSONEncoder):

    def default(self, obj):

            if plugin.canModelBeConstructedFromISD(isd, model_name):

                return plugin.constructModelFromISD(isd, model_name)

@singledispatch

def generate_ground_point(dem, image_pt, camera):

    '''

    Generates a longitude, latitude, and altitude coordinate for a given x, y

    pixel coordinate

    Parameters

    ----------

    dem : float or object

          Either a float that represents the height above the datum or a

          GeoDataset object generated from Plio off of a Digital Elevation

          Model (DEM)

    image_pt : tuple

               Pair of x, y coordinates in pixel space

    camera : object

             USGSCSM camera model object

    max_its : int, optional

              Maximum number of iterations to go through if the height does not

              converge

    tolerance : float, optional

                Number of decimal places to solve to when solving for height

    Returns

    -------

    : object

      CSM EcefCoord containing the newly computed lon, lat, and alt values

      corresponding to the original image_pt coordinates

    '''

    if not isinstance(image_pt, csmapi.ImageCoord):

        # Support a call where image_pt is in the form (x,y)

        image_pt = csmapi.ImageCoord(*image_pt)

    return camera.imageToGround(image_pt, dem)

@generate_ground_point.register(GeoDataset)

def _(dem, image_pt, camera, max_its = 20, tolerance = 0.001):

    if not isinstance(image_pt, csmapi.ImageCoord):

        # Support a call where image_pt is in the form (x,y)

        image_pt = csmapi.ImageCoord(*image_pt)

    intersection = generate_ground_point(0.0, image_pt, camera)

    iterations = 0

    semi_major, semi_minor = get_radii(camera)

    ecef = pyproj.Proj(proj='geocent', a=semi_major, b=semi_minor)

    lla = pyproj.Proj(proj='latlon', a=semi_major, b=semi_minor)

    while iterations != max_its:

        lon, lat, alt = pyproj.transform(ecef, lla, intersection.x, intersection.y, intersection.z)

        px, py = dem.latlon_to_pixel(lon, lat)

        height = dem.read_array(1, [px, py, 1, 1])[0][0]

        next_intersection = camera.imageToGround(image_pt, float(height))

        dist = max(abs(intersection.x - next_intersection.x),

            abs(intersection.y - next_intersection.y),

            abs(intersection.z - next_intersection.z))

        intersection = next_intersection

        iterations += 1

        if dist < tolerance:

            break

    return intersection

def generate_boundary(isize, npoints=10):

    '''

    Generates a bounding box given a camera model

    return boundary

def generate_latlon_boundary(camera, boundary, radii=None):

def generate_latlon_boundary(camera, boundary, dem=0.0, radii=None, **kwargs):

    '''

    Generates a latlon bounding box given a camera model

    alts : list

           List of altitude values

    '''

    if radii is None:

        semi_major, semi_minor = get_radii(camera)

    else:

    for i, b in enumerate(boundary):

        # Could be potential errors or warnings from imageToGround

        try:

            gnd = camera.imageToGround(csmapi.ImageCoord(*b), 0)

            gnd = generate_ground_point(dem, b, camera, **kwargs)

        except:

            pass

             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

    return gcps

def generate_latlon_footprint(camera, boundary, radii=None):

def generate_latlon_footprint(camera, boundary, dem=0.0, radii=None, **kwargs):

    '''

    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

    else:

        semi_major, semi_minor = radii

    lons, lats, _ = generate_latlon_boundary(camera, boundary, radii=radii)

    lons, lats, _ = generate_latlon_boundary(camera, boundary, dem=dem, radii=radii, **kwargs)

    # Transform coords from -180, 180 to 0, 360

    # Makes crossing the meridian easier to identify