Skip to content
Snippets Groups Projects
Commit 8c234ffe authored by Tyler Thatcher's avatar Tyler Thatcher Committed by GitHub
Browse files

Merge pull request #47 from acpaquette/ipf 

isis2socet
parents 44acc278 76450e6e
No related branches found
No related tags found
No related merge requests found
Hide whitespace changes
Inline Side-by-side
Showing
with 581 additions and 360 deletions
bin/isisnet2socet 0 → 100644
+ 102
− 0
View file @ 8c234ffe
#!/usr/bin/env python
import argparse
import os
import pandas as pd
from plio.io.io_bae import save_gpf, save_ipf
from plio.spatial.transformations import apply_isis_transformations
import plio.io.io_controlnetwork as cn
import plio.io.isis_serial_number as sn
def parse_args():
parser = argparse.ArgumentParser()
# Add args here
parser.add_argument('cnet_file', help='Path to an isis control network.')
parser.add_argument('e_radius', type=float, help='The semimajor radius of a given target.')
parser.add_argument('p_radius', type=float, help='The semiminor radius of a given target.')
parser.add_argument('cub_path', help='Path to the cub files associated with a control network.')
parser.add_argument('cub_extension', help='Extension for all cubes.')
parser.add_argument('cub_list', help='Path to a list file of all cubes being used')
parser.add_argument('out_gpf', help='Path to save location of gpf file and new ipf files.')
parser.add_argument('--adjusted', help='Flag for saving apriori values or adjusted values',
default=False, required = False)
return parser.parse_args()
def main(args):
# Create cub dict to map ipf to cub
df = cn.from_isis(args.cnet_file)
e_radius = args.e_radius
p_radius = e_radius * (1 - args.p_radius)
cub_path = args.cub_path
extension = args.cub_extension
with open(args.cub_list, 'r') as f:
lines = f.readlines()
cub_list = [cub.replace('\n', '') for cub in lines]
out_gpf = args.out_gpf
adjusted_flag = args.adjusted
# Create cub dict to map ipf to cub
cub_dict = {i: i + extension for i in cub_list}
# Create serial dict to match serial to ipf
serial_dict = {sn.generate_serial_number(os.path.join(cub_path, i + extension)): i for i in cub_list}
# Remove duplicate columns
# There are better ways to do this but pandas was not having it
columns = []
column_index = []
for i, column in enumerate(list(df.columns)):
if column not in columns:
column_index.append(i)
columns.append(column)
df = df.iloc[:, column_index]
# Begin translation
# Remap the ISIS columns to socet column names
column_map = {'id': 'pt_id', 'line': 'l.', 'sample': 's.',
'lineResidual': 'res_l', 'sampleResidual': 'res_s', 'type': 'known',
'aprioriLatitudeSigma': 'sig0', 'aprioriLongitudeSigma': 'sig1', 'aprioriRadiusSigma': 'sig2',
'linesigma': 'sig_l', 'samplesigma': 'sig_s', 'ignore': 'stat'}
# Depending on the adjusted flag, set the renames for columns appropriately
if adjusted_flag:
column_map['adjustedY'] = 'lat_Y_North'
column_map['adjustedX'] = 'long_X_East'
column_map['adjustedZ'] = 'ht'
else:
column_map['aprioriY'] = 'lat_Y_North'
column_map['aprioriX'] = 'long_X_East'
column_map['aprioriZ'] = 'ht'
df.rename(columns = column_map, inplace=True)
apply_isis_transformations(df, e_radius, p_radius, serial_dict, extension, cub_path)
# Save the ipf(s)
save_ipf(df, os.path.split(out_gpf)[0])
# Get the first record from each group as there all the same, put them
# into a list, and sort it
points = [int(i[1].index[0]) for i in df.groupby('pt_id')]
points.sort()
# Set the gpf_df to only the values we need and do a small rename
gpf_df = df.iloc[points].copy()
gpf_df.rename(columns = {'pt_id': 'point_id'}, inplace=True)
# Save the gpf
save_gpf(gpf_df, out_gpf)
if __name__ == '__main__':
main(parse_args())
bin/socet2isis → bin/socetnet2isis
+ 17
− 26
View file @ 8c234ffe
#!/usr/bin/env python
import argparse
import os
import sys
import argparse
import warnings
import csv
import numpy as np
from plio.examples import get_path
from plio.io.io_bae import read_atf, read_gpf, read_ipf
from plio.spatial.transformations import *
from plio.spatial.transformations import apply_socet_transformations, serial_numbers
import plio.io.io_controlnetwork as cn
import pandas as pd
......@@ -19,9 +16,9 @@ def parse_args():
# Add args here
parser.add_argument('at_file', help='Path to the .atf file for a project.')
parser.add_argument('cub_file_path', help='Path to cube files related to ipf files.')
parser.add_argument('cub_ipf_map', help='Path to map file for all ipfs and cubes.')
parser.add_argument('--outpath', help='Directory for the control network to be output to.',
required = False)
parser.add_argument('extension', help='Extension for all cubes being used.')
parser.add_argument('target_name', help='Name of the target body used in the control net')
parser.add_argument('--outpath', help='Directory for the control network to be output to.')
return parser.parse_args()
......@@ -37,10 +34,6 @@ def main(args):
else:
outpath = os.path.split(at_file)[0]
with open(args.cub_ipf_map) as cub_ipf_map:
reader = csv.reader(cub_ipf_map, delimiter = ',')
image_dict = dict([(row[0], row[1]) for row in reader])
# Read in and setup the atf dict of information
atf_dict = read_atf(at_file)
......@@ -60,33 +53,31 @@ def main(args):
point_diff = ipf_pt_idx.difference(gpf_pt_idx)
if len(point_diff) != 0:
warnings.warn("The following points found in ipf files missing from gpf file: " +
"\n\n{}\n\n".format("\n".join(point_diff)) +
"Continuing, but these points will be missing from the control " +
"network.", stacklevel=3)
warnings.warn("The following points found in ipf files missing from gpf file: \n\n{}. \
\n\nContinuing, but these points will be missing from the control network".format(list(point_diff)))
# Merge the two dataframes on their point id columns
socet_df = ipf_df.merge(gpf_df, left_on='pt_id', right_on='point_id')
# Apply the transformations
apply_transformations(atf_dict, socet_df)
apply_socet_transformations(atf_dict, socet_df)
# Define column remap for socet dataframe
column_remap = {'l.': 'y', 's.': 'x',
'res_l': 'LineResidual', 'res_s': 'SampleResidual', 'known': 'Type',
'lat_Y_North': 'AprioriY', 'long_X_East': 'AprioriX', 'ht': 'AprioriZ',
'sig0': 'AprioriLatitudeSigma', 'sig1': 'AprioriLongitudeSigma',
'sig2': 'AprioriRadiusSigma'}
column_map = {'pt_id': 'id', 'l.': 'y', 's.': 'x',
'res_l': 'lineResidual', 'res_s': 'sampleResidual', 'known': 'Type',
'lat_Y_North': 'aprioriY', 'long_X_East': 'aprioriX', 'ht': 'aprioriZ',
'sig0': 'aprioriLatitudeSigma', 'sig1': 'aprioriLongitudeSigma', 'sig2': 'aprioriRadiusSigma',
'sig_l': 'linesigma', 'sig_s': 'samplesigma'}
# Rename the columns using the column remap above
socet_df.rename(columns = column_remap, inplace=True)
images = pd.unique(socet_df['ipf_file'])
socet_df.rename(columns = column_map, inplace=True)
# Build an image and serial dict assuming the cubes will be named as the IPFs are
image_dict = {i: i + args.extension for i in pd.unique(socet_df['ipf_file'])}
serial_dict = serial_numbers(image_dict, cub_path)
# creates the control network
cn.to_isis(os.path.join(outpath, cnet_out + '.net'), socet_df, serial_dict)
cn.to_isis(os.path.join(outpath, cnet_out + '.net'), socet_df, serial_dict, targetname = args.target_name)
if __name__ == '__main__':
main(parse_args())
notebooks/Socet2ISIS.ipynb
+ 122
− 295
View file @ 8c234ffe
%% Cell type:code id: tags:
``` python
import os
import sys
from functools import singledispatch
import warnings
import pandas as pd
import numpy as np
import math
import pyproj
from plio.examples import get_path
from plio.io.io_bae import read_gpf, read_ipf
from collections import defaultdict
from plio.io.io_bae import read_gpf, read_ipf, read_atf, save_gpf, save_ipf
from plio.utils.utils import find_in_dict
from plio.spatial.transformations import apply_isis_transformations, apply_socet_transformations, serial_numbers
import plio.io.io_controlnetwork as cn
import plio.io.isis_serial_number as sn
```
%% Cell type:code id: tags:
``` python
# Reads a .atf file and outputs all of the
# .ipf, .gpf, .sup, .prj, and path to locate the
# .apf file (should be the same as all others)
def read_atf(atf_file):
with open(atf_file) as f:
# Extensions of files we want
files_ext = ['.prj', '.sup', '.ipf', '.gpf']
files_dict = []
files = defaultdict(list)
for line in f:
ext = os.path.splitext(line)[-1].strip()
# Check is needed for split as all do not have a space
if ext in files_ext:
# If it is the .prj file, it strips the directory away and grabs file name
if ext == '.prj':
files[ext].append(line.strip().split(' ')[1].split('\\')[-1])
# If the ext is in the list of files we care about, it addes to the dict
files[ext].append(line.strip().split(' ')[-1])
else:
# Adds to the dict even if not in files we care about
files[ext.strip()].append(line)
# Gets the base filepath
files['basepath'] = os.path.dirname(os.path.abspath(atf_file))
# Creates a dict out of file lists for GPF, PRJ, IPF, and ATF
files_dict = (dict(files_dict))
# Sets the value of IMAGE_IPF to all IPF images
files_dict['IMAGE_IPF'] = files['.ipf']
# Sets the value of IMAGE_SUP to all SUP images
files_dict['IMAGE_SUP'] = files['.sup']
# Sets value for GPF file
files_dict['GP_FILE'] = files['.gpf'][0]
# Sets value for PRJ file
files_dict['PROJECT'] = files['.prj'][0]
# Sets the value of PATH to the path of the ATF file
files_dict['PATH'] = files['basepath']
return files_dict
# converts columns l. and s. to isis
def line_sample_size(record, path):
with open(os.path.join(path, record['ipf_file'] + '.sup')) as f:
for i, line in enumerate(f):
if i == 2:
img_index = line.split('\\')
img_index = img_index[-1].strip()
img_index = img_index.split('.')[0]
if i == 3:
line_size = line.split(' ')
line_size = line_size[-1].strip()
assert int(line_size) > 0, "Line number {} from {} is a negative number: Invalid Data".format(line_size, record['ipf_file'])
if i == 4:
sample_size = line.split(' ')
sample_size = sample_size[-1].strip()
assert int(sample_size) > 0, "Sample number {} from {} is a negative number: Invalid Data".format(sample_size, record['ipf_file'])
break
line_size = int(line_size)/2.0 + record['l.'] + 1
sample_size = int(sample_size)/2.0 + record['s.'] + 1
return sample_size, line_size, img_index
# converts known to ISIS keywords
def known(record):
if record['known'] == 0:
return 'Free'
elif record['known'] == 1 or record['known'] == 2 or record['known'] == 3:
return 'Constrained'
# converts +/- 180 system to 0 - 360 system
def to_360(num):
return num % 360
# ocentric to ographic latitudes
def oc2og(dlat, dMajorRadius, dMinorRadius):
try:
dlat = math.radians(dlat)
dlat = math.atan(((dMajorRadius / dMinorRadius)**2) * (math.tan(dlat)))
dlat = math.degrees(dlat)
except:
print ("Error in oc2og conversion")
return dlat
# ographic to ocentric latitudes
def og2oc(dlat, dMajorRadius, dMinorRadius):
try:
dlat = math.radians(dlat)
dlat = math.atan((math.tan(dlat) / ((dMajorRadius / dMinorRadius)**2)))
dlat = math.degrees(dlat)
except:
print ("Error in og2oc conversion")
return dlat
# gets eRadius and pRadius from a .prj file
def get_axis(file):
with open(file) as f:
from collections import defaultdict
files = defaultdict(list)
for line in f:
ext = line.strip().split(' ')
files[ext[0]].append(ext[-1])
eRadius = float(files['A_EARTH'][0])
pRadius = eRadius * (1 - float(files['E_EARTH'][0]))
return eRadius, pRadius
# function to convert lat_Y_North to ISIS_lat
def lat_ISIS_coord(record, semi_major, semi_minor):
ocentric_coord = og2oc(record['lat_Y_North'], semi_major, semi_minor)
coord_360 = to_360(ocentric_coord)
return coord_360
# function to convert long_X_East to ISIS_lon
def lon_ISIS_coord(record, semi_major, semi_minor):
ocentric_coord = og2oc(record['long_X_East'], semi_major, semi_minor)
coord_360 = to_360(ocentric_coord)
return coord_360
def body_fix(record, semi_major, semi_minor, inverse=False):
"""
Parameters
----------
record : ndarray
(n,3) where columns are x, y, height or lon, lat, alt
"""
ecef = pyproj.Proj(proj='geocent', a=semi_major, b=semi_minor)
lla = pyproj.Proj(proj='latlon', a=semi_major, b=semi_minor)
if inverse:
lon, lat, height = pyproj.transform(ecef, lla, record[0], record[1], record[2])
return lon, lat, height
else:
y, x, z = pyproj.transform(lla, ecef, record[0], record[1], record[2])
return y, x, z
def socet2isis(at_file, cub_file_path, extension, target_name, outpath=None):
# Setup the at_file, path to cubes, and control network out path
at_file = at_file
cnet_out = os.path.split(os.path.splitext(at_file)[0])[1]
cub_path = cub_file_path
def ignore_toggle(record):
if record['stat'] == 0:
return True
if(outpath):
outpath = outpath
else:
return False
# TODO: Does isis cnet need a convariance matrix for sigmas? Even with a static matrix of 1,1,1,1
def compute_sigma_covariance_matrix(lat, lon, rad, latsigma, lonsigma, radsigma, semimajor_axis):
"""
Given geospatial coordinates, desired accuracy sigmas, and an equitorial radius, compute a 2x3
sigma covariange matrix.
Parameters
----------
lat : float
A point's latitude in degrees
lon : float
A point's longitude in degrees
rad : float
The radius (z-value) of the point in meters
latsigma : float
The desired latitude accuracy in meters (Default 10.0)
lonsigma : float
The desired longitude accuracy in meters (Default 10.0)
radsigma : float
The desired radius accuracy in meters (Defualt: 15.0)
semimajor_axis : float
The semi-major or equitorial radius in meters (Default: 1737400.0 - Moon)
Returns
-------
rectcov : ndarray
(2,3) covariance matrix
"""
lat = math.radians(lat)
lon = math.radians(lon)
# SetSphericalSigmasDistance
scaled_lat_sigma = latsigma / semimajor_axis
# This is specific to each lon.
scaled_lon_sigma = lonsigma * math.cos(lat) / semimajor_axis
# SetSphericalSigmas
cov = np.eye(3,3)
cov[0,0] = scaled_lat_sigma ** 2
cov[1,1] = scaled_lon_sigma ** 2
cov[2,2] = radsigma ** 2
# Approximate the Jacobian
j = np.zeros((3,3))
cosphi = math.cos(lat)
sinphi = math.sin(lat)
coslambda = math.cos(lon)
sinlambda = math.sin(lon)
rcosphi = rad * cosphi
rsinphi = rad * sinphi
j[0,0] = -rsinphi * coslambda
j[0,1] = -rcosphi * sinlambda
j[0,2] = cosphi * coslambda
j[1,0] = -rsinphi * sinlambda
j[1,1] = rcosphi * coslambda
j[1,2] = cosphi * sinlambda
j[2,0] = rcosphi
j[2,1] = 0.
j[2,2] = sinphi
mat = j.dot(cov)
mat = mat.dot(j.T)
rectcov = np.zeros((2,3))
rectcov[0,0] = mat[0,0]
rectcov[0,1] = mat[0,1]
rectcov[0,2] = mat[0,2]
rectcov[1,0] = mat[1,1]
rectcov[1,1] = mat[1,2]
rectcov[1,2] = mat[2,2]
return rectcov
# return np.array([[1.0, 1.0, 1.0], [1.0, 1.0, 1.0]])
def compute_cov_matrix(record, semimajor_axis):
cov_matrix = compute_sigma_covariance_matrix(record['lat_Y_North'], record['long_X_East'], record['ht'], record['sig0'], record['sig1'], record['sig2'], semimajor_axis)
return cov_matrix.ravel().tolist()
# applys transformations to columns
def apply_transformations(atf_dict, df):
prj_file = os.path.join(atf_dict['PATH'], atf_dict['PROJECT'])
eRadius, pRadius = get_axis(prj_file)
lla = np.array([[df['long_X_East']], [df['lat_Y_North']], [df['ht']]])
ecef = body_fix(lla, semi_major = eRadius, semi_minor = pRadius, inverse=False)
df['s.'], df['l.'], df['image_index'] = (zip(*df.apply(line_sample_size, path = atf_dict['PATH'], axis=1)))
df['known'] = df.apply(known, axis=1)
df['long_X_East'] = ecef[0][0]
df['lat_Y_North'] = ecef[1][0]
df['ht'] = ecef[2][0]
df['aprioriCovar'] = df.apply(compute_cov_matrix, semimajor_axis = eRadius, axis=1)
df['ignore'] = df.apply(ignore_toggle, axis=1)
outpath = os.path.split(at_file)[0]
def socet2isis(prj_file):
# Read in and setup the atf dict of information
atf_dict = read_atf(prj_file)
atf_dict = read_atf(at_file)
# Get the gpf and ipf files using atf dict
gpf_file = os.path.join(atf_dict['PATH'], atf_dict['GP_FILE']);
ipf_list = [os.path.join(atf_dict['PATH'], i) for i in atf_dict['IMAGE_IPF']]
# Read in the gpf file and ipf file(s) into seperate dataframes
gpf_df = read_gpf(gpf_file)
ipf_df = read_ipf(ipf_list)
# Check for differences between point ids using each dataframes
# point ids as a reference
gpf_pt_idx = pd.Index(pd.unique(gpf_df['point_id']))
ipf_pt_idx = pd.Index(pd.unique(ipf_df['pt_id']))
point_diff = ipf_pt_idx.difference(gpf_pt_idx)
if len(point_diff) != 0:
warnings.warn("The following points found in ipf files missing from gpf file: \n\n{}. \
\n\nContinuing, but these points will be missing from the control network".format(list(point_diff)))
# Merge the two dataframes on their point id columns
socet_df = ipf_df.merge(gpf_df, left_on='pt_id', right_on='point_id')
# Apply the transformations
apply_transformations(atf_dict, socet_df)
apply_socet_transformations(atf_dict, socet_df)
# Define column remap for socet dataframe
column_remap = {'l.': 'y', 's.': 'x',
'res_l': 'lineResidual', 'res_s': 'sampleResidual', 'known': 'Type',
'lat_Y_North': 'AprioriY', 'long_X_East': 'AprioriX', 'ht': 'AprioriZ',
'sig0': 'AprioriLatitudeSigma', 'sig1': 'AprioriLongitudeSigma', 'sig2': 'AprioriRadiusSigma',
'sig_l': 'linesigma', 'sig_s': 'samplesigma'}
column_map = {'pt_id': 'id', 'l.': 'y', 's.': 'x',
'res_l': 'lineResidual', 'res_s': 'sampleResidual', 'known': 'Type',
'lat_Y_North': 'aprioriY', 'long_X_East': 'aprioriX', 'ht': 'aprioriZ',
'sig0': 'aprioriLatitudeSigma', 'sig1': 'aprioriLongitudeSigma', 'sig2': 'aprioriRadiusSigma',
'sig_l': 'linesigma', 'sig_s': 'samplesigma'}
# Rename the columns using the column remap above
socet_df.rename(columns = column_remap, inplace=True)
socet_df.rename(columns = column_map, inplace=True)
# Build an image and serial dict assuming the cubes will be named as the IPFs are
image_dict = {i: i + extension for i in pd.unique(socet_df['ipf_file'])}
serial_dict = serial_numbers(image_dict, cub_path)
# creates the control network
cn.to_isis(os.path.join(outpath, cnet_out + '.net'), socet_df, serial_dict, targetname = targetname)
```
%% Cell type:code id: tags:
``` python
def isis2socet(cnet_path, eRadius, eccentricity, cub_path, extension, cub_list, out_gpf, adjusted_flag = False):
pRadius = eRadius * math.sqrt(1 - (eccentricity ** 2))
df = cn.from_isis(cnet_path)
# Create cub dict to map ipf to cub
cub_dict = {i: i + extension for i in cub_list}
# Create serial dict to match serial to ipf
serial_dict = {sn.generate_serial_number(os.path.join(cub_path, i + extension)): i for i in cub_list}
# Remove duplicate columns
# There are better ways to do this but pandas was not having it
columns = []
column_index = []
for i, column in enumerate(list(df.columns)):
if column not in columns:
column_index.append(i)
columns.append(column)
df = df.iloc[:, column_index]
# Begin translation
# Remap the ISIS columns to socet column names
column_map = {'id': 'pt_id', 'line': 'l.', 'sample': 's.',
'lineResidual': 'res_l', 'sampleResidual': 'res_s', 'type': 'known',
'aprioriLatitudeSigma': 'sig0', 'aprioriLongitudeSigma': 'sig1', 'aprioriRadiusSigma': 'sig2',
'linesigma': 'sig_l', 'samplesigma': 'sig_s', 'ignore': 'stat'}
# Depending on the adjusted flag, set the renames for columns appropriately
if adjusted_flag:
column_map['adjustedY'] = 'lat_Y_North'
column_map['adjustedX'] = 'long_X_East'
column_map['adjustedZ'] = 'ht'
else:
column_map['aprioriY'] = 'lat_Y_North'
column_map['aprioriX'] = 'long_X_East'
column_map['aprioriZ'] = 'ht'
# Return the socet dataframe to be converted to a control net
return socet_df
df.rename(columns = column_map, inplace=True)
# creates a dict of serial numbers with the cub being the key
def serial_numbers(images, path, extension):
serial_dict = dict()
for image in images:
snum = sn.generate_serial_number(os.path.join(path, image + extension))
snum = snum.replace('Mars_Reconnaissance_Orbiter', 'MRO')
serial_dict[image] = snum
return serial_dict
apply_isis_transformations(df, eRadius, pRadius, serial_dict, extension, cub_path)
# Save the ipf
save_ipf(df, os.path.split(out_gpf)[0])
# Get the first record from each group as there all the same, put them
# into a list, and sort it
points = [int(i[1].index[0]) for i in df.groupby('pt_id')]
points.sort()
# Set the gpf_df to only the values we need and do a small rename
gpf_df = df.iloc[points].copy()
gpf_df.rename(columns = {'pt_id': 'point_id'}, inplace=True)
# Save the gpf
save_gpf(gpf_df, out_gpf)
```
%% Cell type:code id: tags:
``` python
# Setup stuffs for the cub information namely the path and extension
path = '/path/where/cub/files/are/'
cub_path = '/Path/to/cubs'
# Extension of your cub files
extension = '.something.cub'
# Name of the target body
targetname = 'Mars'
extension = 'cub.-->extension<--'
# Path to atf file
atf_file = ('/path/to/atf/file')
atf_file = 'Path/to/socket/set/at_file.atf'
socet2isis(atf_file, cub_path, extension, targetname)
```
%% Cell type:code id: tags:
``` python
# Setup stuffs for the cub information namely the path and extension
# along with eRadius and eccentricity
cnet = "Path/to/control/network.net"
eRadius = 3.39619000000000e+006
eccentricity = 1.08339143554195e-001
cub_path = 'Path/to/cubs'
extension = 'cub.-->extension<--'
socet_df = socet2isis(atf_file)
# List of cubes to use
cub_list = ['D06_029601_1846_XN_04N224W',
'F05_037684_1857_XN_05N224W']
images = pd.unique(socet_df['ipf_file'])
out_gpf = "/Users/adampaquette/Desktop/InSightE09_XW.gpf"
serial_dict = serial_numbers(images, path, extension)
adjusted_flag = False
# creates the control network
cn.to_isis('/path/you/want/the/cnet/to/be/in/cn.net', socet_df, serial_dict)
isis2socet(cnet, eRadius, eccentricity, cub_path, extension, cub_list, out_gpf, adjusted_flag)
```
%% Cell type:code id: tags:
``` python
```
......
plio/examples/SocetSet/cub_map2.csv 0 → 100644
+ 2
− 0
View file @ 8c234ffe
D06_029601_1846_XN_04N224W,D06_029601_1846_XN_04N224W.8bit.cub
F05_037684_1857_XN_05N224W,F05_037684_1857_XN_05N224W.8bit.cub
plio/io/io_bae.py
+ 28
− 23
View file @ 8c234ffe
import json
import re
import os
import json
from collections import defaultdict
from functools import singledispatch
import numpy as np
......@@ -290,44 +291,48 @@ def read_atf(atf_file):
project
"""
with open(atf_file) as f:
files = []
ipf = []
sup = []
# Extensions of files we want
files_ext = ['.prj', '.sup', '.ipf', '.gpf']
files_dict = []
files = defaultdict(list)
# Grabs every PRJ, GPF, SUP, and IPF image from the ATF file
for line in f:
if line[-4:-1] == 'prj' or line[-4:-1] == 'gpf' or line[-4:-1] == 'sup' or line[-4:-1] == 'ipf' or line[-4:-1] == 'atf':
files.append(line)
ext = os.path.splitext(line)[-1].strip()
files = np.array(files)
# Check is needed for split as all do not have a space
if ext in files_ext:
# Creates appropriate arrays for certain files in the right format
for file in files:
file = file.strip()
file = file.split(' ')
# If it is the .prj file, it strips the directory away and grabs file name
if ext == '.prj':
files[ext].append(line.strip().split(' ')[1].split('\\')[-1])
# Grabs all the IPF files
if file[1].endswith('.ipf'):
ipf.append(file[1])
# If the ext is in the list of files we care about, it addes to the dict
files[ext].append(line.strip().split(' ')[-1])
# Grabs all the SUP files
if file[1].endswith('.sup'):
sup.append(file[1])
else:
files_dict.append(file)
# Adds to the dict even if not in files we care about
files[ext.strip()].append(line)
# Gets the base filepath
files['basepath'] = os.path.dirname(os.path.abspath(atf_file))
# Creates a dict out of file lists for GPF, PRJ, IPF, and ATF
files_dict = (dict(files_dict))
# Sets the value of IMAGE_IPF to all IPF images
files_dict['IMAGE_IPF'] = ipf
files_dict['IMAGE_IPF'] = files['.ipf']
# Sets the value of IMAGE_SUP to all SUP images
files_dict['IMAGE_SUP'] = sup
files_dict['IMAGE_SUP'] = files['.sup']
# Sets value for GPF file
files_dict['GP_FILE'] = files['.gpf'][0]
# Sets value for PRJ file
files_dict['PROJECT'] = files['.prj'][0]
# Sets the value of PATH to the path of the ATF file
files_dict['PATH'] = os.path.dirname(os.path.abspath(atf_file))
files_dict['PATH'] = files['basepath']
return files_dict
plio/io/io_controlnetwork.py
+ 6
− 5
View file @ 8c234ffe
......@@ -187,11 +187,12 @@ class IsisStore(object):
header_bytes = find_in_dict(pvl_header, 'HeaderBytes')
point_start_byte = find_in_dict(pvl_header, 'PointsStartByte')
version = find_in_dict(pvl_header, 'Version')
if version == 2:
point_attrs = [i for i in cnf._CONTROLPOINTFILEENTRYV0002.fields_by_name if i != 'measures']
measure_attrs = [i for i in cnf._CONTROLPOINTFILEENTRYV0002_MEASURE.fields_by_name]
self.point_attrs = [i for i in cnf._CONTROLPOINTFILEENTRYV0002.fields_by_name if i != 'measures']
self.measure_attrs = [i for i in cnf._CONTROLPOINTFILEENTRYV0002_MEASURE.fields_by_name]
cols = point_attrs + measure_attrs
cols = self.point_attrs + self.measure_attrs
cp = cnf.ControlPointFileEntryV0002()
self._handle.seek(header_start_byte)
......@@ -203,10 +204,10 @@ class IsisStore(object):
pts = []
for s in pbuf_header.pointMessageSizes:
cp.ParseFromString(self._handle.read(s))
pt = [getattr(cp, i) for i in point_attrs if i != 'measures']
pt = [getattr(cp, i) for i in self.point_attrs if i != 'measures']
for measure in cp.measures:
meas = pt + [getattr(measure, j) for j in measure_attrs]
meas = pt + [getattr(measure, j) for j in self.measure_attrs]
pts.append(meas)
df = IsisControlNetwork(pts, columns=cols)
df.header = pvl_header
......
plio/spatial/transformations.py
+ 302
− 9
View file @ 8c234ffe
import os
import pvl
import math
import pyproj
import numpy as np
import plio.io.isis_serial_number as sn
from plio.utils.utils import find_in_dict
def line_sample_size(record, path):
"""
......@@ -163,7 +167,7 @@ def get_axis(file):
files[ext[0]].append(ext[-1])
eRadius = float(files['A_EARTH'][0])
pRadius = eRadius * (1 - float(files['E_EARTH'][0]))
pRadius = eRadius * math.sqrt(1 - (float(files['E_EARTH'][0]) ** 2))
return eRadius, pRadius
......@@ -217,7 +221,57 @@ def lon_ISIS_coord(record, semi_major, semi_minor):
coord_360 = to_360(ocentric_coord)
return coord_360
def body_fix(record, semi_major, semi_minor):
def lat_socet_coord(record, semi_major, semi_minor):
"""
Function to convert lat_Y_North to ISIS_lat
Parameters
----------
record : object
Pandas series object
semi_major : float
Radius from the center of the body to the equater
semi_minor : float
Radius from the pole to the center of mass
Returns
-------
coord_360 : float
Converted latitude into ocentric space, and mapped
into 0 to 360
"""
ographic_coord = oc2og(record['lat_Y_North'], semi_major, semi_minor)
coord_180 = ((ographic_coord + 180) % 360) - 180
return coord_180
def lon_socet_coord(record, semi_major, semi_minor):
"""
Function to convert long_X_East to ISIS_lon
Parameters
----------
record : object
Pandas series object
semi_major : float
Radius from the center of the body to the equater
semi_minor : float
Radius from the pole to the center of mass
Returns
-------
coord_360 : float
Converted longitude into ocentric space, and mapped
into 0 to 360
"""
ographic_coord = oc2og(record['long_X_East'], semi_major, semi_minor)
coord_180 = ((ographic_coord + 180) % 360) - 180
return coord_180
def body_fix(record, semi_major, semi_minor, inverse = False, **kwargs):
"""
Transforms latitude, longitude, and height of a socet point into
a body fixed point
......@@ -241,10 +295,85 @@ def body_fix(record, semi_major, semi_minor):
"""
ecef = pyproj.Proj(proj='geocent', a=semi_major, b=semi_minor)
lla = pyproj.Proj(proj='latlon', a=semi_major, b=semi_minor)
lon, lat, height = pyproj.transform(lla, ecef, record['long_X_East'], record['lat_Y_North'], record['ht'])
return lon, lat, height
def apply_transformations(atf_dict, df):
if inverse:
lon, lat, height = pyproj.transform(ecef, lla, record[0], record[1], record[2], **kwargs)
return lon, lat, height
else:
y, x, z = pyproj.transform(lla, ecef, record[0], record[1], record[2], **kwargs)
return y, x, z
def stat_toggle(record):
if record['stat'] == 0:
return True
else:
return False
def apply_isis_transformations(df, eRadius, pRadius, serial_dict, extension, cub_path):
"""
Takes a atf dictionary and a socet dataframe and applies the necessary
transformations to convert that dataframe into a isis compatible
dataframe
Parameters
----------
df : object
Pandas dataframe object
eRadius : float
Equitorial radius
pRadius : float
Polar radius
serial_dict : dict
Dictionary mapping serials as keys to images as the values
extension : str
String extension of all cubes being used
cub_path : str
Path to all cubes being used
"""
# Convert from geocentered coords (x, y, z), to lat lon coords (latitude, longitude, alltitude)
ecef = np.array([[df['long_X_East']], [df['lat_Y_North']], [df['ht']]])
lla = body_fix(ecef, semi_major = eRadius, semi_minor = pRadius, inverse=True)
df['long_X_East'], df['lat_Y_North'], df['ht'] = lla[0][0], lla[1][0], lla[2][0]
# df['lat_Y_North'] = df.apply(lat_socet_coord, semi_major = eRadius, semi_minor = pRadius, axis=1)
# df['long_X_East'] = df.apply(lon_socet_coord, semi_major = eRadius, semi_minor = pRadius, axis=1)
# Update the stat fields and add the val field as it is just a clone of stat
df['stat'] = df.apply(ignore_toggle, axis = 1)
df['val'] = df['stat']
# Update the known field, add the ipf_file field for saving, and
# update the line, sample using data from the cubes
df['known'] = df.apply(reverse_known, axis = 1)
df['ipf_file'] = df['serialnumber'].apply(lambda serial_number: serial_dict[serial_number])
df['l.'], df['s.'] = zip(*df.apply(fix_sample_line, serial_dict = serial_dict,
extension = extension,
cub_path = cub_path, axis = 1))
# Add dummy for generic value setting
x_dummy = lambda x: np.full(len(df), x)
df['sig0'] = x_dummy(1)
df['sig1'] = x_dummy(1)
df['sig2'] = x_dummy(1)
df['res0'] = x_dummy(0)
df['res1'] = x_dummy(0)
df['res2'] = x_dummy(0)
df['fid_x'] = x_dummy(0)
df['fid_y'] = x_dummy(0)
df['no_obs'] = x_dummy(1)
df['fid_val'] = x_dummy(0)
def apply_socet_transformations(atf_dict, df):
"""
Takes a atf dictionary and a socet dataframe and applies the necessary
transformations to convert that dataframe into a isis compatible
......@@ -259,15 +388,24 @@ def apply_transformations(atf_dict, df):
Pandas dataframe object
"""
prj_file = os.path.join(atf_dict['PATH'], atf_dict['PROJECT'].split('\\')[-1])
prj_file = os.path.join(atf_dict['PATH'], atf_dict['PROJECT'])
eRadius, pRadius = get_axis(prj_file)
# df['lat_Y_North'] = df.apply(lat_ISIS_coord, semi_major = eRadius, semi_minor = pRadius, axis=1)
# df['long_X_East'] = df.apply(lon_ISIS_coord, semi_major = eRadius, semi_minor = pRadius, axis=1)
lla = np.array([[df['long_X_East']], [df['lat_Y_North']], [df['ht']]])
ecef = body_fix(lla, semi_major = eRadius, semi_minor = pRadius, inverse=False)
df['s.'], df['l.'], df['image_index'] = (zip(*df.apply(line_sample_size, path = atf_dict['PATH'], axis=1)))
df['known'] = df.apply(known, axis=1)
df['lat_Y_North'] = df.apply(lat_ISIS_coord, semi_major = eRadius, semi_minor = pRadius, axis=1)
df['long_X_East'] = df.apply(lon_ISIS_coord, semi_major = eRadius, semi_minor = pRadius, axis=1)
df['long_X_East'], df['lat_Y_North'], df['ht'] = zip(*df.apply(body_fix, semi_major = eRadius, semi_minor = pRadius, axis = 1))
df['long_X_East'] = ecef[0][0]
df['lat_Y_North'] = ecef[1][0]
df['ht'] = ecef[2][0]
df['aprioriCovar'] = df.apply(compute_cov_matrix, semimajor_axis = eRadius, axis=1)
df['stat'] = df.apply(stat_toggle, axis=1)
def serial_numbers(image_dict, path):
"""
......@@ -290,3 +428,158 @@ def serial_numbers(image_dict, path):
for key in image_dict:
serial_dict[key] = sn.generate_serial_number(os.path.join(path, image_dict[key]))
return serial_dict
# TODO: Does isis cnet need a convariance matrix for sigmas? Even with a static matrix of 1,1,1,1
def compute_sigma_covariance_matrix(lat, lon, rad, latsigma, lonsigma, radsigma, semimajor_axis):
"""
Given geospatial coordinates, desired accuracy sigmas, and an equitorial radius, compute a 2x3
sigma covariange matrix.
Parameters
----------
lat : float
A point's latitude in degrees
lon : float
A point's longitude in degrees
rad : float
The radius (z-value) of the point in meters
latsigma : float
The desired latitude accuracy in meters (Default 10.0)
lonsigma : float
The desired longitude accuracy in meters (Default 10.0)
radsigma : float
The desired radius accuracy in meters (Defualt: 15.0)
semimajor_axis : float
The semi-major or equitorial radius in meters (Default: 1737400.0 - Moon)
Returns
-------
rectcov : ndarray
(2,3) covariance matrix
"""
lat = math.radians(lat)
lon = math.radians(lon)
# SetSphericalSigmasDistance
scaled_lat_sigma = latsigma / semimajor_axis
# This is specific to each lon.
scaled_lon_sigma = lonsigma * math.cos(lat) / semimajor_axis
# SetSphericalSigmas
cov = np.eye(3,3)
cov[0,0] = math.radians(scaled_lat_sigma) ** 2
cov[1,1] = math.radians(scaled_lon_sigma) ** 2
cov[2,2] = radsigma ** 2
# Approximate the Jacobian
j = np.zeros((3,3))
cosphi = math.cos(lat)
sinphi = math.sin(lat)
cos_lmbda = math.cos(lon)
sin_lmbda = math.sin(lon)
rcosphi = rad * cosphi
rsinphi = rad * sinphi
j[0,0] = -rsinphi * cos_lmbda
j[0,1] = -rcosphi * sin_lmbda
j[0,2] = cosphi * cos_lmbda
j[1,0] = -rsinphi * sin_lmbda
j[1,1] = rcosphi * cos_lmbda
j[1,2] = cosphi * sin_lmbda
j[2,0] = rcosphi
j[2,1] = 0.
j[2,2] = sinphi
mat = j.dot(cov)
mat = mat.dot(j.T)
rectcov = np.zeros((2,3))
rectcov[0,0] = mat[0,0]
rectcov[0,1] = mat[0,1]
rectcov[0,2] = mat[0,2]
rectcov[1,0] = mat[1,1]
rectcov[1,1] = mat[1,2]
rectcov[1,2] = mat[2,2]
return rectcov
def compute_cov_matrix(record, semimajor_axis):
cov_matrix = compute_sigma_covariance_matrix(record['lat_Y_North'], record['long_X_East'], record['ht'], record['sig0'], record['sig1'], record['sig2'], semimajor_axis)
return cov_matrix.ravel().tolist()
def reverse_known(record):
"""
Converts the known field from an isis dataframe into the
socet known column
Parameters
----------
record : object
Pandas series object
Returns
-------
: str
String representation of a known field
"""
record_type = record['known']
if record_type == 0 or record_type == 2:
return 0
elif record_type == 1 or record_type == 3 or record_type == 4:
return 3
def fix_sample_line(record, serial_dict, extension, cub_path):
"""
Extracts the sample, line data from a cube and computes deviation from the
center of the image
Parameters
----------
record : dict
Dict containing the key serialnumber, l., and s.
serial_dict : dict
Maps serial numbers to images
extension : str
Extension for cube being looked at
cub_path : str
Path to a given cube being looked at
Returns
-------
new_line : int
new line deviation from the center
new_sample : int
new sample deviation from the center
"""
# Cube location to load
cube = pvl.load(os.path.join(cub_path, serial_dict[record['serialnumber']] + extension))
line_size = find_in_dict(cube, 'Lines')
sample_size = find_in_dict(cube, 'Samples')
new_line = record['l.'] - (int(line_size/2.0)) - 1
new_sample = record['s.'] - (int(sample_size/2.0)) - 1
return new_line, new_sample
def ignore_toggle(record):
"""
Maps the stat column in a record to 0 or 1 based on True or False
Parameters
----------
record : dict
Dict containing the key stat
"""
if record['stat'] == True:
return 0
else:
return 1
setup.py
+ 2
− 2
View file @ 8c234ffe
......@@ -36,7 +36,7 @@ def setup_package():
package_data={'plio' : list(examples) + ['data/*.db', 'data/*.py'] +\
['sqlalchemy_json/*.py', 'sqlalchemy_json/LICENSE']},
zip_safe=False,
scripts=['bin/socet2isis'],
scripts=['bin/socet2isis', 'bin/isis2socet'],
install_requires=[
'gdal',
'numpy',
......@@ -46,7 +46,7 @@ def setup_package():
'pandas',
'sqlalchemy',
'pyyaml',
'networkx',
'networkx',
'affine',
'scipy'],
classifiers=[
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment