0

I have a pandas dataframe where each row is a bounding box. i'd like to crop my raster by each feature.

       xmin       xmax        ymin        ymax  \
1  398296.6  398301.46  3280997.15  3280998.96   

                                                 box  
1  NEON_D03_OSBS_DP1_398000_3280000_classified_po...  

I first created a function to create a geojson polygon from each row.

def data2geojson(df):
    features = []
    insert_features = lambda X: features.append(
            {"type": "polygon",
                 "coordinates": 
                 [[(float(X["xmin"]),float(X["ymin"])),
                     (float(X["xmax"]),float(X["ymin"])),
                     (float(X["xmax"]),float(X["ymax"])),
                     (float(X["xmin"]),float(X["ymax"])),
                     (float(X["xmin"]),float(X["ymin"]))]]}
        )

    df.apply(insert_features, axis=1)
    return geojson.FeatureCollection(features)

such that the output of each row looks like

{'type': 'polygon', 'coordinates': [[(398534.16, 3280995.59), (398539.84, 3280995.59), (398539.84, 3280998.99), (398534.16, 3280998.99), (398534.16, 3280995.59)]]}

The idea would be then to iterate through features to create the crops

out_image,out_transform=mask(img.tile,[features[0]],crop=True)

This yields no error message, but just kills the python kernel.

My initial thought is that my format of the geojson is incorrect? I was following rasterio masking with GeoJSON from file issue

System info

3.6.5 | packaged by conda-forge | (default, Mar 30 2018, 00:00:55) 
[GCC 4.2.1 Compatible Apple LLVM 6.1.0 (clang-602.0.53)]

rasterio.__version__
'0.36.0'

Edit #1: Looking at other answer, it appears we are looking for a geojson polygon, not feature. Same behavior, simplified the function.

Edit #2: Closed the polygon by providing the initial point as the last coordinate.

Edit #3, in case there was worry that my installation was invalid. I cloned the rasterio repo and ran this example. No issues. From that example, it appears that the features should be a dict, not a geojson polygon directly, which makes sense looking at the source. Revised. No improvement.

  • Shouldn't polygon be closed? – user30184 May 15 '18 at 19:54
  • no change, but edited to reflect that try. – bw4sz May 15 '18 at 20:04
1

It was only evident to me after cloning the rasterio repo and trying itself. The docs say

    GeoJSON-like dict representation of polygons that will be used to
    create the mask.

My mistake was literally trying to make a geojson object. I succeeded when I just made a dict that looked like a geojson object...just as the docs suggested. Taken literally, it works.

def data2geojson(df):
    features = []
    insert_features = lambda X: features.append(
            {"type": "Polygon",
                 "coordinates": 
                 [[(float(X["xmin"]),float(X["ymin"])),
                     (float(X["xmax"]),float(X["ymin"])),
                     (float(X["xmax"]),float(X["ymax"])),
                     (float(X["xmin"]),float(X["ymax"])),
                     (float(X["xmin"]),float(X["ymin"]))]]}
        )

Note that case matters! "Polygon"

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.