Commit 3ef72c90 authored by Lauren Adoram-Kershner's avatar Lauren Adoram-Kershner Committed by Kelvin Rodriguez
Browse files

Bug fixes in Demo notebook and for Demo notebook (#346)

* general bug fixes for Demo

* Demo notebook bug fix

* updated autocnet version num

* tests updated to reflect data pvl format of nodes

* oops
parent 09ba2960
Loading
Loading
Loading
Loading
+9 −16
Original line number Diff line number Diff line
@@ -1534,7 +1534,7 @@ WHERE points.active = True AND measures.active=TRUE AND measures.jigreject=FALSE
        elif os.path.exists(filelist):
            filelist = io_utils.file_to_list(filelist)
        else:
            warning.warn('Unable to parse the passed filelist')
            warnings.warn('Unable to parse the passed filelist')

        for f in filelist:
            # Create the nodes in the graph. Really, this is creating the
@@ -1572,23 +1572,21 @@ WHERE points.active = True AND measures.active=TRUE AND measures.jigreject=FALSE
        intersect the user provided polygon (the LINESTRING) in the given spatial reference system
        (SRID), 949900.

        "SELECT * FROM Images WHERE ST_INTERSECTS(footprint_latlon, ST_Polygon(ST_GeomFromText('LINESTRING(159 10, 159 11, 160 11, 160 10, 159 10)'),949900)) = TRUE"
from_database
        SELECT * FROM Images WHERE ST_INTERSECTS(footprint_latlon, ST_Polygon(ST_GeomFromText('LINESTRING(159 10, 159 11, 160 11, 160 10, 159 10)'),949900)) = TRUE

        ## Select from a specific orbit
        This example selects those images that are from a particular orbit. In this case,
        the regex string pulls all P##_* orbits and creates a graph from them. This method
        does not guarantee that the graph is fully connected.

        "SELECT * FROM Images WHERE (split_part(path, '/', 6) ~ 'P[0-9]+_.+') = True"

        SELECT * FROM Images WHERE (split_part(path, '/', 6) ~ 'P[0-9]+_.+') = True
        """
        composite_query = """WITH
	i as ({})
SELECT i1.id as i1_id,i1.path as i1_path, i2.id as i2_id, i2.path as i2_path
FROM
	i as i1, i as i2

        composite_query = '''WITH i as ({}) SELECT i1.id
        as i1_id,i1.path as i1_path, i2.id as i2_id, i2.path as i2_path
        FROM i  as i1, i as i2
        WHERE ST_INTERSECTS(i1.footprint_latlon, i2.footprint_latlon) = TRUE
AND i1.id < i2.id""".format(query_string)
        AND i1.id < i2.id'''.format(query_string)

        session = Session()
        res = session.execute(composite_query)
@@ -1634,8 +1632,3 @@ AND i1.id < i2.id""".format(query_string)
                session.rollback()
        session.commit()
        session.close()




            
+9 −9
Original line number Diff line number Diff line
@@ -185,12 +185,13 @@ def place_points_in_overlap(nodes, geom, dem=None, cam_type="csm",

    valid = compgeom.distribute_points_in_geom(geom)
    if not valid:
        raise ValueError('Failed to distribute points in overlap')
        warnings.warn('Failed to distribute points in overlap')
        return []

    # Grab the source image. This is just the node with the lowest ID, nothing smart.
    source = nodes[0]
    nodes.remove(source)
    source_camera = source.camera
    source_camera = source["data"].camera
    for v in valid:
        lon = v[0]
        lat = v[1]
@@ -217,28 +218,27 @@ def place_points_in_overlap(nodes, geom, dem=None, cam_type="csm",

        point.measures.append(Measures(sample=ssample,
                                       line=sline,
                                       imageid=source['node_id'],
                                       serial=source.isis_serial,
                                       imageid=source['data']['node_id'],
                                       serial=source['data'].isis_serial,
                                       measuretype=3))


        for i, dest in enumerate(nodes):
            if cam_type == "csm":
                dic = dest.camera.groundToImage(gnd)
                dic = dest['data'].camera.groundToImage(gnd)
                dline, dsample = dic.line, dic.samp
            if cam_type == "isis":
                dline, dsample = isis.groud_to_image(dest["data"]["image_path"], lat, lon)

            dx, dy, _ = iterative_phase(ssample, sline, dsample, dline,
                                        source.geodata, dest.geodata,
                                        source['data'].geodata, dest['data'].geodata,
                                        **iterative_phase_kwargs)
            if dx is not None or dy is not None:
                point.measures.append(Measures(sample=dx,
                                               line=dy,
                                               imageid=dest['node_id'],
                                               serial=dest.isis_serial,
                                               imageid=dest['data']['node_id'],
                                               serial=dest['data'].isis_serial,
                                               measuretype=3))
        if len(point.measures) >= 2:
            points.append(point)
    return points
+27 −27
Original line number Diff line number Diff line
@@ -8,26 +8,26 @@ import csmapi
@patch('autocnet.cg.cg.distribute_points_in_geom', return_value=[(0, 0), (5, 5), (10, 10)])
def test_place_points_in_overlap(point_distributer, phase_matcher):
    # Mock setup
    first_node = MagicMock()
    first_node.camera = MagicMock()
    first_node.camera.groundToImage.return_value = csmapi.ImageCoord(1.0, 0.0)
    first_node.isis_serial = '1'
    first_node.__getitem__.return_value = 1
    second_node = MagicMock()
    second_node.camera = MagicMock()
    second_node.camera.groundToImage.return_value = csmapi.ImageCoord(1.0, 1.0)
    second_node.isis_serial = '2'
    second_node.__getitem__.return_value = 2
    third_node = MagicMock()
    third_node.camera = MagicMock()
    third_node.camera.groundToImage.return_value = csmapi.ImageCoord(0.0, 1.0)
    third_node.isis_serial = '3'
    third_node.__getitem__.return_value = 3
    fourth_node = MagicMock()
    fourth_node.camera = MagicMock()
    fourth_node.camera.groundToImage.return_value = csmapi.ImageCoord(0.0, 0.0)
    fourth_node.isis_serial = '4'
    fourth_node.__getitem__.return_value = 4
    first_node = {'data':MagicMock()}
    first_node['data'].camera = MagicMock()
    first_node['data'].camera.groundToImage.return_value = csmapi.ImageCoord(1.0, 0.0)
    first_node['data'].isis_serial = '1'
    first_node['data'].__getitem__.return_value = 1
    second_node = {'data':MagicMock()}
    second_node['data'].camera = MagicMock()
    second_node['data'].camera.groundToImage.return_value = csmapi.ImageCoord(1.0, 1.0)
    second_node['data'].isis_serial = '2'
    second_node['data'].__getitem__.return_value = 2
    third_node = {'data':MagicMock()}
    third_node['data'].camera = MagicMock()
    third_node['data'].camera.groundToImage.return_value = csmapi.ImageCoord(0.0, 1.0)
    third_node['data'].isis_serial = '3'
    third_node['data'].__getitem__.return_value = 3
    fourth_node = {'data':MagicMock()}
    fourth_node['data'].camera = MagicMock()
    fourth_node['data'].camera.groundToImage.return_value = csmapi.ImageCoord(0.0, 0.0)
    fourth_node['data'].isis_serial = '4'
    fourth_node['data'].__getitem__.return_value = 4
    dem = MagicMock()
    dem.latlon_to_pixel.return_value = (1.0, 1.0)
    dem.read_array.return_value = [[0.0]]
@@ -48,13 +48,13 @@ def test_place_points_in_overlap(point_distributer, phase_matcher):
    point_distributer.assert_called_with(Polygon([(0, 0), (0, 10), (10, 10), (10, 0)]))
    dem.latlon_to_pixel.assert_called()
    dem.read_array.assert_called()
    first_node.camera.groundToImage.assert_called()
    second_node.camera.groundToImage.assert_called()
    third_node.camera.groundToImage.assert_called()
    fourth_node.camera.groundToImage.assert_called()
    first_node['data'].camera.groundToImage.assert_called()
    second_node['data'].camera.groundToImage.assert_called()
    third_node['data'].camera.groundToImage.assert_called()
    fourth_node['data'].camera.groundToImage.assert_called()
    phase_matcher.assert_any_call(0.0, 1.0, 1.0, 1.0,
                                  first_node.geodata, second_node.geodata, size=71)
                                  first_node['data'].geodata, second_node['data'].geodata, size=71)
    phase_matcher.assert_any_call(0.0, 1.0, 1.0, 0.0,
                                  first_node.geodata, third_node.geodata, size=71)
                                  first_node['data'].geodata, third_node['data'].geodata, size=71)
    phase_matcher.assert_any_call(0.0, 1.0, 0.0, 0.0,
                                  first_node.geodata, fourth_node.geodata, size=71)
                                  first_node['data'].geodata, fourth_node['data'].geodata, size=71)
+2 −4
Original line number Diff line number Diff line
%% Cell type:markdown id: tags:

# AutoCNet Demo Notebook

%% Cell type:markdown id: tags:

This first cell is largely book keeping and environment setup. The second line places a configuration file into the environment that provides URLs, paths, and login information for the services that AutoCNet uses, as well as information about the spatial reference system the project is going to use.

Lines 4-6 get the USGS Community Sensor Model plugin loaded and ready for use.

%% Cell type:code id: tags:

``` python
import os
os.environ['autocnet_config'] = 'config/sample.yml'
os.environ['PROJ_LIB'] = '/home/ladoramkershner/miniconda3/env/autocnet_local/share/proj' #point to you local environment path
os.environ['PROJ_LIB'] = '/home/ladoramkershner/miniconda3/envs/autocnet_local/share/proj' #point to you local environment path

import ctypes
from ctypes.util import find_library
ctypes.CDLL(find_library('usgscsm'))

from autocnet.graph.network import NetworkCandidateGraph
```

%% Cell type:markdown id: tags:

The primary way to use AutoCNet is through the CandidateGraph object. In this demo, the derived, NetworkCandidateGraph is used. This object has an identical interface to the CandidateGraph. The difference lies in where the data are stored. In the CandidateGraph everything is stored in memory and all processing occurs in serial. On the NetworkCandidateGraph, data are stored in a database and processing occurs either in serial or on a compute cluster.

Below, the `from_filelist` method is used to perform an initial database populate. A few things are happening here behind the scenes:

* The database (name specified in the configuration file) is being created if it does not already exist. This creates all of the tables, relationships, and triggers.
* The `images` table is being populated with metadata about the images in the file list including a lat/lon footprint.
* The `cameras` table is being populated with a CSM compliant state string.

For a one-and-done style project, this cell should be run once.

%% Cell type:code id: tags:

``` python
# First run
import glob
ncg = NetworkCandidateGraph.from_filelist(glob.glob('/scratch/jlaura/elysium_subset/cal/*.cub'))
```

%% Cell type:markdown id: tags:

The cell below is the primary mechanism for accessing an existing project. Here, the `from_database` class method is used. This method makes a spatial query to the `images` table and builds a graph object where the nodes are images and the edges linking nodes indicate that the footprints overlap. This call makes no modifications to the database.

%% Cell type:code id: tags:

``` python
# On subsequent runs
ncg = NetworkCandidateGraph().from_database()
```

%% Cell type:markdown id: tags:

The database currently contains a populated `images` table and a populated `cameras` table. In the current pipeline style flow, we now need to compute polygons generated by the overlaps of footprints. It is within those polygons that correspondences can be found. We run two SQL queries (found in the `sql` directory).

The first query computes the n-wise overlapping polygons and the second query identifies those polygons that contributed to the overlapping geometry. The picture below illustrates what the result looks like.

![overlap](overlap.png)

Each polygon (A, B, AB, ABC, BC) will be an independent row in the database. The `intersections` column will be populated with an array where the values in the array are the `id`s of the image that contributed to the particular overlap.

%% Cell type:code id: tags:

``` python
# Run the 2 sql commands in the sql directory of the autocnet repo to compute overlaps and get the overlap arrays populates
```

%% Cell type:markdown id: tags:

Once overlaps have been computed, points are placed into the overlapping geometries. This process adds points to the `points` table and measures to the `measures` table. These points and measures are synonymous with the points and measures in an ISIS control network. This code is making use of the `images` and `cameras` tables that have been previously populated.

%% Cell type:code id: tags:

``` python
# This block is used to compute the overlapping polygon components and then place points into them.
from autocnet.spatial.overlap import place_points_in_overlaps

# Place points
place_points_in_overlaps(ncg, height=-3000) # This value is good for elysium, but needs to be more granularly parameterizable
# A bad height value results in very poor results... The height is height above (below) the sphere (the aeroid).
# To generalize this, we would spawn a new cluster job for each geomety and pull a height dynamically from from reference
place_points_in_overlaps(ncg)
```

%% Cell type:markdown id: tags:

Next, the points/measures are converted to be pairwise matches between images. Instead of a point-centric representation where a point has many measures, this representation is image-to-image centric where a pair of images has some set of shared correspondences. This command results in the `matches` table being populated with the pairwise matches between edges. Using the above example, all of the points/measures placed in the overlapping ABC polygon are decomposed into matches between AB, AC, and BC. The data are identical, but the representation has changed.

We do this because now we want to use classic computer vision techniques that generally operate best with pairwise representations.

%% Cell type:code id: tags:

``` python
# This block converts the points into matches
for s, d, e in ncg.edges(data='data'):  # intentionally in a loop so this doesn't spawn a cluster job
    e.network_to_matches()

```

%% Cell type:markdown id: tags:

Now that a pairwise representation exists, we apply a standard CV technique (computation of the fundamental matrix (F)). This operation updates the `edges` table to add the F matrix and a mask to the `masks` column indicating whether a particular match has been flagged as a blunder by the ransac procedure. This code is making use of the `matches` table that has been previously populated.

%% Cell type:code id: tags:

``` python
# This block computes the fundamental matrices
ncg.compute_fundamental_matrices(method='ransac', maskname='fundamental')
```

%% Output

    /home/jlaura/autocnet/autocnet/transformation/fundamental_matrix.py:310: UserWarning: F Computation Failed.
      warnings.warn("F Computation Failed.")

%% Cell type:markdown id: tags:

Finally, the next three cells are a janky heuristic designed to use the F matrix to remove blunders in the `points` and `measures` tables. This is an aggregation step whereby we seek to determine if a measure should be made inactive. This code makes use of the `edges` table (`masks` column) and updates the `measures` table (`active` column).

These three cells are not well integrated into the NetworkCandidateGraph and provide the best view of the type of operations happening under the hood in the previous steps.

%% Cell type:code id: tags:

``` python
# This block converts the points into matches
counters = []
for s, d, e in ncg.edges(data='data'):  # intentionally in a loop so this doesn't spawn a cluster job
    counters.append(e.mask_to_counter('fundamental'))
```

%% Cell type:code id: tags:

``` python
from collections import Counter
from autocnet import Session
from autocnet.io.db.model import Measures

aggregate = sum(counters, Counter())

# Now I need to take the output here and then look in qnet to see wtf is going on. Do we have a threshold here
# for blowing away bad stuff? If so, where? I should probably normalize all of these too based on the number of other
# images that they exist in. In other words, count/n-images
to_pop = []
session = Session()
for k, v in aggregate.items():
    pid = session.query(Measures).filter(Measures.id == k).first().pointid
    nimages = len(session.query(Measures).filter(Measures.pointid == pid).all())
    outlier_ratio = v / nimages  # This is the metric to test on (maybe?) - the ratio of the # of times the measure is
                                  # flagged as bad to the number of measures associated with the point.
    # These are rules that are going to need testing / vetting. Are these appropriate values?
    if outlier_ratio <= 0.5 or (outlier_ratio <= 0.5 and nimages == 2):
        to_pop.append(k)
    else:
        aggregate[k] = v / len(session.query(Measures).filter(Measures.pointid == pid).all())
for k in to_pop:
    aggregate.pop(k)
```

%% Cell type:code id: tags:

``` python
session = Session()
make_inactive = list(aggregate.keys())
session.query(Measures).filter(Measures.id.in_(make_inactive)).update({'active':False}, synchronize_session='fetch')
session.commit()
```

%% Cell type:code id: tags:

``` python
ncg.to_isis('/scratch/jlaura/elysium_subset/demo.net')
```

%% Cell type:code id: tags:

``` python
```
+1 −1
Original line number Diff line number Diff line
@@ -20,7 +20,7 @@ def setup_package():

    setup(
        name = "autocnet",
        version = '0.2.5',
        version = '0.2.6',
        author = "Jay Laura",
        author_email = "jlaura@usgs.gov",
        description = ("I/O API to support planetary data formats."),