sgillies

import rasterio

with rasterio.open(raster_filename, 'r') as src:
    image = src.read_band(1) # read

features = rasterio.features.polygonize(image) # no I/O!

with fiona.open(vector_filename, 'w', **meta) as dst:
    dst.writerecords(features) # write

sgillies

/******************************************************************************
 * $Id$
 * Project:  GDAL
 * Purpose:  Raster to Polygon Converter
 * Author:   Frank Warmerdam, [email protected]
 *
 ******************************************************************************
 * Copyright (c) 2008, Frank Warmerdam
 *
 * Permission is hereby granted, free of charge, to any person obtaining a

sgillies

sgillies

>>> import numpy
>>> from scipy.ndimage.morphology import grey_dilation
>>> a = numpy.zeros((7, 7))
>>> a[3,3] = 1
>>> a
array([[ 0.,  0.,  0.,  0.,  0.,  0.,  0.],
       [ 0.,  0.,  0.,  0.,  0.,  0.,  0.],
       [ 0.,  0.,  0.,  0.,  0.,  0.,  0.],
       [ 0.,  0.,  0.,  1.,  0.,  0.,  0.],
       [ 0.,  0.,  0.,  0.,  0.,  0.,  0.],

sgillies

import rasterio

with rasterio.open('input.tif') as src:
    source = src.read_band(1)

def window_transform(transform, offset):
    result = transform[:]
    result[3] = transform[3] + offset[0]*transform[5]
    result[0] = transform[0] + offset[1]*transform[1]
    return result

sgillies

sgillies

sgillies

...
#GPX 	GPX 	Yes 	Yes 	Yes (read support needs libexpat)
    ("GPX", "raw"),
...

sgillies

import rasterio
import numpy as np

src = rasterio.open('/Users/sean/code/decloud/envs/1_16/output/0_q.tif')

# Fraction of elements with value > 100 at full resolution.
q = src.read_band(1)
ones = np.ones(q.shape)
print "Full, %d pixels" % (q.shape[0]*q.shape[1])
print sum(ones[q>100])/(q.shape[0]*q.shape[1])

sgillies