Skip to content

Commit

Permalink
Merge branch 'master' into Dianne_Sept28
Browse files Browse the repository at this point in the history
  • Loading branch information
fourndo committed Jan 28, 2019
2 parents 0cebf85 + 6a745ed commit e3f5c79
Show file tree
Hide file tree
Showing 31 changed files with 1,081 additions and 1,642 deletions.
2 changes: 1 addition & 1 deletion .gitignore
Expand Up @@ -64,6 +64,6 @@ target/
# Windows

# notebooks
Notebooks/Output/*
# Notebooks/Output/*

!somefolder/.gitkeep
2 changes: 1 addition & 1 deletion .travis.yml
Expand Up @@ -21,7 +21,7 @@ before_install:

install:
- conda env create -f environment.yml
- source activate Toolkit-environment
- source activate geotoolkit
- export MPLBACKEND="agg"
- pip install -r requirements_dev.txt
- python setup.py install
Expand Down
72 changes: 58 additions & 14 deletions GeoToolkit/Mag/DataIO.py
Expand Up @@ -2,7 +2,7 @@
import scipy as sp
from . import (Simulator, MathUtils)
from scipy.spatial import cKDTree
from SimPEG.Utils import mkvc
from matplotlib.contour import QuadContourSet
import matplotlib.pyplot as plt
import gdal
import osr
Expand Down Expand Up @@ -81,7 +81,7 @@ def gridCC(self):

X, Y = np.meshgrid(self.hx, self.hy)

self._gridCC = np.c_[mkvc(X), mkvc(Y)]
self._gridCC = np.c_[X.flatten(order='F'), Y.flatten(order='F')]

return self._gridCC

Expand All @@ -108,7 +108,7 @@ def gridPadded(self):
else:

if np.any(np.isnan(self.values)):
self.indNan = np.isnan(mkvc(self.values))
self.indNan = np.isnan(self.values.flatten(order='F'))
grid = self.valuesFilled
else:
grid = self.values
Expand Down Expand Up @@ -160,18 +160,33 @@ def values(self):
def valuesFilled(self):

if getattr(self, '_valuesFilled', None) is None:
values = mkvc(self.values)
indVal = np.where(~self.indNan)[0]
values = self.values.flatten(order='F')

# Do a minimum curvature extrapolation
isNan = np.isnan(values)

# Get the real values on the outer edge
indVal = np.where(~isNan)[0]
tree = cKDTree(self.gridCC[indVal, :])
dists, indexes = tree.query(self.gridCC[self.indNan, :])
dists, indexes = tree.query(self.gridCC[isNan, :])

uInd = np.unique(indVal[indexes])

xyz = self.gridCC[uInd, :]

_, values[self.indNan] = MathUtils.minCurvatureInterp(
xyz, values[uInd], xyzOut=self.gridCC[self.indNan, :])
_, values[isNan] = MathUtils.minCurvatureInterp(
xyz, values[uInd],
xyzOut=self.gridCC[isNan, :])

# If there are still NDV, just do a neareest neighbour

while np.isnan(values).sum() > 0:
isNan = np.isnan(values)
# Get the real values on the outer edge
indVal = np.where(~isNan)[0]
tree = cKDTree(self.gridCC[indVal, :])
dists, indexes = tree.query(self.gridCC[isNan, :])
values[isNan] = values[indVal[indexes]]

self._valuesFilled = values.reshape(self.values.shape, order='F')

Expand Down Expand Up @@ -236,7 +251,7 @@ def derivativeX(self):
fhxd_pad[self.npady:-self.npady, self.npadx:-self.npadx])

if self.indNan is not None:
derivX = mkvc(derivX)
derivX = derivX.flatten(order='F')

derivX[self.indNan] = np.nan
derivX = derivX.reshape(self.values.shape, order='F')
Expand All @@ -257,7 +272,7 @@ def derivativeY(self):
fhyd_pad[self.npady:-self.npady, self.npadx:-self.npadx])

if self.indNan is not None:
derivY = mkvc(derivY)
derivY = derivY.flatten(order='F')

derivY[self.indNan] = np.nan
derivY = derivY.reshape(self.values.shape, order='F')
Expand All @@ -277,7 +292,7 @@ def firstVertical(self):
fhzd_pad[self.npady:-self.npady, self.npadx:-self.npadx])

if self.indNan is not None:
firstVD = mkvc(firstVD)
firstVD = firstVD.flatten(order='F')

firstVD[self.indNan] = np.nan
firstVD = firstVD.reshape(self.values.shape, order='F')
Expand Down Expand Up @@ -338,7 +353,7 @@ def setRTP(self, isRTP):
rtp_pad[self.npady:-self.npady, self.npadx:-self.npadx])

if self.indNan is not None:
rtp = mkvc(rtp)
rtp = rtp.flatten(order='F')

rtp[self.indNan] = np.nan
rtp = rtp.reshape(self.values.shape, order='F')
Expand Down Expand Up @@ -374,7 +389,7 @@ def upwardContinuation(self, z=0):
self._valuesFilledUC = zUpw.copy()

if self.indNan is not None:
zUpw = mkvc(zUpw)
zUpw = zUpw.flatten(order='F')

zUpw[self.indNan] = np.nan
zUpw = zUpw.reshape(self.values.shape, order='F')
Expand Down Expand Up @@ -560,8 +575,20 @@ def exportShapefile(

}

if isinstance(polylines, QuadContourSet):

temp, attributes = [], []
for segs, level in zip(polylines.allsegs, polylines.levels):

for poly in segs:

temp += [poly]
attributes += [level]

polylines = temp

with fiona.open(
saveAs + '.shp', 'w', 'ESRI Shapefile', schema, crs=crs
saveAs + '.shp', 'w', driver='ESRI Shapefile', schema=schema, crs=crs
) as c:

# If there are multiple geometries, put the "for" loop here
Expand All @@ -577,6 +604,23 @@ def exportShapefile(
res['geometry'] = mapping(pline)
c.write(res)

# function to generate .prj file information using spatialreference.org
def getWKT_PRJ(epsg_code):
import urllib
# access projection information
wkt = urllib.request.urlopen("http://spatialreference.org/ref/epsg/{0}/prettywkt/".format(epsg_code))
# remove spaces between charachters
remove_spaces = wkt.read().replace(b" ", b"")
# place all the text on one line
output = remove_spaces.replace(b"\n", b"")
return output

# create the .prj file
prj = open(saveAs + ".prj", "w")
# call the function and supply the epsg code
epsg = getWKT_PRJ(str(int(EPSGcode)))
prj.write(epsg.decode("utf-8"))
prj.close()

def fetchData(
path="./assets/Search/", checkDir=False, file="",
Expand Down

0 comments on commit e3f5c79

Please sign in to comment.