Merge pull request #1 from worldbank/andres-update

Resolving conflicts

README.md

@@ -16,7 +16,7 @@ Future releases can be built from source, but pip will contain the most recent s
 
 ## Contributing
 
-Please refer to the World Bank's Github [Contributing](docs/Contributing.md) guidelines.
+Please refer to the World Bank's Github [Contributing](docs/CONTRIBUTING.md) guidelines.
 
 ## Code of Conduct
 

src/GOSTrocks/ghslMisc.py

@@ -6,7 +6,7 @@
 import xarray as xr
 import matplotlib.pyplot as plt
 
-from GOSTRocks.misc import tPrint
+from GOSTrocks.misc import tPrint
 
 curPath = os.path.realpath(os.path.abspath(os.path.split(inspect.getfile( inspect.currentframe() ))[0]))
 if not curPath in sys.path:

src/GOSTrocks/misc.py

@@ -5,9 +5,11 @@
 import numpy as np
 
 from math import ceil
-from shapely.geometry import Point, Polygon
+from shapely.geometry import Point, Polygon, box
 
-wgs84 = {'init':'epsg:4326'}
+from pyproj import CRS
+from pyproj.aoi import AreaOfInterest
+from pyproj.database import query_utm_crs_info
 
 def loggingInfo(lvl = logging.INFO):
     """ Set logging settings to info (default) and print useful information
@@ -17,7 +19,6 @@ def loggingInfo(lvl = logging.INFO):
     """
     logging.basicConfig(format='%(asctime)s:%(levelname)s: %(message)s', level=lvl)
 
-
 def tPrint(s):
     '''prints the time along with the message'''
     print("%s\t%s" % (time.strftime("%H:%M:%S"), s))
@@ -212,4 +213,22 @@ def project_UTM(inD):
         letter = '7'
     outUTM = '32%s%s' % (letter, ll_utm[2])
     return(inD.to_crs({'init': 'epsg:%s' % outUTM}))
+
+def get_utm(gdf, bbox=True):
+
+    if bbox:
+        center = box(*gdf.total_bounds).centroid
+    else:
+        center = gdf.unary_union.centroid
+    utm_crs_list = query_utm_crs_info(
+    datum_name="WGS 84",
+    area_of_interest=AreaOfInterest(
+        west_lon_degree=center.x,
+        south_lat_degree=center.y,
+        east_lon_degree=center.x,
+        north_lat_degree=center.y
+        ),
+    )
+    utm_crs = CRS.from_epsg(utm_crs_list[0].code)
+    return(utm_crs)
 

src/GOSTrocks/osmMisc.py

@@ -5,7 +5,7 @@
 import pandas as pd
 
 from shapely.geometry import box
-from GOSTRocks import misc
+from GOSTrocks import misc
 
 # Highway features are reclassified to 4 OSMLR classes for simplification and standardization
 #   https://mapzen.com/blog/osmlr-2nd-technical-preview/

src/GOSTrocks/rasterMisc.py

@@ -16,9 +16,6 @@
 from rasterio.io import MemoryFile
 from contextlib import contextmanager
 
-import seaborn as sns
-sns.set(font_scale=1.5, style="whitegrid")
-
 curPath = os.path.realpath(os.path.abspath(os.path.split(inspect.getfile( inspect.currentframe() ))[0]))
 if not curPath in sys.path:
     sys.path.append(curPath)
@@ -139,7 +136,7 @@ def project_raster(srcRst, dstCrs, output_raster=''):
 
     return([dstRst, kwargs])
 
-def clipRaster(inR, inD, outFile='', crop=True):
+def clipRaster(inR, inD, outFile=None, crop=True):
     ''' Clip input raster
     
     :param inR: rasterio object to clip
@@ -166,17 +163,18 @@ def getFeatures(gdf):
         tD = gpd.GeoDataFrame([[1]], geometry=[box(*inD.total_bounds)])
 
     coords = getFeatures(tD)
-    out_img, out_transform = mask(inR, shapes=coords, crop=True)
+    out_img, out_transform = mask(inR, shapes=coords, crop=True, all_touched=True)
     out_meta.update({"driver": "GTiff",
                      "height": out_img.shape[1],
                      "width": out_img.shape[2],
                      "transform": out_transform})
-    if outFile != '':
+    if outFile:
         with rasterio.open(outFile, "w", **out_meta) as dest:
             dest.write(out_img)
     return([out_img, out_meta])
 
-def rasterizeDataFrame(inD, outFile='', idField='', templateRaster='', templateMeta = '', nCells=0, res=0, mergeAlg="REPLACE", re_proj=False):
+def rasterizeDataFrame(inD, outFile, idField='', templateRaster='', templateMeta = '',
+                       nCells=0, res=0, mergeAlg="REPLACE", re_proj=False, nodata=np.nan):
     """ Convert input geopandas dataframe into a raster file
 
     :param inD: input data frame to rasterize
@@ -252,16 +250,17 @@ def rasterizeDataFrame(inD, outFile='', idField='', templateRaster='', templateM
         else:
             crs = inD.crs
         cMeta = {'count':1, 'crs': crs, 'dtype':inD['VALUE'].dtype, 'driver':'GTiff',
-                 'transform':nTransform, 'height':height, 'width':width}
+                 'transform':nTransform, 'height':height, 'width':width, 'nodata':nodata}
     shapes = ((row.geometry,row.VALUE) for idx, row in inD.iterrows())
-    burned = features.rasterize(shapes=shapes, out_shape=(cMeta['height'], cMeta['width']), transform=nTransform, dtype=cMeta['dtype'], merge_alg=mAlg)
-    try:
-        with rasterio.open(outFile, 'w', **cMeta) as out:
-            out.write_band(1, burned)
-        return({'meta':cMeta, 'vals': burned})
-    except:
-        print("Error writing raster")
-        return({'meta':cMeta, 'vals': burned})
+    burned = features.rasterize(shapes=shapes, out_shape=(cMeta['height'], cMeta['width']),
+                                transform=nTransform, dtype=cMeta['dtype'], merge_alg=mAlg, fill=nodata)
+    if outFile:
+        try:
+            with rasterio.open(outFile, 'w', **cMeta) as out:
+                out.write_band(1, burned)
+        except:
+            print("Error writing raster")
+    return({'meta':cMeta, 'vals': burned})
 
 def polygonizeArray(data, curRaster):
     """ Convert input array (data) to a geodataframe