@@ -1419,33 +1419,35 @@ <h3>Data Pre-Processing</h3>
14191419< p > These files were download from < a
14201420href ="https://dataverse.harvard.edu/dataset.xhtml?persistentId=doi:10.7910/DVN/YGIVCD "> Harvard
14211421Dataverse</ a > and converted to Cloud-Optimized GeoTIFFs using GDAL.</ p >
1422- < p > Example command
1423- < code > gdalwarp -of COG -co COMPRESS=DEFLATE -co PREDICTOR=2 -co NUM_THREADS=ALL_CPUS 2021_HasMask/LongNTL_2021.tif 2021.tif -te -180 -90 180 90 -dstnodata 0</ code > </ p >
1422+ < p > Example command</ p >
1423+ < pre > < code > gdalwarp -of COG 2021_HasMask/LongNTL_2021.tif 2021.tif \
1424+ -te -180 -90 180 90 -dstnodata 0 \
1425+ -co COMPRESS=DEFLATE -co PREDICTOR=2 -co NUM_THREADS=ALL_CPUS</ code > </ pre >
14241426< p > The resulting files are now hosted on a Google Cloud Storage
14251427bucket.</ p >
1426- < div class ="sourceCode " id ="cb82 "> < pre
1427- class ="sourceCode python "> < code class ="sourceCode python "> < span id ="cb82-1 "> < a href ="#cb82-1 " tabindex ="-1 "> </ a > start_year < span class ="op "> =</ span > < span class ="dv "> 2015</ span > </ span >
1428- < span id ="cb82-2 "> < a href ="#cb82-2 " tabindex ="-1 "> </ a > end_year < span class ="op "> =</ span > < span class ="dv "> 2020</ span > </ span > </ code > </ pre > </ div >
14291428< div class ="sourceCode " id ="cb83 "> < pre
1430- class ="sourceCode python "> < code class ="sourceCode python "> < span id ="cb83-1 "> < a href ="#cb83-1 " tabindex ="-1 "> </ a > </ span >
1431- < span id ="cb83-2 "> < a href ="#cb83-2 " tabindex ="-1 "> </ a > ntl_folder < span class ="op "> =</ span > < span class ="st "> 'https://storage.googleapis.com/spatialthoughts-public-data/ntl/npp_viirs_ntl'</ span > </ span >
1432- < span id ="cb83-3 "> < a href ="#cb83-3 " tabindex ="-1 "> </ a > </ span >
1433- < span id ="cb83-4 "> < a href ="#cb83-4 " tabindex ="-1 "> </ a > da_list < span class ="op "> =</ span > []</ span >
1434- < span id ="cb83-5 "> < a href ="#cb83-5 " tabindex ="-1 "> </ a > </ span >
1435- < span id ="cb83-6 "> < a href ="#cb83-6 " tabindex ="-1 "> </ a > < span class ="cf "> for</ span > year < span class ="kw "> in</ span > < span class ="bu "> range</ span > (start_year, end_year < span class ="op "> +</ span > < span class ="dv "> 1</ span > ):</ span >
1436- < span id ="cb83-7 "> < a href ="#cb83-7 " tabindex ="-1 "> </ a > cog_url < span class ="op "> =</ span > < span class ="ss "> f'</ span > < span class ="sc "> {</ span > ntl_folder< span class ="sc "> }</ span > < span class ="ss "> /</ span > < span class ="sc "> {</ span > year< span class ="sc "> }</ span > < span class ="ss "> .tif'</ span > </ span >
1437- < span id ="cb83-8 "> < a href ="#cb83-8 " tabindex ="-1 "> </ a > da < span class ="op "> =</ span > rxr.open_rasterio(</ span >
1438- < span id ="cb83-9 "> < a href ="#cb83-9 " tabindex ="-1 "> </ a > cog_url,</ span >
1439- < span id ="cb83-10 "> < a href ="#cb83-10 " tabindex ="-1 "> </ a > chunks< span class ="op "> =</ span > < span class ="va "> True</ span > ).rio.clip_box(< span class ="op "> *</ span > bbox)</ span >
1440- < span id ="cb83-11 "> < a href ="#cb83-11 " tabindex ="-1 "> </ a > dt < span class ="op "> =</ span > pd.to_datetime(year, < span class ="bu "> format</ span > < span class ="op "> =</ span > < span class ="st "> '%Y'</ span > )</ span >
1441- < span id ="cb83-12 "> < a href ="#cb83-12 " tabindex ="-1 "> </ a > da < span class ="op "> =</ span > da.assign_coords(time < span class ="op "> =</ span > dt)</ span >
1442- < span id ="cb83-13 "> < a href ="#cb83-13 " tabindex ="-1 "> </ a > da < span class ="op "> =</ span > da.expand_dims(dim< span class ="op "> =</ span > < span class ="st "> "time"</ span > )</ span >
1443- < span id ="cb83-14 "> < a href ="#cb83-14 " tabindex ="-1 "> </ a > da_list.append(da)</ span > </ code > </ pre > </ div >
1429+ class ="sourceCode python "> < code class ="sourceCode python "> < span id ="cb83-1 "> < a href ="#cb83-1 " tabindex ="-1 "> </ a > start_year < span class ="op "> =</ span > < span class ="dv "> 2015</ span > </ span >
1430+ < span id ="cb83-2 "> < a href ="#cb83-2 " tabindex ="-1 "> </ a > end_year < span class ="op "> =</ span > < span class ="dv "> 2020</ span > </ span > </ code > </ pre > </ div >
14441431< div class ="sourceCode " id ="cb84 "> < pre
1445- class ="sourceCode python "> < code class ="sourceCode python "> < span id ="cb84-1 "> < a href ="#cb84-1 " tabindex ="-1 "> </ a > ntl_datacube < span class ="op "> =</ span > xr.concat(da_list, dim< span class ="op "> =</ span > < span class ="st "> 'time'</ span > ).chunk(< span class ="st "> 'auto'</ span > )</ span >
1446- < span id ="cb84-2 "> < a href ="#cb84-2 " tabindex ="-1 "> </ a > ntl_datacube</ span > </ code > </ pre > </ div >
1432+ class ="sourceCode python "> < code class ="sourceCode python "> < span id ="cb84-1 "> < a href ="#cb84-1 " tabindex ="-1 "> </ a > </ span >
1433+ < span id ="cb84-2 "> < a href ="#cb84-2 " tabindex ="-1 "> </ a > ntl_folder < span class ="op "> =</ span > < span class ="st "> 'https://storage.googleapis.com/spatialthoughts-public-data/ntl/npp_viirs_ntl'</ span > </ span >
1434+ < span id ="cb84-3 "> < a href ="#cb84-3 " tabindex ="-1 "> </ a > </ span >
1435+ < span id ="cb84-4 "> < a href ="#cb84-4 " tabindex ="-1 "> </ a > da_list < span class ="op "> =</ span > []</ span >
1436+ < span id ="cb84-5 "> < a href ="#cb84-5 " tabindex ="-1 "> </ a > </ span >
1437+ < span id ="cb84-6 "> < a href ="#cb84-6 " tabindex ="-1 "> </ a > < span class ="cf "> for</ span > year < span class ="kw "> in</ span > < span class ="bu "> range</ span > (start_year, end_year < span class ="op "> +</ span > < span class ="dv "> 1</ span > ):</ span >
1438+ < span id ="cb84-7 "> < a href ="#cb84-7 " tabindex ="-1 "> </ a > cog_url < span class ="op "> =</ span > < span class ="ss "> f'</ span > < span class ="sc "> {</ span > ntl_folder< span class ="sc "> }</ span > < span class ="ss "> /</ span > < span class ="sc "> {</ span > year< span class ="sc "> }</ span > < span class ="ss "> .tif'</ span > </ span >
1439+ < span id ="cb84-8 "> < a href ="#cb84-8 " tabindex ="-1 "> </ a > da < span class ="op "> =</ span > rxr.open_rasterio(</ span >
1440+ < span id ="cb84-9 "> < a href ="#cb84-9 " tabindex ="-1 "> </ a > cog_url,</ span >
1441+ < span id ="cb84-10 "> < a href ="#cb84-10 " tabindex ="-1 "> </ a > chunks< span class ="op "> =</ span > < span class ="va "> True</ span > ).rio.clip_box(< span class ="op "> *</ span > bbox)</ span >
1442+ < span id ="cb84-11 "> < a href ="#cb84-11 " tabindex ="-1 "> </ a > dt < span class ="op "> =</ span > pd.to_datetime(year, < span class ="bu "> format</ span > < span class ="op "> =</ span > < span class ="st "> '%Y'</ span > )</ span >
1443+ < span id ="cb84-12 "> < a href ="#cb84-12 " tabindex ="-1 "> </ a > da < span class ="op "> =</ span > da.assign_coords(time < span class ="op "> =</ span > dt)</ span >
1444+ < span id ="cb84-13 "> < a href ="#cb84-13 " tabindex ="-1 "> </ a > da < span class ="op "> =</ span > da.expand_dims(dim< span class ="op "> =</ span > < span class ="st "> "time"</ span > )</ span >
1445+ < span id ="cb84-14 "> < a href ="#cb84-14 " tabindex ="-1 "> </ a > da_list.append(da)</ span > </ code > </ pre > </ div >
14471446< div class ="sourceCode " id ="cb85 "> < pre
1448- class ="sourceCode python "> < code class ="sourceCode python "> < span id ="cb85-1 "> < a href ="#cb85-1 " tabindex ="-1 "> </ a > ntl_datacube < span class ="op "> =</ span > ntl_datacube.sel(band< span class ="op "> =</ span > < span class ="dv "> 1</ span > , drop< span class ="op "> =</ span > < span class ="va "> True</ span > )</ span > </ code > </ pre > </ div >
1447+ class ="sourceCode python "> < code class ="sourceCode python "> < span id ="cb85-1 "> < a href ="#cb85-1 " tabindex ="-1 "> </ a > ntl_datacube < span class ="op "> =</ span > xr.concat(da_list, dim< span class ="op "> =</ span > < span class ="st "> 'time'</ span > ).chunk(< span class ="st "> 'auto'</ span > )</ span >
1448+ < span id ="cb85-2 "> < a href ="#cb85-2 " tabindex ="-1 "> </ a > ntl_datacube</ span > </ code > </ pre > </ div >
1449+ < div class ="sourceCode " id ="cb86 "> < pre
1450+ class ="sourceCode python "> < code class ="sourceCode python "> < span id ="cb86-1 "> < a href ="#cb86-1 " tabindex ="-1 "> </ a > ntl_datacube < span class ="op "> =</ span > ntl_datacube.sel(band< span class ="op "> =</ span > < span class ="dv "> 1</ span > , drop< span class ="op "> =</ span > < span class ="va "> True</ span > )</ span > </ code > </ pre > </ div >
14491451</ div >
14501452< div id ="zonal-stats " class ="section level3 ">
14511453< h3 > Zonal Stats</ h3 >
@@ -1456,44 +1458,44 @@ <h3>Zonal Stats</h3>
14561458geodataframe, it’s unique value as integer and converts the geodataframe
14571459into a < code > xarray dataset</ code > having dimension and coordinates same
14581460as of given input xarray dataset</ p >
1459- < div class ="sourceCode " id ="cb86 "> < pre
1460- class ="sourceCode python "> < code class ="sourceCode python "> < span id ="cb86-1 "> < a href ="#cb86-1 " tabindex ="-1 "> </ a > < span class ="co "> # Create mask of multiple regions from shapefile</ span > </ span >
1461- < span id ="cb86-2 "> < a href ="#cb86-2 " tabindex ="-1 "> </ a > mask < span class ="op "> =</ span > regionmask.mask_3D_geopandas(</ span >
1462- < span id ="cb86-3 "> < a href ="#cb86-3 " tabindex ="-1 "> </ a > zones,</ span >
1463- < span id ="cb86-4 "> < a href ="#cb86-4 " tabindex ="-1 "> </ a > ntl_datacube.x,</ span >
1464- < span id ="cb86-5 "> < a href ="#cb86-5 " tabindex ="-1 "> </ a > ntl_datacube.y,</ span >
1465- < span id ="cb86-6 "> < a href ="#cb86-6 " tabindex ="-1 "> </ a > drop< span class ="op "> =</ span > < span class ="va "> True</ span > ,</ span >
1466- < span id ="cb86-7 "> < a href ="#cb86-7 " tabindex ="-1 "> </ a > numbers< span class ="op "> =</ span > < span class ="st "> "id"</ span > ,</ span >
1467- < span id ="cb86-8 "> < a href ="#cb86-8 " tabindex ="-1 "> </ a > overlap< span class ="op "> =</ span > < span class ="va "> True</ span > </ span >
1468- < span id ="cb86-9 "> < a href ="#cb86-9 " tabindex ="-1 "> </ a > )</ span > </ code > </ pre > </ div >
14691461< div class ="sourceCode " id ="cb87 "> < pre
1470- class ="sourceCode python "> < code class ="sourceCode python "> < span id ="cb87-1 "> < a href ="#cb87-1 " tabindex ="-1 "> </ a > ntl_datacube < span class ="op "> =</ span > ntl_datacube.where(mask).chunk(< span class ="st "> 'auto'</ span > )</ span >
1471- < span id ="cb87-2 "> < a href ="#cb87-2 " tabindex ="-1 "> </ a > ntl_datacube</ span > </ code > </ pre > </ div >
1462+ class ="sourceCode python "> < code class ="sourceCode python "> < span id ="cb87-1 "> < a href ="#cb87-1 " tabindex ="-1 "> </ a > < span class ="co "> # Create mask of multiple regions from shapefile</ span > </ span >
1463+ < span id ="cb87-2 "> < a href ="#cb87-2 " tabindex ="-1 "> </ a > mask < span class ="op "> =</ span > regionmask.mask_3D_geopandas(</ span >
1464+ < span id ="cb87-3 "> < a href ="#cb87-3 " tabindex ="-1 "> </ a > zones,</ span >
1465+ < span id ="cb87-4 "> < a href ="#cb87-4 " tabindex ="-1 "> </ a > ntl_datacube.x,</ span >
1466+ < span id ="cb87-5 "> < a href ="#cb87-5 " tabindex ="-1 "> </ a > ntl_datacube.y,</ span >
1467+ < span id ="cb87-6 "> < a href ="#cb87-6 " tabindex ="-1 "> </ a > drop< span class ="op "> =</ span > < span class ="va "> True</ span > ,</ span >
1468+ < span id ="cb87-7 "> < a href ="#cb87-7 " tabindex ="-1 "> </ a > numbers< span class ="op "> =</ span > < span class ="st "> "id"</ span > ,</ span >
1469+ < span id ="cb87-8 "> < a href ="#cb87-8 " tabindex ="-1 "> </ a > overlap< span class ="op "> =</ span > < span class ="va "> True</ span > </ span >
1470+ < span id ="cb87-9 "> < a href ="#cb87-9 " tabindex ="-1 "> </ a > )</ span > </ code > </ pre > </ div >
14721471< div class ="sourceCode " id ="cb88 "> < pre
1473- class ="sourceCode python "> < code class ="sourceCode python "> < span id ="cb88-1 "> < a href ="#cb88-1 " tabindex ="-1 "> </ a > grouped < span class ="op "> =</ span > ntl_datacube.groupby( < span class =" st " > 'time' </ span > ). < span class =" bu " > sum </ span > ([ < span class ="st "> 'x '</ span > , < span class =" st " > 'y' </ span > ] )</ span >
1474- < span id ="cb88-2 "> < a href ="#cb88-2 " tabindex ="-1 "> </ a > grouped </ span > </ code > </ pre > </ div >
1472+ class ="sourceCode python "> < code class ="sourceCode python "> < span id ="cb88-1 "> < a href ="#cb88-1 " tabindex ="-1 "> </ a > ntl_datacube < span class ="op "> =</ span > ntl_datacube.where(mask).chunk( < span class ="st "> 'auto '</ span > )</ span >
1473+ < span id ="cb88-2 "> < a href ="#cb88-2 " tabindex ="-1 "> </ a > ntl_datacube </ span > </ code > </ pre > </ div >
14751474< div class ="sourceCode " id ="cb89 "> < pre
1476- class ="sourceCode python "> < code class ="sourceCode python "> < span id ="cb89-1 "> < a href ="#cb89-1 " tabindex ="-1 "> </ a > < span class ="op "> %</ span > time grouped < span class ="op "> =</ span > grouped.compute()</ span > </ code > </ pre > </ div >
1475+ class ="sourceCode python "> < code class ="sourceCode python "> < span id ="cb89-1 "> < a href ="#cb89-1 " tabindex ="-1 "> </ a > grouped < span class ="op "> =</ span > ntl_datacube.groupby(< span class ="st "> 'time'</ span > ).< span class ="bu "> sum</ span > ([< span class ="st "> 'x'</ span > ,< span class ="st "> 'y'</ span > ])</ span >
1476+ < span id ="cb89-2 "> < a href ="#cb89-2 " tabindex ="-1 "> </ a > grouped</ span > </ code > </ pre > </ div >
14771477< div class ="sourceCode " id ="cb90 "> < pre
1478- class ="sourceCode python "> < code class ="sourceCode python "> < span id ="cb90-1 "> < a href ="#cb90-1 " tabindex ="-1 "> </ a > stats < span class ="op "> =</ span > grouped.drop(< span class ="st "> 'spatial_ref'</ span > ).to_dataframe(< span class ="st "> 'sum'</ span > ).reset_index()</ span >
1479- < span id ="cb90-2 "> < a href ="#cb90-2 " tabindex ="-1 "> </ a > stats</ span > </ code > </ pre > </ div >
1478+ class ="sourceCode python "> < code class ="sourceCode python "> < span id ="cb90-1 "> < a href ="#cb90-1 " tabindex ="-1 "> </ a > < span class ="op "> %</ span > time grouped < span class ="op "> =</ span > grouped.compute()</ span > </ code > </ pre > </ div >
14801479< div class ="sourceCode " id ="cb91 "> < pre
1481- class ="sourceCode python "> < code class ="sourceCode python "> < span id ="cb91-1 "> < a href ="#cb91-1 " tabindex ="-1 "> </ a > stats[< span class ="st "> 'year'</ span > ] < span class ="op "> =</ span > stats.time.dt.year</ span >
1482- < span id ="cb91-2 "> < a href ="#cb91-2 " tabindex ="-1 "> </ a > stats < span class ="op "> =</ span > stats.rename(columns< span class ="op "> =</ span > {< span class ="st "> 'region'</ span > : < span class ="st "> 'id'</ span > })</ span >
1483- < span id ="cb91-3 "> < a href ="#cb91-3 " tabindex ="-1 "> </ a > stats</ span > </ code > </ pre > </ div >
1480+ class ="sourceCode python "> < code class ="sourceCode python "> < span id ="cb91-1 "> < a href ="#cb91-1 " tabindex ="-1 "> </ a > stats < span class ="op "> =</ span > grouped.drop(< span class ="st "> 'spatial_ref'</ span > ).to_dataframe(< span class ="st "> 'sum'</ span > ).reset_index()</ span >
1481+ < span id ="cb91-2 "> < a href ="#cb91-2 " tabindex ="-1 "> </ a > stats</ span > </ code > </ pre > </ div >
14841482< div class ="sourceCode " id ="cb92 "> < pre
1485- class ="sourceCode python "> < code class ="sourceCode python "> < span id ="cb92-1 "> < a href ="#cb92-1 " tabindex ="-1 "> </ a > results < span class ="op " > = </ span > zones.merge(stats, on < span class ="op "> =</ span > < span class =" st " > 'id' </ span > ) </ span >
1486- < span id ="cb92-2 "> < a href ="#cb92-2 " tabindex ="-1 "> </ a > results < /span > </ code > </ pre > </ div >
1487- < p > Finally, we save the result to disk. </ p >
1483+ class ="sourceCode python "> < code class ="sourceCode python "> < span id ="cb92-1 "> < a href ="#cb92-1 " tabindex ="-1 "> </ a > stats[ < span class ="st " > 'year' </ span > ] < span class ="op "> =</ span > stats.time.dt.year </ span >
1484+ < span id ="cb92-2 "> < a href ="#cb92-2 " tabindex ="-1 "> </ a > stats < span class =" op " > = < /span > stats.rename(columns < span class =" op " > = </ span > { < span class =" st " > 'region' </ span > : < span class =" st " > 'id' </ span > }) </ span >
1485+ < span id =" cb92-3 " > < a href =" #cb92-3 " tabindex =" -1 " > </ a > stats </ span > </ code > </ pre > </ div >
14881486< div class ="sourceCode " id ="cb93 "> < pre
1489- class ="sourceCode python "> < code class ="sourceCode python "> < span id ="cb93-1 "> < a href ="#cb93-1 " tabindex ="-1 "> </ a > output < span class ="op "> =</ span > results[[< span class ="st "> 'adm1_code'</ span > , < span class ="st "> 'name'</ span > , < span class ="st "> 'iso_a2'</ span > , < span class ="st "> 'year'</ span > , < span class ="st "> 'sum'</ span > ]]</ span >
1490- < span id ="cb93-2 "> < a href ="#cb93-2 " tabindex ="-1 "> </ a > output</ span > </ code > </ pre > </ div >
1487+ class ="sourceCode python "> < code class ="sourceCode python "> < span id ="cb93-1 "> < a href ="#cb93-1 " tabindex ="-1 "> </ a > results < span class ="op "> =</ span > zones.merge(stats, on< span class ="op "> =</ span > < span class ="st "> 'id'</ span > )</ span >
1488+ < span id ="cb93-2 "> < a href ="#cb93-2 " tabindex ="-1 "> </ a > results</ span > </ code > </ pre > </ div >
1489+ < p > Finally, we save the result to disk.</ p >
14911490< div class ="sourceCode " id ="cb94 "> < pre
1492- class ="sourceCode python "> < code class ="sourceCode python "> < span id ="cb94-1 "> < a href ="#cb94-1 " tabindex ="-1 "> </ a > output_file < span class ="op "> =</ span > < span class ="st "> 'output.csv'</ span > </ span >
1493- < span id ="cb94-2 "> < a href ="#cb94-2 " tabindex ="-1 "> </ a > output_path < span class ="op "> =</ span > os.path.join(output_folder, output_file)</ span >
1494- < span id ="cb94-3 "> < a href ="#cb94-3 " tabindex ="-1 "> </ a > </ span >
1495- < span id ="cb94-4 "> < a href ="#cb94-4 " tabindex ="-1 "> </ a > output.to_csv(output_path, index< span class ="op "> =</ span > < span class ="va "> False</ span > )</ span >
1496- < span id ="cb94-5 "> < a href ="#cb94-5 " tabindex ="-1 "> </ a > < span class ="bu "> print</ span > (< span class ="st "> 'Successfully written output file at </ span > < span class ="sc "> {}</ span > < span class ="st "> '</ span > .< span class ="bu "> format</ span > (output_path))</ span > </ code > </ pre > </ div >
1491+ class ="sourceCode python "> < code class ="sourceCode python "> < span id ="cb94-1 "> < a href ="#cb94-1 " tabindex ="-1 "> </ a > output < span class ="op "> =</ span > results[[< span class ="st "> 'adm1_code'</ span > , < span class ="st "> 'name'</ span > , < span class ="st "> 'iso_a2'</ span > , < span class ="st "> 'year'</ span > , < span class ="st "> 'sum'</ span > ]]</ span >
1492+ < span id ="cb94-2 "> < a href ="#cb94-2 " tabindex ="-1 "> </ a > output</ span > </ code > </ pre > </ div >
1493+ < div class ="sourceCode " id ="cb95 "> < pre
1494+ class ="sourceCode python "> < code class ="sourceCode python "> < span id ="cb95-1 "> < a href ="#cb95-1 " tabindex ="-1 "> </ a > output_file < span class ="op "> =</ span > < span class ="st "> 'output.csv'</ span > </ span >
1495+ < span id ="cb95-2 "> < a href ="#cb95-2 " tabindex ="-1 "> </ a > output_path < span class ="op "> =</ span > os.path.join(output_folder, output_file)</ span >
1496+ < span id ="cb95-3 "> < a href ="#cb95-3 " tabindex ="-1 "> </ a > </ span >
1497+ < span id ="cb95-4 "> < a href ="#cb95-4 " tabindex ="-1 "> </ a > output.to_csv(output_path, index< span class ="op "> =</ span > < span class ="va "> False</ span > )</ span >
1498+ < span id ="cb95-5 "> < a href ="#cb95-5 " tabindex ="-1 "> </ a > < span class ="bu "> print</ span > (< span class ="st "> 'Successfully written output file at </ span > < span class ="sc "> {}</ span > < span class ="st "> '</ span > .< span class ="bu "> format</ span > (output_path))</ span > </ code > </ pre > </ div >
14971499</ div >
14981500< div id ="exercise-3 " class ="section level3 ">
14991501< h3 > Exercise</ h3 >
0 commit comments