In [65]:
show(src.read(1))
Out[65]:
<AxesSubplot:>
Okay--I admit this is slightly complicated, but doing GIS in the wild can be that way sometimes!
There's a bunch of ways you may want to clip something, so consider this a primer on what it really takes to clip a raster.
Essentially, you need:
Let's do it!
import rasterio as rio
from rasterio.plot import show
from rasterio.mask import mask #this is the masking aka clip module
import matplotlib.pyplot as plt
import geopandas as gpd #to read the input shapefile
import shapely #to transform extent into coordinate geometries
import json
import os
Important: only run this once or you'll wind up in the wrong directory.
os.chdir('../')
os.getcwd()
First, let's read in the image:
src = rio.open('workshopdata/NAIP_Campus.tif')
src = rio.open('workshopdata/NAIP_Campus.tif')
Let's take a peek:
show(src.read(1))
show(src.read(1))
<AxesSubplot:>
Let's clip so that we're only working with the quad area...
Read in the campus shapefile:
shape = gpd.read_file('workshopdata/UCB_MainCampus_Boundaries.shp')
shape = gpd.read_file('workshopdata/UCB_MainCampus_Boundaries.shp')
Let's take a look:
shape.plot(figsize = (12,12))
shape[22:23].plot(figsize = (12,12))
<AxesSubplot:>
Okay, but if we're going to clip with this, it needs to be in the same CRS as the image, this is easy to change using geopandas:
src.crs
src.crs
CRS.from_epsg(26913)
Change it to EPSG 26913
shape = shape.to_crs(26913)
shape = shape.to_crs(26913)
Now, let's get the geometry of the southwest part of main campus:
We'll use GeoPandas to access the geometry of that particular feature.
swQuad = shape.geometry[22]
swQuad = shape.geometry[22]
print(swQuad)
POLYGON ((477547.8159370111 4428256.238502387, 477603.2191063287 4428311.991142652, 477610.260153695 4428313.372457464, 477663.4512336894 4428360.636259837, 477676.9024611069 4428375.471382552, 477687.7264066462 4428392.210929826, 477695.4980979763 4428422.068033623, 477694.6469493901 4428437.489914686, 477693.0019507 4428695.097440739, 477835.9787919673 4428694.515272426, 477916.8683550697 4428683.440503193, 477921.0539890607 4428429.005324838, 477916.7407354645 4428311.39728243, 477915.3693515466 4428244.320595865, 477916.9353471698 4428134.864784565, 477901.0624630528 4428038.386414838, 477886.135831317 4427999.586153438, 477877.0739491967 4427987.72881689, 477862.4903698746 4427970.25427815, 477820.7270960809 4427936.344001541, 477819.520040813 4427935.282557529, 477818.3290782706 4427934.203029216, 477817.154509462 4427933.105614842, 477815.9966359601 4427931.990612612, 477814.8557593446 4427930.858320715, 477813.7319806968 4427929.708938523, 477812.6258015168 4427928.542763111, 477811.5372229281 4427927.359994397, 477810.4666464725 4427926.160930035, 477809.4142737513 4427924.945868777, 477808.3804063456 4427923.715108819, 477807.3651453439 4427922.468849537, 477806.3687928739 4427921.207489084, 477805.3915499914 4427919.931226264, 477804.4337182654 4427918.640359277, 477803.4954993092 4427917.335186888, 477802.5769947746 4427916.016008425, 477801.6785062311 4427914.683122078, 477800.8002352882 4427913.336826622, 477799.9423841297 4427911.977520775, 477799.1050538309 4427910.605403904, 477798.288645932 4427909.220773641, 477797.4930621551 4427907.823930446, 477796.7187051666 4427906.41537187, 477795.9656754816 4427904.99519733, 477795.2341758478 4427903.563905496, 477794.524307339 4427902.121695743, 477793.8362721362 4427900.668966787, 477793.1702718479 4427899.206017399, 477792.5263087277 4427897.733247428, 477791.9047837424 4427896.250854544, 477791.3054997821 4427894.759339692, 477790.7289577812 4427893.258899976, 477790.1749600645 4427891.749936378, 477789.6438082064 4427890.232747095, 477789.1357049404 4427888.707830824, 477788.6505508677 4427887.175288074, 477788.1885492817 4427885.635717503, 477787.7499012398 4427884.089317914, 477787.3346089792 4427882.536489155, 477786.9427741635 4427880.977530549, 477786.5743990258 4427879.412841951, 477786.2297857127 4427877.842821515, 477785.9088359431 4427876.267769698, 477785.6116519319 4427874.688085789, 477785.3384352833 4427873.104068547, 477785.0889883284 4427871.516118938, 477784.8636126365 4427869.924535177, 477784.6622110528 4427868.329817632, 477784.4848846761 4427866.73216567, 477784.3317362638 4427865.132078532, 477784.2026675509 4427863.529856669, 477784.0978801451 4427861.925798854, 477784.0171769317 4427860.320406016, 477764.0059635132 4427843.964513178, 477691.1781734904 4427836.376911336, 477557.1093193134 4427822.112503104, 477538.7085833196 4427820.357116409, 477524.0411849411 4427821.059182097, 477511.666475973 4427825.259645761, 477500.959943642 4427832.135803921, 477489.6752551041 4427846.450752284, 477461.6442778364 4427895.972177809, 477412.4114870515 4427981.081018358, 477400.8901550272 4428003.800306964, 477383.0070603646 4428033.443299215, 477353.9488704999 4428071.291930779, 477352.2769082641 4428073.062926424, 477350.588161236 4428074.817822033, 477348.8825288863 4428076.5565182, 477347.1604099455 4428078.27881276, 477345.4219032535 4428079.984505221, 477343.6671082046 4428081.673495064, 477341.8962236026 4428083.345581235, 477340.1093488453 4428085.000663205, 477338.3067826954 4428086.638539372, 477336.4886239966 4428088.259009245, 477334.6549727068 4428089.862072263, 477332.8061275925 4428091.447526808, 477330.9420869607 4428093.015073, 477329.0631512661 4428094.564809112, 477327.1695187447 4428096.09643413, 477325.261289362 4428097.609947495, 477323.3385613879 4428099.105048764, 477321.4017341178 4428100.58163571, 477319.4508069815 4428102.039608387, 477317.4859793524 4428103.478865695, 477315.5074505804 4428104.899306553, 477313.515418912 4428106.300629951, 477311.5100842722 4428107.682834766, 477309.4916460204 4428109.045819901, 477307.4602029956 4428110.389384883, 477305.4159545616 4428111.713428614, 477303.3592000392 4428113.017849445, 477301.2899383128 4428114.302447462, 477299.2085692268 4428115.567220406, 477297.1150922154 4428116.812068317, 477295.0099054495 4428118.036689062, 477292.8930094854 4428119.241182595, 477290.7647030932 4428120.425347318, 477288.6250851009 4428121.588982739, 477286.4744559654 4428122.732187124, 477284.3130139228 4428123.854659474, 477282.1409594635 4428124.956498613, 477279.9583908561 4428126.037404109, 477277.7655080289 4428127.097374827, 477275.5627097142 4428128.136208594, 477273.3499964655 4428129.154005374, 477271.1276670465 4428130.150563559, 477268.8958208594 4428131.125782616, 477266.6548571929 4428132.079560339, 477209.4152747923 4428156.313806572, 477187.7892743901 4428166.421779409, 477177.194418984 4428179.985711107, 477176.8040512021 4428180.660874734, 477176.4346424642 4428181.347816408, 477176.0866880798 4428182.04573362, 477175.7604851356 4428182.754124873, 477175.4562296178 4428183.472289305, 477175.1744173985 4428184.199524378, 477174.9151445121 4428184.935129776, 477174.6787063504 4428185.678304122, 477174.4653994193 4428186.428445954, 477174.2754197185 4428187.18485441, 477174.1087633058 4428187.946829758, 477173.9658255418 4428188.713570041, 477173.8467018818 4428189.484275002, 477173.7512890003 4428190.258345429, 477173.6799822467 4428191.034979374, 477173.6325777611 4428191.813478222, 477173.6093709322 4428192.593040587, 477173.6101579024 4428193.372967857, 477173.6351341018 4428194.152459215, 477173.6840956695 4428194.930816034, 477173.7571386293 4428195.707338025, 477173.8541585258 4428196.481226042, 477173.9750508963 4428197.251680962, 477174.1196124483 4428198.01810412, 477174.2878381203 4428198.779595856, 477174.4795251755 4428199.535657492, 477174.6944686298 4428200.285390479, 477174.9324646249 4428201.028096208, 477175.1933093024 4428201.763076073, 477175.4766988347 4428202.489632011, 477175.7824299322 4428203.20716538, 477176.1101982062 4428203.914878164, 477176.4595998718 4428204.612072872, 477176.8304321958 4428205.298250816, 477177.2222902677 4428205.972514586, 477177.6348719501 4428206.634466018, 477178.0675729743 4428207.283308783, 477178.5201900377 4428207.91844423, 477178.9922199606 4428208.539275399, 477179.4831601199 4428209.145305285, 477179.9927066957 4428209.735835834, 477180.5201565699 4428210.310271215, 477181.0651076504 4428210.868213817, 477181.6270567491 4428211.409066685, 477182.2054007218 4428211.932233408, 477182.7996375024 4428212.437316949, 477183.4092644748 4428212.923820315, 477184.0335790899 4428213.391247617, 477184.6720787196 4428213.839101858, 477185.3240619444 4428214.267087089, 477185.9890261393 4428214.674706299, 477186.6662693243 4428215.061563584, 477187.3550900713 4428215.427362978, 477188.0548857894 4428215.771608061, 477188.7649561894 4428216.094202786, 477189.4846981218 4428216.394550761, 477190.2134118468 4428216.672655917, 477190.9502948502 4428216.928022946, 477191.6947468037 4428217.160555249, 477192.4460662788 4428217.369956874, 477193.203452451 4428217.556032399, 477193.9662044613 4428217.718585831, 477194.7336220084 4428217.857521146, 477195.5050042305 4428217.972642345, 477196.2795514342 4428218.063953933, 477197.0565627486 4428218.131259918, 477387.9783805712 4428217.047942437, 477395.0443542381 4428216.497783144, 477402.1120582443 4428215.97090608, 477409.181492589 4428215.467311241, 477416.2525567608 4428214.986899227, 477423.3250513832 4428214.52977112, 477430.3991763868 4428214.095925801, 477437.4745319249 4428213.685365513, 477444.5513179202 4428213.298089139, 477451.6292339248 4428212.933998397, 477458.708381023 4428212.593292645, 477465.7886581301 4428212.275772531, 477472.8698664483 4428211.981639099, 477479.9520048493 4428211.710692424, 477480.9932983483 4428211.70213766, 477482.0345240208 4428211.717075197, 477483.0751820586 4428211.755507837, 477484.1146732539 4428211.817538918, 477485.1524972225 4428211.90307129, 477486.1881535951 4428212.012007798, 477487.2210420331 4428212.144251862, 477488.2507632548 4428212.299905678, 477489.2766163896 4428212.478773261, 477490.2982015976 4428212.68085687, 477491.3149174095 4428212.905859977, 477492.3262645806 4428213.153885357, 477493.3318421338 4428213.424735343, 477494.3309492086 4428213.718213939, 477495.323184836 4428214.034123464, 477496.3080491976 4428214.372466736, 477497.2850413595 4428214.733046638, 477498.2535598584 4428215.115566659, 477499.2133036837 4428215.519828551, 477500.1635725404 4428215.945736295, 477501.1039648892 4428216.392992257, 477502.033980354 4428216.861499269, 477502.953216275 4428217.350759781, 477544.6755198449 4428252.896918824, 477547.8159370111 4428256.238502387))
Perfect!
Now finally, we can run the mask function to get the clipped image. mask() will return two things: the image, and it's affine transformation. You must set crop to True--without it, the image will be the same size but with no values beyond our clip area. Setting crop to True gets rid of those areas:
out_img, out_transform = mask(src, [swQuad], crop=True)
** One other caveat: wrap the swQuad variable in square brackets so rasterio knows to treat it like an iterable.
out_img, out_transform = mask(src, [swQuad], crop=True)
Let's take a look...
show(out_img)
show(out_img)
<AxesSubplot:>
Okay, maybe you don't want it clipped to the exact geometry... maybe you would rather have it by the extent of the geometry?
Okay, let's do it!
What we'll need is the bounding box geometry of our swQuad shape. We can use geopandas' envelope property:
envelope = swQuad.envelope
envelope
envelope = swQuad.envelope
envelope
print it...
print(envelope)
POLYGON ((477173.6093709322 4427820.357116409, 477921.0539890607 4427820.357116409, 477921.0539890607 4428695.097440739, 477173.6093709322 4428695.097440739, 477173.6093709322 4427820.357116409))
Now, let's do it again, this time with the envelope as the geometry object:
out_img, out_transform = mask(src, [envelope], crop=True)
out_img, out_transform = mask(src,[envelope], crop=True)
Check it out:
show(out_img)
show(out_img)
<AxesSubplot:>
Now, let's get the original profile
src.profile
src.profile
{'driver': 'GTiff', 'dtype': 'uint8', 'nodata': None, 'width': 9590, 'height': 12260, 'count': 4, 'crs': CRS.from_epsg(26913), 'transform': Affine(0.6, 0.0, 473124.0,
0.0, -0.6, 4434942.0), 'tiled': False, 'interleave': 'pixel'}
Copy it to a variable:
profile = src.profile.copy()
(Note, this is the same as the kwargs variable we used in the last notebook)
profile = src.profile.copy()
Grab the height and the width of the output image:
height = out_img.shape[1]
width = out_img.shape[2]
height = out_img.shape[1]
width = out_img.shape[2]
Update the profile. Since we have a new transform, a new width & height, we need to use the .update method to take the old profile and update certain elements:
profile.update(transform = out_transform,
width = width,
height = height,
photometric = 'rgb',
alpha = 'no')
(more on photometric & alpha in the next notebook)
profile.update(transform = out_transform, width = width, height = height, photometric = 'rgb', alpha = 'no')
Open a new empty dataset:
new = rio.open('NAIP_swQuad.tif', 'w', **profile)
new = rio.open('NAIP_swQuad.tif', 'w', **profile)
Write it to disc:
new.write(out_img)
new.write(out_img)
Close it...
new.close()
new.close()