A regional view of the past with applications to the present. Python spatial and tabular analysis using multiple variables to look at a single time slice. This code can be a combination of cloud and local based analysis. A purpose for this exercise is to use use public research data for derived conclusions view the impacts climate change has on agriculture a look at what farmers can do and why we need to support them.
+
NASA Acres has defined the Essential Agricultural Variables or EAVs. These were designed by NASA Acres to define the capabilities of the top satellite data scientists and practitioners who make up NASA Acres Research, Development, and Extension partners, and the needs of decision-making-collaborators already in our network, to identify an initial set of focus.
+
+
Cropland and Crop Type Mapping
+
Crop and Crop Type Area Estimation
+
+
Derived results will be the county level statistics of impacts of harvestable acreage. To achieve this, a basic workflow we can
+
+
Use api to call in the CDL dataset to map crop types.
+
+
Then
+
+
Access the Sea level rise (elevation dataset) to identify new areas of potential areas that are at risk of future flooding. Clip to county of interest
+
Compare the CDL with the SLR mask and with out to identify the NOAA estimated loss of land.
+
+
Then using MODIS create a time series for NDVI of the masked crop land, derive insights about trends in NDVI
+
+
Use NDVI to create a time series looking backwards in time at areas that have experienced flooding to visualize the moving from productive farms to moderate quality.
+
+
This all together would allow us to make a predictive analysis for Maryland in the future under the projections of sea level rise. Given the current conditions, subtracting the sea level rise inundated areas. Then given the remaining areas, making predictions about the trends of NDVI given the trends in historical areas that are near impacted areas. We cannot say that “sea level rise is driving the decrease” because there are many other factor and decisions that farmers make around how much to plant, how much to harvest, chemicals applied, growing degree days, soil and soil moisture conditions. But this module should provide the ability to make insights about potential causes and impacts.
+
The data story that we have derived is about sea level rise in Maryland and the impact that it has on the production levels with in the state. From this module, we can provide students with the ability to draw from multiple data sources, as well as derive insights using historical and future viewing data sets.
+
The U.S. Mid-Atlantic has seen higher rates of sea level rise, marshes may be especially vulnerable. “In the Chesapeake Bay, sea level rise has already contributed to the degradation of over 80,000 ha (70%) of tidal marsh” (Taylor et al. 2020). This view can help us to understand the impacts that small levels of sea level rise have on land.
+
We can prompt the user to think about future impacts outside the direct sea level rise projections, allowing them to include a more full picture, and finally using that picture to identify economic impacts that action or inaction causes. This begs the question, what can the public do to enact changes, rather than putting the pressure on farmers to change?
To understand the nature of Maryland agriculture we can begin by using the CDL. This raster dataset comes at a 30 meter spatial resolution. We can access this data through their api. The current configuration allows data to be pulled at the county level.
for condition in ['corn_below_water', 'corn_above_water']:
+
+print(f"\n{condition.replace('_', ' ').title()} Corn Analysis per MODIS:")
+
+for stat in ['mean', 'max']:
+
+print(f"{stat.capitalize()} NDVI: {ndvi_analysis[condition][stat]}")
What are the insights that we can derive from this module?
+
What are the impacts of land being close to the water currently?
+
What is the potential loss of productive land?
+
Where could there be errors given the sampling methods?
+
+
+
+
What to do next
+
+
Clip the MODIS data to the CDL data of only the agricultural values of interest
+
Find the average NDVI of these values
+
Clip this to the SLR area
+
Find the average NDVI of these values
+
Compare the changes for the year of 2018
+
+
This will tell us the impacted area that could be effected by the floods
+
As we look into the future we should see this change become more drastic
+
+
+
ABOUT NASA ACRES AND THE CLIMATE RESILIENCE NETWORK
+
+
NASA Acres
+
NASA Acres is NASA’s U.S.-focused agriculture and food security Consortium. NASA Acres is commissioned under NASA’s Applied Sciences Program and led by the University of Maryland, with over 30 partner institutions and growing. The consortium approach brings together public and private sectors, and allows for flexible partnerships and rapid action in delivering the benefits of NASA data and tools to U.S. agriculture decision makers. We work with people across the agriculture sector to develop Earth observation (EO)-based data and tools that strengthen U.S. agriculture.
+
+
+
Climate Resilience Network
+
Climate Resilience Network has the goal of unifying the University of Maryland’s groundbreaking research with the immediate requirements of Maryland’s decision-makers to address the complexities of climate change as it affects our home. GEOG’s Grand Challenges project is aimed at promoting climate-resilient agriculture in Maryland, creating interactive tools designed to support research based decision-making.italicized text
+Taylor, Lotem, David Curson, Gregory M. Verutes, and Chad Wilsey. 2020. “Mapping Sea Level Rise Impacts to Identify Climate Change Adaptation Opportunities in the Chesapeake and Delaware Bays, USA.”Wetlands Ecology and Management 28 (3): 527–41. https://doi.org/10.1007/s11273-020-09729-w.
+
+
+
+
+
+
+
+
\ No newline at end of file
diff --git a/docs/index.html b/docs/index.html
new file mode 100644
index 0000000..bade81e
--- /dev/null
+++ b/docs/index.html
@@ -0,0 +1,4129 @@
+
+
+
+
+
+
+
+
+
+Welcome SCHOOL Module 2: Climate and Agriculture
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
+
Welcome to the fourth module of the SCHOOL curriculum!
+
The Science Core Heuristics for Open Science Outcomes in Learning (SCHOOL) is part of the NASA Transform to Open Science (TOPS) Training (TOPST) initiative, designed to teach the data science lifecycle using data from the NASA Earth Sciences division and to foster an inclusive culture of open science. You can learn more about the SCHOOL Project and other modules on the SCHOOL Project home page.
+
This TOPS SCHOOL module on the theme of ’Climate and Agriculture” provides a foundational understanding of the complex interplay between Earth’s ecosystems, climate and human-induced climate change, and environmental justice. We will explore how we can use spatial data for the study of climate, agriculture, and their impacts on vulnerable communities. You will also gain critical insights into how a changing climate affect agricultural resources people rely on for food. As you engage with the lessons, you will develop a comprehensive understanding of key concepts of earth’s climate system, the significance of social vulnerability to climate risks, and the tools available for analyzing environmental and agricultural data. Ultimately, these lessons will empower you to contribute meaningfully to discussions and actions pertaining to earth’s climate, agricultural systems, and environmental justice.
+
+
+
+
+
+
+Tip
+
+
+
+
This module is designed to provide undergraduate students and early-career researchers with an introduction to the data science “life cycle” and sample code describing how Open Science principles can be effectively applied to earth sciences, particularly in the context of air pollution and its impact on health.
+
Teaching earth science topics and coding skills in depth is beyond the scope of the SCHOOL Modules. Instead, the SCHOOL Project aims to provide users with the skills to adapt the skills learned in SCHOOL lessons to the users’ own Open Science workflow. To learn more about Open Science, explore NASA’s TOPS Open Science 101 Curriculum. To explore other themes in the SCHOOL project, visit our Modules Page.
+
+
+
+
Module 4: Climate and Agriculture datasets and use cases cover:
+
+
NASA ACRES Climate Resilience Network: This lesson will use public research data for derived conclusions and to view the impacts climate change has on agriculture.
+
This course was made possible thanks to the work of our NASA Transform to Open Science (TOPS) team, our SCHOOL Open Science team, open science Subject Matter Experts (SMEs), and the SCHOOL Development team!
+
+
+
+
+
+
+
+
+
+
+
\ No newline at end of file
diff --git a/docs/search.json b/docs/search.json
new file mode 100644
index 0000000..af9c176
--- /dev/null
+++ b/docs/search.json
@@ -0,0 +1,44 @@
+[
+ {
+ "objectID": "NASA-Acres-Climate-Resilience-Network.html",
+ "href": "NASA-Acres-Climate-Resilience-Network.html",
+ "title": "NASA Acres and Climate Resilience Network",
+ "section": "",
+ "text": "A regional view of the past with applications to the present. Python spatial and tabular analysis using multiple variables to look at a single time slice. This code can be a combination of cloud and local based analysis. A purpose for this exercise is to use use public research data for derived conclusions view the impacts climate change has on agriculture a look at what farmers can do and why we need to support them.\nNASA Acres has defined the Essential Agricultural Variables or EAVs. These were designed by NASA Acres to define the capabilities of the top satellite data scientists and practitioners who make up NASA Acres Research, Development, and Extension partners, and the needs of decision-making-collaborators already in our network, to identify an initial set of focus.\n\nCropland and Crop Type Mapping\nCrop and Crop Type Area Estimation\n\nDerived results will be the county level statistics of impacts of harvestable acreage. To achieve this, a basic workflow we can\n\nUse api to call in the CDL dataset to map crop types.\n\nThen\n\nAccess the Sea level rise (elevation dataset) to identify new areas of potential areas that are at risk of future flooding. Clip to county of interest\nCompare the CDL with the SLR mask and with out to identify the NOAA estimated loss of land.\n\nThen using MODIS create a time series for NDVI of the masked crop land, derive insights about trends in NDVI\n\nUse NDVI to create a time series looking backwards in time at areas that have experienced flooding to visualize the moving from productive farms to moderate quality.\n\nThis all together would allow us to make a predictive analysis for Maryland in the future under the projections of sea level rise. Given the current conditions, subtracting the sea level rise inundated areas. Then given the remaining areas, making predictions about the trends of NDVI given the trends in historical areas that are near impacted areas. We cannot say that “sea level rise is driving the decrease” because there are many other factor and decisions that farmers make around how much to plant, how much to harvest, chemicals applied, growing degree days, soil and soil moisture conditions. But this module should provide the ability to make insights about potential causes and impacts.\nThe data story that we have derived is about sea level rise in Maryland and the impact that it has on the production levels with in the state. From this module, we can provide students with the ability to draw from multiple data sources, as well as derive insights using historical and future viewing data sets.\nThe U.S. Mid-Atlantic has seen higher rates of sea level rise, marshes may be especially vulnerable. “In the Chesapeake Bay, sea level rise has already contributed to the degradation of over 80,000 ha (70%) of tidal marsh” (Taylor et al. 2020). This view can help us to understand the impacts that small levels of sea level rise have on land.\nWe can prompt the user to think about future impacts outside the direct sea level rise projections, allowing them to include a more full picture, and finally using that picture to identify economic impacts that action or inaction causes. This begs the question, what can the public do to enact changes, rather than putting the pressure on farmers to change?\n\n!pip install arcgis arcgis-mapping rasterio earthaccess\n\nRequirement already satisfied: arcgis in c:\\users\\jmartine\\appdata\\local\\esri\\conda\\envs\\agsdev\\lib\\site-packages (2.2.0)\n\n\n\n\n\n\n#import the required packages\nimport geopandas as gpd\nimport rasterio\nimport numpy as np\nimport pandas as pd\nimport matplotlib.pyplot as plt\nfrom shapely.geometry import box\nimport numpy as np\nimport requests\nimport xml.etree.ElementTree as ET\nimport matplotlib.pyplot as plt\nfrom io import BytesIO\nimport os\nfrom rasterio.warp import calculate_default_transform, reproject, Resampling\nfrom datetime import datetime, timedelta\n\n\n#define fuctions\n\ndef reproject_raster(src_dataset, src_crs, src_transform, dst_crs='EPSG:4326'):\n dst_crs = rasterio.crs.CRS.from_string(dst_crs)\n dst_transform, dst_width, dst_height = calculate_default_transform(\n src_crs, dst_crs,\n src_dataset.width, src_dataset.height,\n *src_dataset.bounds)\n dst_kwargs = src_dataset.meta.copy()\n dst_kwargs.update({\n 'crs': dst_crs,\n 'transform': dst_transform,\n 'width': dst_width,\n 'height': dst_height})\n dst_data = np.zeros((src_dataset.count, dst_height, dst_width), dtype=src_dataset.dtypes[0])\n for i in range(1, src_dataset.count + 1):\n reproject(\n source=rasterio.band(src_dataset, i),\n destination=dst_data[i-1],\n src_transform=src_transform,\n src_crs=src_crs,\n dst_transform=dst_transform,\n dst_crs=dst_crs,\n resampling=Resampling.nearest)\n return dst_data, dst_kwargs\n\ndef format_cdl_url(fips,year):\n base_url = \"https://nassgeodata.gmu.edu/axis2/services/CDLService/GetCDLFile\"\n url = f\"{base_url}?year={year}&fips={fips}\"\n return url\n\ndef get_srl_raster():\n slr_path = 'data/MD_East_slr_depth_3_5ft.tif'\n with rasterio.open(slr_path) as slr_raster:\n slr_meta = slr_raster.meta\n slr_meta['nodata'] = 0\n slr_reprojected, slr_meta_reprojected = reproject_raster(src_dataset=slr_raster, src_crs=slr_meta['crs'], src_transform=slr_meta['transform'])\n return slr_reprojected, slr_meta_reprojected"
+ },
+ {
+ "objectID": "NASA-Acres-Climate-Resilience-Network.html#combine-and-compare-the-results",
+ "href": "NASA-Acres-Climate-Resilience-Network.html#combine-and-compare-the-results",
+ "title": "NASA Acres and Climate Resilience Network",
+ "section": "Combine and compare the results",
+ "text": "Combine and compare the results\n\nndvi_resampled, ndvi_resampled_meta = resample_using_cdl(ndvi, meta, cdl_meta_reprojected)\n\nslr_resampled, slr_resampled_meta = resample_using_cdl(slr_reprojected, slr_meta_reprojected, cdl_meta_reprojected)\n\nprint(\"CDL shape:\", cdl_reprojected.shape)\n\nprint(\"Original NDVI shape:\", ndvi.shape)\n\nprint(\"Resampled NDVI shape:\", ndvi_resampled.shape)\n\nprint(\"Original SLR shape:\", slr_reprojected.shape)\n\nprint(\"Resampled SLR shape:\", slr_resampled.shape)\n\nCDL shape: (1, 1681, 2153)\nOriginal NDVI shape: (9, 166, 402)\nResampled NDVI shape: (9, 1681, 2153)\nOriginal SLR shape: (1, 44790, 28460)\nResampled SLR shape: (1, 1681, 2153)\n\n\n\ndef mask_rasters(cdl_raster, slr_raster, ndvi_rasters):\n\n corn_mask = cdl_raster == 1\n\n ABOVE_water_mask = slr_raster <= 1\n\n corn_above_water_mask = corn_mask & ~ABOVE_water_mask\n\n corn_below_water_mask = corn_mask & ABOVE_water_mask\n\n corn_above_water_mask = corn_above_water_mask[0,:,:]\n\n corn_below_water_mask = corn_below_water_mask[0,:,:]\n\n masked_ndvi_below_water = ndvi_rasters[:, corn_below_water_mask]\n\n masked_ndvi_above_water = ndvi_rasters[:, corn_above_water_mask]\n\n return {\n\n 'corn_below_water': {\n\n 'mask': corn_below_water_mask,\n\n 'ndvi': masked_ndvi_below_water,\n\n 'mean': np.mean(masked_ndvi_below_water, axis=1) if masked_ndvi_below_water.size > 0 else None,\n\n 'min': np.min(masked_ndvi_below_water, axis=1) if masked_ndvi_below_water.size > 0 else None,\n\n 'max': np.max(masked_ndvi_below_water, axis=1) if masked_ndvi_below_water.size > 0 else None\n\n },\n\n 'corn_above_water': {\n\n 'mask': corn_above_water_mask,\n\n 'ndvi': masked_ndvi_above_water,\n\n 'mean': np.mean(masked_ndvi_above_water, axis=1) if masked_ndvi_above_water.size > 0 else None,\n\n 'min': np.min(masked_ndvi_above_water, axis=1) if masked_ndvi_above_water.size > 0 else None,\n\n 'max': np.max(masked_ndvi_above_water, axis=1) if masked_ndvi_above_water.size > 0 else None\n\n }\n\n }\n\n\ndef visualize_ndvi_analysis(ndvi_analysis, days):\n\n plt.figure(figsize=(12, 6))\n\n plt.subplot(1, 2, 1)\n\n plt.plot(days, ndvi_analysis['corn_above_water']['mean'], label='Corn Above Water', color='green')\n\n plt.plot(days, ndvi_analysis['corn_below_water']['mean'], label='Corn Below Water', color='blue')\n\n plt.title('Mean NDVI')\n\n plt.xlabel('Days')\n\n plt.xticks(rotation=90)\n\n plt.ylabel('NDVI')\n\n plt.legend()\n\n plt.subplot(1, 2, 2)\n\n plt.plot(days, ndvi_analysis['corn_above_water']['max'], label='Corn Above Water', color='green')\n\n plt.plot(days, ndvi_analysis['corn_below_water']['max'], label='Corn Below Water', color='blue')\n\n plt.title('Maximum NDVI')\n\n plt.xlabel('Days')\n\n plt.xticks(rotation=90)\n\n plt.ylabel('NDVI')\n\n plt.legend()\n\n plt.tight_layout()\n\n plt.show()\n\n\nndvi_analysis = mask_rasters(cdl_reprojected, slr_resampled, ndvi_resampled)\n\nvisualize_ndvi_analysis(ndvi_analysis, day)\n\n\n\n\n\nfor condition in ['corn_below_water', 'corn_above_water']:\n\n print(f\"\\n{condition.replace('_', ' ').title()} Corn Analysis per MODIS:\")\n\n for stat in ['mean', 'max']:\n\n print(f\"{stat.capitalize()} NDVI: {ndvi_analysis[condition][stat]}\")\n\n\nCorn Below Water Corn Analysis per MODIS:\nMean NDVI: [21.55594865 20.51568452 7.4782058 26.24400004 26.60275338 23.33767715\n 13.16068297 11.80940537 12.10314245]\nMax NDVI: [81.11699055 76.91399092 80.58758897 82.75593908 89.70690955 91.41303851\n 88.49125289 85.55169137 90.92957006]\n\nCorn Above Water Corn Analysis per MODIS:\nMean NDVI: [54.61076598 59.69546623 65.22443398 65.90975493 78.52357012 76.01061572\n 73.06511618 54.97502392 54.67252283]\nMax NDVI: [95.95130762 89.55483178 90.90481424 95.45510286 96.85160514 99.7260847\n 93.69302615 98.97050876 96.44831104]\n\n\n\ndef visualize_ndvi_histogram(ndvi_analysis):\n\n plt.figure(figsize=(12, 6))\n\n plt.hist(ndvi_analysis['corn_below_water']['ndvi'].flatten(), bins=30, alpha=0.5, label='Corn Below Water', color='blue')\n\n plt.hist(ndvi_analysis['corn_above_water']['ndvi'].flatten(), bins=30, alpha=0.5, label='Corn Above Water', color='green')\n\n plt.xlabel(\"NDVI values\")\n\n plt.ylabel(\"Frequency\")\n\n plt.title(\"Histogram of NDVI values for corn\")\n\n plt.legend(loc='upper right')\n\nvisualize_ndvi_histogram(ndvi_analysis)\n\n\n\n\nWhat are the insights that we can derive from this module?\nWhat are the impacts of land being close to the water currently?\nWhat is the potential loss of productive land?\nWhere could there be errors given the sampling methods?"
+ },
+ {
+ "objectID": "NASA-Acres-Climate-Resilience-Network.html#nasa-acres",
+ "href": "NASA-Acres-Climate-Resilience-Network.html#nasa-acres",
+ "title": "NASA Acres and Climate Resilience Network",
+ "section": "NASA Acres",
+ "text": "NASA Acres\nNASA Acres is NASA’s U.S.-focused agriculture and food security Consortium. NASA Acres is commissioned under NASA’s Applied Sciences Program and led by the University of Maryland, with over 30 partner institutions and growing. The consortium approach brings together public and private sectors, and allows for flexible partnerships and rapid action in delivering the benefits of NASA data and tools to U.S. agriculture decision makers. We work with people across the agriculture sector to develop Earth observation (EO)-based data and tools that strengthen U.S. agriculture."
+ },
+ {
+ "objectID": "NASA-Acres-Climate-Resilience-Network.html#climate-resilience-network",
+ "href": "NASA-Acres-Climate-Resilience-Network.html#climate-resilience-network",
+ "title": "NASA Acres and Climate Resilience Network",
+ "section": "Climate Resilience Network",
+ "text": "Climate Resilience Network\nClimate Resilience Network has the goal of unifying the University of Maryland’s groundbreaking research with the immediate requirements of Maryland’s decision-makers to address the complexities of climate change as it affects our home. GEOG’s Grand Challenges project is aimed at promoting climate-resilient agriculture in Maryland, creating interactive tools designed to support research based decision-making.italicized text\nCOUNTY LEVEL SHAPEFILE https://resilience.climate.gov/datasets/nationalclimate::u-s-climate-and-coastal-inundation-projections-by-geography/explore?layer=0&location=2.619718%2C0.314282%2C1.74"
+ },
+ {
+ "objectID": "index.html",
+ "href": "index.html",
+ "title": "Welcome SCHOOL Module 2: Climate and Agriculture",
+ "section": "",
+ "text": "Welcome to the fourth module of the SCHOOL curriculum!\nThe Science Core Heuristics for Open Science Outcomes in Learning (SCHOOL) is part of the NASA Transform to Open Science (TOPS) Training (TOPST) initiative, designed to teach the data science lifecycle using data from the NASA Earth Sciences division and to foster an inclusive culture of open science. You can learn more about the SCHOOL Project and other modules on the SCHOOL Project home page.\nThis TOPS SCHOOL module on the theme of ’Climate and Agriculture” provides a foundational understanding of the complex interplay between Earth’s ecosystems, climate and human-induced climate change, and environmental justice. We will explore how we can use spatial data for the study of climate, agriculture, and their impacts on vulnerable communities. You will also gain critical insights into how a changing climate affect agricultural resources people rely on for food. As you engage with the lessons, you will develop a comprehensive understanding of key concepts of earth’s climate system, the significance of social vulnerability to climate risks, and the tools available for analyzing environmental and agricultural data. Ultimately, these lessons will empower you to contribute meaningfully to discussions and actions pertaining to earth’s climate, agricultural systems, and environmental justice."
+ },
+ {
+ "objectID": "index.html#module-4-climate-and-agriculture-datasets-and-use-cases-cover",
+ "href": "index.html#module-4-climate-and-agriculture-datasets-and-use-cases-cover",
+ "title": "Welcome SCHOOL Module 2: Climate and Agriculture",
+ "section": "Module 4: Climate and Agriculture datasets and use cases cover:",
+ "text": "Module 4: Climate and Agriculture datasets and use cases cover:\n\nNASA ACRES Climate Resilience Network: This lesson will use public research data for derived conclusions and to view the impacts climate change has on agriculture.\n\nLesson 1: Acquiring, Pre-Processing, and Visualizing Student-Monitored Data for New York City schools\n\n\n\nStart Lesson 1\n\nThis course was made possible thanks to the work of our NASA Transform to Open Science (TOPS) team, our SCHOOL Open Science team, open science Subject Matter Experts (SMEs), and the SCHOOL Development team!"
+ }
+]
\ No newline at end of file
diff --git a/index.qmd b/index.qmd
new file mode 100644
index 0000000..b96dc9c
--- /dev/null
+++ b/index.qmd
@@ -0,0 +1,32 @@
+---
+title: "Welcome SCHOOL Module 2: Climate and Agriculture"
+format: html
+---
+
+
+
+Welcome to the fourth module of the SCHOOL curriculum!
+
+The Science Core Heuristics for Open Science Outcomes in Learning (SCHOOL) is part of the [NASA Transform to Open Science (TOPS)](https://nasa.github.io/Transform-to-Open-Science/) Training (TOPST) initiative, designed to teach the data science lifecycle using data from the NASA Earth Sciences division and to foster an inclusive culture of open science. You can learn more about the SCHOOL Project and other modules on the [SCHOOL Project home page.](https://ciesin-geospatial.github.io/TOPSTSCHOOL/)
+
+This TOPS SCHOOL module on the theme of 'Climate and Agriculture" provides a foundational understanding of the complex interplay between Earth's ecosystems, climate and human-induced climate change, and environmental justice. We will explore how we can use spatial data for the study of climate, agriculture, and their impacts on vulnerable communities. You will also gain critical insights into how a changing climate affect agricultural resources people rely on for food. As you engage with the lessons, you will develop a comprehensive understanding of key concepts of earth's climate system, the significance of social vulnerability to climate risks, and the tools available for analyzing environmental and agricultural data. Ultimately, these lessons will empower you to contribute meaningfully to discussions and actions pertaining to earth's climate, agricultural systems, and environmental justice.
+
+
+
+::: {.callout-tip style="color: #5a7a2b;"}
+**This module is designed to provide undergraduate students and early-career researchers with an introduction to the data science “life cycle” and sample code describing how Open Science principles can be effectively applied to earth sciences, particularly in the context of air pollution and its impact on health.**
+
+Teaching earth science topics and coding skills in depth is beyond the scope of the SCHOOL Modules. Instead, the SCHOOL Project aims to provide users with the skills to adapt the skills learned in SCHOOL lessons to the users’ own Open Science workflow. To learn more about Open Science, explore [NASA's TOPS Open Science 101 Curriculum](https://openscience101.org/). To explore other themes in the SCHOOL project, visit our [Modules Page](https://ciesin-geospatial.github.io/TOPSTSCHOOL/modules.html).
+:::
+
+## Module 4: Climate and Agriculture datasets and use cases cover:
+
+- **NASA ACRES Climate Resilience Network:** This lesson will use public research data for derived conclusions and to view the impacts climate change has on agriculture.
+ - [Lesson 1: Acquiring, Pre-Processing, and Visualizing Student-Monitored Data for New York City schools](https://ciesin-geospatial.github.io/TOPSTSCHOOL-climate-agriculture/m401-student-led-monitoring-nyc.html)
+
+
+::: {style="text-align: right;"}
+[Start Lesson 1](https://ciesin-geospatial.github.io/TOPSTSCHOOL-climate-agriculture/m401-student-led-monitoring-nyc.html){.btn .btn-primary .btn role="button"}
+:::
+
+This course was made possible thanks to the work of our NASA Transform to Open Science (TOPS) team, our SCHOOL Open Science team, open science Subject Matter Experts (SMEs), and the SCHOOL Development team!
\ No newline at end of file
diff --git a/m401-nasa-acres-alimate-resilience-network.qmd b/m401-nasa-acres-alimate-resilience-network.qmd
new file mode 100644
index 0000000..5bb53b4
--- /dev/null
+++ b/m401-nasa-acres-alimate-resilience-network.qmd
@@ -0,0 +1,672 @@
+---
+title: "NASA Acres and Climate Resilience Network"
+author: "Jacom Orser"
+format: html
+bibliography: acres-references.bib
+---
+
+# Introduction
+
+## CLIMATE SMART AGRICULTURE IN MARYLAND
+
+A regional view of the past with applications to the present. Python spatial and tabular analysis using multiple variables to look at a single time slice. This code can be a combination of cloud and local based analysis. A purpose for this exercise is to use public research data for derived conclusions, view the impacts climate change has on agriculture, and look at what farmers can do and why we need to support them.
+
+NASA Acres has defined the **Essential Agricultural Variables (EAVs)**. These were designed by NASA Acres to define the capabilities of the top satellite data scientists and practitioners who make up NASA Acres Research, Development, and Extension partners, and the needs of decision-making-collaborators already in our network, to identify an initial set of focus.
+
+1. Cropland and Crop Type Mapping
+
+2. Crop and Crop Type Area Estimation
+
+Derived results will be the county level statistics of impacts of harvestable acreage. To achieve this, a basic workflow we can
+
+1. Use api to call in the CDL dataset to map crop types.
+
+Then
+
+1. Access the Sea level rise (elevation dataset) to identify new areas of potential areas that are at risk of future flooding. Clip to county of interest
+
+2. Compare the CDL with the SLR mask and with out to identify the NOAA estimated loss of land.
+
+Then using MODIS create a time series for NDVI of the masked crop land, derive insights about trends in NDVI
+
+1. Use NDVI to create a time series looking backwards in time at areas that have experienced flooding to visualize the moving from productive farms to moderate quality.
+
+This all together would allow us to make a predictive analysis for Maryland in the future under the projections of sea level rise. Given the current conditions, subtracting the sea level rise inundated areas. Then given the remaining areas, making predictions about the trends of NDVI given the trends in historical areas that are near impacted areas. We cannot say that "sea level rise is driving the decrease" because there are many other factor and decisions that farmers make around how much to plant, how much to harvest, chemicals applied, growing degree days, soil and soil moisture conditions. But this module should provide the ability to make insights about potential causes and impacts.
+
+The data story that we have derived is about sea level rise in Maryland and the impact that it has on the production levels with in the state. From this module, we can provide students with the ability to draw from multiple data sources, as well as derive insights using historical and future viewing data sets.
+
+The U.S. Mid-Atlantic has seen higher rates of sea level rise, marshes may be especially vulnerable. "In the Chesapeake Bay, sea level rise has already contributed to the degradation of over 80,000 ha (70%) of tidal marsh" [@taylor2020]. This view can help us to understand the impacts that small levels of sea level rise have on land.
+
+We can prompt the user to think about future impacts outside the direct sea level rise projections, allowing them to include a more full picture, and finally using that picture to identify economic impacts that action or inaction causes. This begs the question, what can the public do to enact changes, rather than putting the pressure on farmers to change?
+
+```{python}
+!pip install arcgis arcgis-mapping rasterio earthaccess
+```
+
+```{python}
+#import the required packages
+import geopandas as gpd
+import rasterio
+import numpy as np
+import pandas as pd
+import matplotlib.pyplot as plt
+from shapely.geometry import box
+import numpy as np
+import requests
+import xml.etree.ElementTree as ET
+import matplotlib.pyplot as plt
+from io import BytesIO
+import os
+from rasterio.warp import calculate_default_transform, reproject, Resampling
+from datetime import datetime, timedelta
+
+```
+
+```{python}
+#define fuctions
+
+def reproject_raster(src_dataset, src_crs, src_transform, dst_crs='EPSG:4326'):
+ dst_crs = rasterio.crs.CRS.from_string(dst_crs)
+ dst_transform, dst_width, dst_height = calculate_default_transform(
+ src_crs, dst_crs,
+ src_dataset.width, src_dataset.height,
+ *src_dataset.bounds)
+ dst_kwargs = src_dataset.meta.copy()
+ dst_kwargs.update({
+ 'crs': dst_crs,
+ 'transform': dst_transform,
+ 'width': dst_width,
+ 'height': dst_height})
+ dst_data = np.zeros((src_dataset.count, dst_height, dst_width), dtype=src_dataset.dtypes[0])
+ for i in range(1, src_dataset.count + 1):
+ reproject(
+ source=rasterio.band(src_dataset, i),
+ destination=dst_data[i-1],
+ src_transform=src_transform,
+ src_crs=src_crs,
+ dst_transform=dst_transform,
+ dst_crs=dst_crs,
+ resampling=Resampling.nearest)
+ return dst_data, dst_kwargs
+
+def format_cdl_url(fips,year):
+ base_url = "https://nassgeodata.gmu.edu/axis2/services/CDLService/GetCDLFile"
+ url = f"{base_url}?year={year}&fips={fips}"
+ return url
+
+def get_srl_raster():
+ slr_path = 'data/MD_East_slr_depth_3_5ft.tif'
+ with rasterio.open(slr_path) as slr_raster:
+ slr_meta = slr_raster.meta
+ slr_meta['nodata'] = 0
+ slr_reprojected, slr_meta_reprojected = reproject_raster(src_dataset=slr_raster, src_crs=slr_meta['crs'], src_transform=slr_meta['transform'])
+ return slr_reprojected, slr_meta_reprojected
+```
+
+# CDL
+
+To understand the nature of Maryland agriculture we can begin by using the CDL. This raster dataset comes at a 30 meter spatial resolution. We can access this data through their api. The current configuration allows data to be pulled at the county level.
+
+COUNTY LEVEL CDL FOR FIPS
+
+This view provides insight into the common commodities that are grown in Maryland in the desired year.
+
+```{python}
+
+crop_mapping = {
+ 0: "Background", 1: "Corn", 2: "Cotton", 3: "Rice",
+ 4: "Sorghum", 5: "Soybeans", 6: "Sunflower", 10: "Peanuts",
+ 11: "Tobacco", 12: "Sweet Corn", 13: "Pop or Orn Corn",
+ 14: "Mint", 21: "Barley", 22: "Durum Wheat",23: "Spring Wheat",
+ 24: "Winter Wheat", 25: "Other Small Grains", 26: "Dbl Crop WinWht/Soybeans",
+ 27: "Rye", 28: "Oats", 29: "Millet", 30: "Speltz", 31: "Canola",
+ 32: "Flaxseed", 33: "Safflower", 34: "Rape Seed",35: "Mustard",
+ 36: "Alfalfa", 37: "Other Hay/Non Alfalfa", 38: "Camelina", 39: "Buckwheat",
+ 41: "Sugarbeets", 42: "Dry Beans", 43: "Potatoes", 44: "Other Crops",
+ 45: "Sugarcane", 46: "Sweet Potatoes", 47: "Misc Vegs & Fruits",
+ 48: "Watermelons", 49: "Onions", 50: "Cucumbers", 51: "Chick Peas",
+ 52: "Lentils", 53: "Peas", 54: "Tomatoes", 55: "Caneberries", 56: "Hops",
+ 57: "Herbs", 58: "Clover/Wildflowers", 59: "Sod/Grass Seed", 60: "Switchgrass",
+ 61: "Fallow/Idle Cropland", 63: "Forest", 64: "Shrubland", 65: "Barren",
+ 66: "Cherries", 67: "Peaches", 68: "Apples", 69: "Grapes",
+ 70: "Christmas Trees", 71: "Other Tree Crops", 72: "Citrus", 74: "Pecans",
+ 75: "Almonds", 76: "Walnuts", 77: "Pears", 81: "Clouds/No Data",
+ 82: "Developed", 83: "Water", 87: "Wetlands", 88: "Nonag/Undefined",
+ 92: "Aquaculture", 111: "Open Water", 112: "Perennial Ice/Snow",
+ 121: "Developed/Open Space",122: "Developed/Low Intensity",
+ 123: "Developed/Med Intensity", 124: "Developed/High Intensity",
+ 131: "Barren", 141: "Deciduous Forest", 142: "Evergreen Forest",
+ 143: "Mixed Forest", 152: "Shrubland", 176: "Grass/Pasture",
+ 190: "Woody Wetlands", 195: "Herbaceous Wetlands", 204: "Pistachios",
+ 205: "Triticale", 206: "Carrots", 207: "Asparagus", 208: "Garlic",
+ 209: "Cantaloupes", 210: "Prunes", 211: "Olives", 212: "Oranges",
+ 213: "Honeydew Melons", 214: "Broccoli", 215: "Avocados", 216: "Peppers",
+ 217: "Pomegranates", 218: "Nectarines", 219: "Greens", 220: "Plums",
+ 221: "Strawberries", 222: "Squash", 223: "Apricots", 224: "Vetch",
+ 225: "Dbl Crop WinWht/Corn", 226: "Dbl Crop Oats/Corn", 227: "Lettuce",
+ 228: "Dbl Crop Triticale/Corn", 229: "Pumpkins",
+ 230: "Dbl Crop Lettuce/Durum Wht", 231: "Dbl Crop Lettuce/Cantaloupe",
+ 232: "Dbl Crop Lettuce/Cotton",233: "Dbl Crop Lettuce/Barley",
+ 234: "Dbl Crop Durum Wht/Sorghum", 235: "Dbl Crop Barley/Sorghum",
+ 236: "Dbl Crop WinWht/Sorghum", 237: "Dbl Crop Barley/Corn",
+ 238: "Dbl Crop WinWht/Cotton",239: "Dbl Crop Soybeans/Cotton",
+ 240: "Dbl Crop Soybeans/Oats", 241: "Dbl Crop Corn/Soybeans",
+ 242: "Blueberries", 243: "Cabbage", 244: "Cauliflower", 245: "Celery",
+ 246: "Radishes", 247: "Turnips", 248: "Eggplants",249: "Gourds",
+ 250: "Cranberries", 254: "Dbl Crop Barley/Soybeans", 255: ""
+ }
+
+```
+
+```{python}
+
+county = 'Talbot'
+code = "24041"
+year = 2018
+
+```
+
+```{python}
+service_url = format_cdl_url(code, year)
+response = requests.get(service_url)
+root = ET.fromstring(response.content)
+tif_url = root.find('.//returnURL').text
+cdl_data = rasterio.open(tif_url)
+cdl_meta = cdl_data.meta
+cdl_reprojected, cdl_meta_reprojected = reproject_raster(src_dataset=cdl_data, src_crs=cdl_meta['crs'], src_transform=cdl_meta['transform'])
+
+```
+
+The eastern seaboard has seen changes in the sea level along with increase flooding and the decrease in predictability of water flow from rivers. This daa comes from NOAA (learn more here )
+
+This raster is the projection of inudated areas by the year of 2050. Data access can be found here
+
+```{python}
+slr_reprojected, slr_meta_reprojected = get_srl_raster()
+
+```
+
+# Gather NDVI values
+
+```{python}
+def get_appeears_token(username, password):
+
+ response = requests.post('https://appeears.earthdatacloud.nasa.gov/api/login',auth=(username, password))
+
+ return response.json()['token']
+
+def submit_appeears_task(token, task):
+
+ response = requests.post('https://appeears.earthdatacloud.nasa.gov/api/task',json=task,headers={'Authorization': f'Bearer {token}'})
+
+ return response.json()['task_id']
+
+def get_appeears_bundle(token, task_id):
+
+ response = requests.get(f'https://appeears.earthdatacloud.nasa.gov/api/bundle/{task_id}',headers={'Authorization': f'Bearer {token}'})
+
+ return response.json()
+
+def load_county_data():
+
+ sea_level_path = 'data/Talbot_County.geojson'
+
+ sea_level_gdf = gpd.read_file(sea_level_path)
+
+ return sea_level_gdf
+
+def day_of_year_to_date2(day_of_year_str):
+
+ year = int(day_of_year_str[:4])
+
+ day_of_year = int(day_of_year_str[4:])
+
+ start_date = datetime(year, 1, 1)
+
+ target_date = start_date + timedelta(days=day_of_year - 1)
+
+ return target_date.strftime('%m-%d-%Y')
+
+def resample_using_cdl(data, data_meta, cdl_meta_reprojected):
+
+ resampled_data = np.zeros(
+
+ (data.shape[0], cdl_meta_reprojected['height'],
+ cdl_meta_reprojected['width']),
+
+ dtype=data.dtype
+ )
+
+ for i in range(data.shape[0]):
+
+ reproject(
+
+ source=data[i],
+
+ destination=resampled_data[i],
+
+ src_transform=data_meta['transform'],
+
+ src_crs=data_meta['crs'],
+
+ dst_transform=cdl_meta_reprojected['transform'],
+
+ dst_crs=cdl_meta_reprojected['crs'],
+
+ resampling=Resampling.bilinear
+
+ )
+ updated_meta = data_meta.copy()
+
+ updated_meta.update({
+
+ 'transform': cdl_meta_reprojected['transform'],
+
+ 'width': cdl_meta_reprojected['width'],
+
+ 'height': cdl_meta_reprojected['height'],
+
+ 'crs': cdl_meta_reprojected['crs']
+
+ })
+
+ return resampled_data, updated_meta
+```
+
+```{python}
+county_data = load_county_data()
+
+county_reprojected = county_data.to_crs('EPSG:4326')
+
+us = county_reprojected.geometry.iloc[0]
+
+minx, miny, maxx, maxy = us.bounds
+```
+
+```{python}
+bbox = [[maxx, miny],
+
+ [maxx, maxy],
+
+ [minx, maxy],
+
+ [minx, maxy],
+
+ [minx, miny]]
+```
+
+```{python}
+start = f'06-01-{year}'
+
+end = f'09-30-{year}'
+
+```
+
+```{python}
+products = [{'layer': '_250m_16_days_NDVI', 'product':'MOD13Q1.061'}]
+
+
+
+```
+
+
+```{python}
+USERNAME = 'USERNAME'
+PASSWORD = 'PASSWORD'
+
+```
+
+```{python}
+#| eval: true
+#| echo: false
+
+USERNAME = 'jfmartinez4'
+PASSWORD = 'Potatoaway1!'
+```
+```{python}
+token = get_appeears_token(f'{USERNAME}', f'{PASSWORD}')
+```
+
+```{python}
+task = {
+ 'task_type': 'area',
+
+ 'task_name': 'MODIS',
+
+ 'params': {
+ "geo": {
+ "type": "FeatureCollection",
+ "features": [{
+ "type": "Feature",
+ "geometry": {
+ "type": "MultiPolygon",
+ "coordinates": [[bbox]]},
+ "properties": {}}],
+ "fileName": "Polygon"},
+
+ 'dates': [{
+ 'startDate': start, 'endDate': end}],
+
+ 'layers': products,
+
+ 'coordinates': bbox,
+
+ "output": {
+ "format": {
+ "type": "geotiff"},
+ "projection": "native"}
+
+ }}
+
+```
+
+```{python}
+task_id = submit_appeears_task(token, task)
+```
+```{python}
+#| eval: true
+#| echo: false
+
+task_id = 'd09681c4-2da0-4513-a721-a041285438bf'
+```
+
+```{python}
+bundle = get_appeears_bundle(token, task_id)
+bundle
+```
+```{python}
+data = {}
+
+for file in bundle['files']:
+
+ file_id = file['file_id']
+
+ if 'NDVI' in file['file_name']:
+
+ if '_doy' in file['file_name']:
+
+ datesy = file['file_name'].split('_doy')[1][:7]
+
+ doy = day_of_year_to_date2(datesy)
+
+ print(file['file_name'])
+
+ else:
+
+ continue
+
+ file_download = requests.get(
+
+ 'https://appeears.earthdatacloud.nasa.gov/api/bundle/{0}/{1}'.format(task_id, file_id),
+
+ headers={'Authorization': 'Bearer {0}'.format(token)},
+
+ allow_redirects=True,
+
+ stream=True)
+
+ file_download.raise_for_status()
+
+ if not file_download.content:
+
+ print(f"Warning: Empty file downloaded for {file['file_name']}")
+
+ continue
+
+ file_content = BytesIO(file_download.content)
+
+ with rasterio.open(file_content) as src_initial:
+
+ src = src_initial.read(1, masked=True)
+
+ src_meta = src_initial.meta
+
+ dst_crs = cdl_data.crs
+
+ transform, width, height = calculate_default_transform(
+
+ src_meta['crs'], dst_crs, src.shape[1], src.shape[0], *src_initial.bounds)
+
+ kwargs = src_meta.copy()
+
+ kwargs.update({
+
+ 'crs': dst_crs,
+
+ 'transform': transform,
+
+ 'width': width,
+
+ 'height': height
+
+ })
+
+ dst = np.zeros((height, width), dtype=src_meta['dtype'])
+
+ reproject(
+
+ source=src,
+
+ destination=dst,
+
+ src_transform=src_meta['transform'],
+
+ src_crs=src_meta['crs'],
+
+ dst_transform=transform,
+
+ dst_crs=dst_crs,
+
+ resampling=Resampling.nearest)
+
+ kwargs.update({'nodata': 0})
+
+ key = 'MODIS_NDVI'
+
+ if key not in data:
+
+ data[key] = {'data': [], 'meta': kwargs, 'doy': []}
+
+ data[key]['data'].append(dst)
+
+ data[key]['doy'].append(doy)
+
+ else:
+
+ continue
+
+```
+
+```{python}
+ndvi = np.stack(data['MODIS_NDVI']['data']) / 100
+
+day = np.stack(data['MODIS_NDVI']['doy'])
+
+meta = data['MODIS_NDVI']['meta']
+
+```
+
+## Combine and compare the results
+
+```{python}
+ndvi_resampled, ndvi_resampled_meta = resample_using_cdl(ndvi, meta, cdl_meta_reprojected)
+
+slr_resampled, slr_resampled_meta = resample_using_cdl(slr_reprojected, slr_meta_reprojected, cdl_meta_reprojected)
+
+print("CDL shape:", cdl_reprojected.shape)
+
+print("Original NDVI shape:", ndvi.shape)
+
+print("Resampled NDVI shape:", ndvi_resampled.shape)
+
+print("Original SLR shape:", slr_reprojected.shape)
+
+print("Resampled SLR shape:", slr_resampled.shape)
+
+```
+
+```{python}
+def mask_rasters(cdl_raster, slr_raster, ndvi_rasters):
+
+ corn_mask = cdl_raster == 1
+
+ ABOVE_water_mask = slr_raster <= 1
+
+ corn_above_water_mask = corn_mask & ~ABOVE_water_mask
+
+ corn_below_water_mask = corn_mask & ABOVE_water_mask
+
+ corn_above_water_mask = corn_above_water_mask[0,:,:]
+
+ corn_below_water_mask = corn_below_water_mask[0,:,:]
+
+ masked_ndvi_below_water = ndvi_rasters[:, corn_below_water_mask]
+
+ masked_ndvi_above_water = ndvi_rasters[:, corn_above_water_mask]
+
+ return {
+
+ 'corn_below_water': {
+
+ 'mask': corn_below_water_mask,
+
+ 'ndvi': masked_ndvi_below_water,
+
+ 'mean': np.mean(masked_ndvi_below_water, axis=1) if masked_ndvi_below_water.size > 0 else None,
+
+ 'min': np.min(masked_ndvi_below_water, axis=1) if masked_ndvi_below_water.size > 0 else None,
+
+ 'max': np.max(masked_ndvi_below_water, axis=1) if masked_ndvi_below_water.size > 0 else None
+
+ },
+
+ 'corn_above_water': {
+
+ 'mask': corn_above_water_mask,
+
+ 'ndvi': masked_ndvi_above_water,
+
+ 'mean': np.mean(masked_ndvi_above_water, axis=1) if masked_ndvi_above_water.size > 0 else None,
+
+ 'min': np.min(masked_ndvi_above_water, axis=1) if masked_ndvi_above_water.size > 0 else None,
+
+ 'max': np.max(masked_ndvi_above_water, axis=1) if masked_ndvi_above_water.size > 0 else None
+
+ }
+
+ }
+
+```
+
+```{python}
+def visualize_ndvi_analysis(ndvi_analysis, days):
+
+ plt.figure(figsize=(12, 6))
+
+ plt.subplot(1, 2, 1)
+
+ plt.plot(days, ndvi_analysis['corn_above_water']['mean'], label='Corn Above Water', color='green')
+
+ plt.plot(days, ndvi_analysis['corn_below_water']['mean'], label='Corn Below Water', color='blue')
+
+ plt.title('Mean NDVI')
+
+ plt.xlabel('Days')
+
+ plt.xticks(rotation=90)
+
+ plt.ylabel('NDVI')
+
+ plt.legend()
+
+ plt.subplot(1, 2, 2)
+
+ plt.plot(days, ndvi_analysis['corn_above_water']['max'], label='Corn Above Water', color='green')
+
+ plt.plot(days, ndvi_analysis['corn_below_water']['max'], label='Corn Below Water', color='blue')
+
+ plt.title('Maximum NDVI')
+
+ plt.xlabel('Days')
+
+ plt.xticks(rotation=90)
+
+ plt.ylabel('NDVI')
+
+ plt.legend()
+
+ plt.tight_layout()
+
+ plt.show()
+
+```
+
+```{python}
+ndvi_analysis = mask_rasters(cdl_reprojected, slr_resampled, ndvi_resampled)
+
+visualize_ndvi_analysis(ndvi_analysis, day)
+```
+
+```{python}
+for condition in ['corn_below_water', 'corn_above_water']:
+
+ print(f"\n{condition.replace('_', ' ').title()} Corn Analysis per MODIS:")
+
+ for stat in ['mean', 'max']:
+
+ print(f"{stat.capitalize()} NDVI: {ndvi_analysis[condition][stat]}")
+
+```
+
+```{python}
+def visualize_ndvi_histogram(ndvi_analysis):
+
+ plt.figure(figsize=(12, 6))
+
+ plt.hist(ndvi_analysis['corn_below_water']['ndvi'].flatten(), bins=30, alpha=0.5, label='Corn Below Water', color='blue')
+
+ plt.hist(ndvi_analysis['corn_above_water']['ndvi'].flatten(), bins=30, alpha=0.5, label='Corn Above Water', color='green')
+
+ plt.xlabel("NDVI values")
+
+ plt.ylabel("Frequency")
+
+ plt.title("Histogram of NDVI values for corn")
+
+ plt.legend(loc='upper right')
+
+visualize_ndvi_histogram(ndvi_analysis)
+
+```
+
+What are the insights that we can derive from this module?
+
+What are the impacts of land being close to the water currently?
+
+What is the potential loss of productive land?
+
+Where could there be errors given the sampling methods?
+
+# What to do next
+
+1. Clip the MODIS data to the CDL data of only the agricultural values of interest
+2. Find the average NDVI of these values
+3. Clip this to the SLR area
+4. Find the average NDVI of these values
+5. Compare the changes for the year of 2018
+
+This will tell us the impacted area that could be effected by the floods
+
+As we look into the future we should see this change become more drastic
+
+# ABOUT NASA ACRES AND THE CLIMATE RESILIENCE NETWORK
+
+## NASA Acres
+
+NASA Acres is NASA's U.S.-focused agriculture and food security Consortium. NASA Acres is commissioned under NASA's Applied Sciences Program and led by the University of Maryland, with over 30 partner institutions and growing. The consortium approach brings together public and private sectors, and allows for flexible partnerships and rapid action in delivering the benefits of NASA data and tools to U.S. agriculture decision makers. We work with people across the agriculture sector to develop Earth observation (EO)-based data and tools that strengthen U.S. agriculture.
+
+## Climate Resilience Network
+
+Climate Resilience Network has the goal of unifying the University of Maryland's groundbreaking research with the immediate requirements of Maryland's decision-makers to address the complexities of climate change as it affects our home. GEOG's Grand Challenges project is aimed at promoting climate-resilient agriculture in Maryland, creating interactive tools designed to support research based decision-making.italicized text
+
+COUNTY LEVEL SHAPEFILE
diff --git a/m401-nasa-acres-climate-resilience-network.ipynb b/m401-nasa-acres-climate-resilience-network.ipynb
new file mode 100644
index 0000000..81b9229
--- /dev/null
+++ b/m401-nasa-acres-climate-resilience-network.ipynb
@@ -0,0 +1,872 @@
+{
+ "cells": [
+ {
+ "cell_type": "raw",
+ "metadata": {},
+ "source": [
+ "---\n",
+ "title: \"NASA Acres and Climate Resilience Network\"\n",
+ "author: \"Jacom Orser\"\n",
+ "format: html\n",
+ "bibliography: acres-references.bib\n",
+ "---"
+ ],
+ "id": "e8554d2a"
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "# Introduction\n",
+ "\n",
+ "## CLIMATE SMART AGRICULTURE IN MARYLAND\n",
+ "\n",
+ "A regional view of the past with applications to the present. Python spatial and tabular analysis using multiple variables to look at a single time slice. This code can be a combination of cloud and local based analysis. A purpose for this exercise is to use use public research data for derived conclusions view the impacts climate change has on agriculture a look at what farmers can do and why we need to support them.\n",
+ "\n",
+ "NASA Acres has defined the Essential Agricultural Variables or EAVs. These were designed by NASA Acres to define the capabilities of the top satellite data scientists and practitioners who make up NASA Acres Research, Development, and Extension partners, and the needs of decision-making-collaborators already in our network, to identify an initial set of focus.\n",
+ "\n",
+ "1. Cropland and Crop Type Mapping\n",
+ "\n",
+ "2. Crop and Crop Type Area Estimation\n",
+ "\n",
+ "Derived results will be the county level statistics of impacts of harvestable acreage. To achieve this, a basic workflow we can\n",
+ "\n",
+ "1. Use api to call in the CDL dataset to map crop types.\n",
+ "\n",
+ "Then\n",
+ "\n",
+ "1. Access the Sea level rise (elevation dataset) to identify new areas of potential areas that are at risk of future flooding. Clip to county of interest\n",
+ "\n",
+ "2. Compare the CDL with the SLR mask and with out to identify the NOAA estimated loss of land.\n",
+ "\n",
+ "Then using MODIS create a time series for NDVI of the masked crop land, derive insights about trends in NDVI\n",
+ "\n",
+ "1. Use NDVI to create a time series looking backwards in time at areas that have experienced flooding to visualize the moving from productive farms to moderate quality.\n",
+ "\n",
+ "This all together would allow us to make a predictive analysis for Maryland in the future under the projections of sea level rise. Given the current conditions, subtracting the sea level rise inundated areas. Then given the remaining areas, making predictions about the trends of NDVI given the trends in historical areas that are near impacted areas. We cannot say that \"sea level rise is driving the decrease\" because there are many other factor and decisions that farmers make around how much to plant, how much to harvest, chemicals applied, growing degree days, soil and soil moisture conditions. But this module should provide the ability to make insights about potential causes and impacts.\n",
+ "\n",
+ "The data story that we have derived is about sea level rise in Maryland and the impact that it has on the production levels with in the state. From this module, we can provide students with the ability to draw from multiple data sources, as well as derive insights using historical and future viewing data sets.\n",
+ "\n",
+ "The U.S. Mid-Atlantic has seen higher rates of sea level rise, marshes may be especially vulnerable. \"In the Chesapeake Bay, sea level rise has already contributed to the degradation of over 80,000 ha (70%) of tidal marsh\" [@taylor2020]. This view can help us to understand the impacts that small levels of sea level rise have on land.\n",
+ "\n",
+ "We can prompt the user to think about future impacts outside the direct sea level rise projections, allowing them to include a more full picture, and finally using that picture to identify economic impacts that action or inaction causes. This begs the question, what can the public do to enact changes, rather than putting the pressure on farmers to change?\n"
+ ],
+ "id": "d264226d"
+ },
+ {
+ "cell_type": "code",
+ "metadata": {},
+ "source": [
+ "!pip install arcgis arcgis-mapping rasterio earthaccess"
+ ],
+ "id": "99fb08e5",
+ "execution_count": null,
+ "outputs": []
+ },
+ {
+ "cell_type": "code",
+ "metadata": {},
+ "source": [
+ "#import the required packages\n",
+ "import geopandas as gpd\n",
+ "import rasterio\n",
+ "import numpy as np\n",
+ "import pandas as pd\n",
+ "import matplotlib.pyplot as plt\n",
+ "from shapely.geometry import box\n",
+ "import numpy as np\n",
+ "import requests\n",
+ "import xml.etree.ElementTree as ET\n",
+ "import matplotlib.pyplot as plt\n",
+ "from io import BytesIO\n",
+ "import os\n",
+ "from rasterio.warp import calculate_default_transform, reproject, Resampling\n",
+ "from datetime import datetime, timedelta"
+ ],
+ "id": "0f48fb31",
+ "execution_count": null,
+ "outputs": []
+ },
+ {
+ "cell_type": "code",
+ "metadata": {},
+ "source": [
+ "#define fuctions\n",
+ "\n",
+ "def reproject_raster(src_dataset, src_crs, src_transform, dst_crs='EPSG:4326'):\n",
+ " dst_crs = rasterio.crs.CRS.from_string(dst_crs)\n",
+ " dst_transform, dst_width, dst_height = calculate_default_transform(\n",
+ " src_crs, dst_crs,\n",
+ " src_dataset.width, src_dataset.height,\n",
+ " *src_dataset.bounds)\n",
+ " dst_kwargs = src_dataset.meta.copy()\n",
+ " dst_kwargs.update({\n",
+ " 'crs': dst_crs,\n",
+ " 'transform': dst_transform,\n",
+ " 'width': dst_width,\n",
+ " 'height': dst_height})\n",
+ " dst_data = np.zeros((src_dataset.count, dst_height, dst_width), dtype=src_dataset.dtypes[0])\n",
+ " for i in range(1, src_dataset.count + 1):\n",
+ " reproject(\n",
+ " source=rasterio.band(src_dataset, i),\n",
+ " destination=dst_data[i-1],\n",
+ " src_transform=src_transform,\n",
+ " src_crs=src_crs,\n",
+ " dst_transform=dst_transform,\n",
+ " dst_crs=dst_crs,\n",
+ " resampling=Resampling.nearest)\n",
+ " return dst_data, dst_kwargs\n",
+ "\n",
+ "def format_cdl_url(fips,year):\n",
+ " base_url = \"https://nassgeodata.gmu.edu/axis2/services/CDLService/GetCDLFile\"\n",
+ " url = f\"{base_url}?year={year}&fips={fips}\"\n",
+ " return url\n",
+ "\n",
+ "def get_srl_raster():\n",
+ " slr_path = 'data/MD_East_slr_depth_3_5ft.tif'\n",
+ " with rasterio.open(slr_path) as slr_raster:\n",
+ " slr_meta = slr_raster.meta\n",
+ " slr_meta['nodata'] = 0\n",
+ " slr_reprojected, slr_meta_reprojected = reproject_raster(src_dataset=slr_raster, src_crs=slr_meta['crs'], src_transform=slr_meta['transform'])\n",
+ " return slr_reprojected, slr_meta_reprojected"
+ ],
+ "id": "06b56597",
+ "execution_count": null,
+ "outputs": []
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "# CDL \n",
+ "\n",
+ "To understand the nature of Maryland agriculture we can begin by using the CDL. This raster dataset comes at a 30 meter spatial resolution. We can access this data through their api. The current configuration allows data to be pulled at the county level.\n",
+ "\n",
+ "COUNTY LEVEL CDL FOR FIPS \n",
+ "\n",
+ "This view provides insight into the common commodities that are grown in Maryland in the desired year.\n"
+ ],
+ "id": "c7767154"
+ },
+ {
+ "cell_type": "code",
+ "metadata": {},
+ "source": [
+ "crop_mapping = {\n",
+ " 0: \"Background\", 1: \"Corn\", 2: \"Cotton\", 3: \"Rice\", \n",
+ " 4: \"Sorghum\", 5: \"Soybeans\", 6: \"Sunflower\", 10: \"Peanuts\",\n",
+ " 11: \"Tobacco\", 12: \"Sweet Corn\", 13: \"Pop or Orn Corn\", \n",
+ " 14: \"Mint\", 21: \"Barley\", 22: \"Durum Wheat\",23: \"Spring Wheat\", \n",
+ " 24: \"Winter Wheat\", 25: \"Other Small Grains\", 26: \"Dbl Crop WinWht/Soybeans\", \n",
+ " 27: \"Rye\", 28: \"Oats\", 29: \"Millet\", 30: \"Speltz\", 31: \"Canola\", \n",
+ " 32: \"Flaxseed\", 33: \"Safflower\", 34: \"Rape Seed\",35: \"Mustard\", \n",
+ " 36: \"Alfalfa\", 37: \"Other Hay/Non Alfalfa\", 38: \"Camelina\", 39: \"Buckwheat\", \n",
+ " 41: \"Sugarbeets\", 42: \"Dry Beans\", 43: \"Potatoes\", 44: \"Other Crops\", \n",
+ " 45: \"Sugarcane\", 46: \"Sweet Potatoes\", 47: \"Misc Vegs & Fruits\",\n",
+ " 48: \"Watermelons\", 49: \"Onions\", 50: \"Cucumbers\", 51: \"Chick Peas\", \n",
+ " 52: \"Lentils\", 53: \"Peas\", 54: \"Tomatoes\", 55: \"Caneberries\", 56: \"Hops\", \n",
+ " 57: \"Herbs\", 58: \"Clover/Wildflowers\", 59: \"Sod/Grass Seed\", 60: \"Switchgrass\",\n",
+ " 61: \"Fallow/Idle Cropland\", 63: \"Forest\", 64: \"Shrubland\", 65: \"Barren\", \n",
+ " 66: \"Cherries\", 67: \"Peaches\", 68: \"Apples\", 69: \"Grapes\", \n",
+ " 70: \"Christmas Trees\", 71: \"Other Tree Crops\", 72: \"Citrus\", 74: \"Pecans\",\n",
+ " 75: \"Almonds\", 76: \"Walnuts\", 77: \"Pears\", 81: \"Clouds/No Data\", \n",
+ " 82: \"Developed\", 83: \"Water\", 87: \"Wetlands\", 88: \"Nonag/Undefined\", \n",
+ " 92: \"Aquaculture\", 111: \"Open Water\", 112: \"Perennial Ice/Snow\", \n",
+ " 121: \"Developed/Open Space\",122: \"Developed/Low Intensity\",\n",
+ " 123: \"Developed/Med Intensity\", 124: \"Developed/High Intensity\", \n",
+ " 131: \"Barren\", 141: \"Deciduous Forest\", 142: \"Evergreen Forest\", \n",
+ " 143: \"Mixed Forest\", 152: \"Shrubland\", 176: \"Grass/Pasture\",\n",
+ " 190: \"Woody Wetlands\", 195: \"Herbaceous Wetlands\", 204: \"Pistachios\", \n",
+ " 205: \"Triticale\", 206: \"Carrots\", 207: \"Asparagus\", 208: \"Garlic\", \n",
+ " 209: \"Cantaloupes\", 210: \"Prunes\", 211: \"Olives\", 212: \"Oranges\",\n",
+ " 213: \"Honeydew Melons\", 214: \"Broccoli\", 215: \"Avocados\", 216: \"Peppers\", \n",
+ " 217: \"Pomegranates\", 218: \"Nectarines\", 219: \"Greens\", 220: \"Plums\", \n",
+ " 221: \"Strawberries\", 222: \"Squash\", 223: \"Apricots\", 224: \"Vetch\",\n",
+ " 225: \"Dbl Crop WinWht/Corn\", 226: \"Dbl Crop Oats/Corn\", 227: \"Lettuce\", \n",
+ " 228: \"Dbl Crop Triticale/Corn\", 229: \"Pumpkins\", \n",
+ " 230: \"Dbl Crop Lettuce/Durum Wht\", 231: \"Dbl Crop Lettuce/Cantaloupe\", \n",
+ " 232: \"Dbl Crop Lettuce/Cotton\",233: \"Dbl Crop Lettuce/Barley\", \n",
+ " 234: \"Dbl Crop Durum Wht/Sorghum\", 235: \"Dbl Crop Barley/Sorghum\", \n",
+ " 236: \"Dbl Crop WinWht/Sorghum\", 237: \"Dbl Crop Barley/Corn\", \n",
+ " 238: \"Dbl Crop WinWht/Cotton\",239: \"Dbl Crop Soybeans/Cotton\", \n",
+ " 240: \"Dbl Crop Soybeans/Oats\", 241: \"Dbl Crop Corn/Soybeans\", \n",
+ " 242: \"Blueberries\", 243: \"Cabbage\", 244: \"Cauliflower\", 245: \"Celery\", \n",
+ " 246: \"Radishes\", 247: \"Turnips\", 248: \"Eggplants\",249: \"Gourds\", \n",
+ " 250: \"Cranberries\", 254: \"Dbl Crop Barley/Soybeans\", 255: \"\"\n",
+ " }"
+ ],
+ "id": "e9eab384",
+ "execution_count": null,
+ "outputs": []
+ },
+ {
+ "cell_type": "code",
+ "metadata": {},
+ "source": [
+ "county = 'Talbot'\n",
+ "code = \"24041\"\n",
+ "year = 2018"
+ ],
+ "id": "a45f2697",
+ "execution_count": null,
+ "outputs": []
+ },
+ {
+ "cell_type": "code",
+ "metadata": {},
+ "source": [
+ "service_url = format_cdl_url(code, year)\n",
+ "response = requests.get(service_url)\n",
+ "root = ET.fromstring(response.content)\n",
+ "tif_url = root.find('.//returnURL').text\n",
+ "cdl_data = rasterio.open(tif_url)\n",
+ "cdl_meta = cdl_data.meta\n",
+ "cdl_reprojected, cdl_meta_reprojected = reproject_raster(src_dataset=cdl_data, src_crs=cdl_meta['crs'], src_transform=cdl_meta['transform'])"
+ ],
+ "id": "3411f259",
+ "execution_count": null,
+ "outputs": []
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "The eastern seaboard has seen changes in the sea level along with increase flooding and the decrease in predictability of water flow from rivers. This daa comes from NOAA (learn more here )\n",
+ "\n",
+ "This raster is the projection of inudated areas by the year of 2050. Data access can be found here \n"
+ ],
+ "id": "6b3d4753"
+ },
+ {
+ "cell_type": "code",
+ "metadata": {},
+ "source": [
+ "slr_reprojected, slr_meta_reprojected = get_srl_raster()"
+ ],
+ "id": "dbd4c962",
+ "execution_count": null,
+ "outputs": []
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "# Gather NDVI values\n"
+ ],
+ "id": "1369bdc7"
+ },
+ {
+ "cell_type": "code",
+ "metadata": {},
+ "source": [
+ "def get_appeears_token(username, password):\n",
+ "\n",
+ " response = requests.post('https://appeears.earthdatacloud.nasa.gov/api/login',auth=(username, password))\n",
+ "\n",
+ " return response.json()['token']\n",
+ "\n",
+ "def submit_appeears_task(token, task):\n",
+ "\n",
+ " response = requests.post('https://appeears.earthdatacloud.nasa.gov/api/task',json=task,headers={'Authorization': f'Bearer {token}'})\n",
+ "\n",
+ " return response.json()['task_id']\n",
+ "\n",
+ "def get_appeears_bundle(token, task_id):\n",
+ "\n",
+ " response = requests.get(f'https://appeears.earthdatacloud.nasa.gov/api/bundle/{task_id}',headers={'Authorization': f'Bearer {token}'})\n",
+ "\n",
+ " return response.json()\n",
+ "\n",
+ "def load_county_data():\n",
+ "\n",
+ " sea_level_path = 'data/Talbot_County.geojson'\n",
+ "\n",
+ " sea_level_gdf = gpd.read_file(sea_level_path)\n",
+ "\n",
+ " return sea_level_gdf\n",
+ "\n",
+ "def day_of_year_to_date2(day_of_year_str):\n",
+ "\n",
+ " year = int(day_of_year_str[:4])\n",
+ "\n",
+ " day_of_year = int(day_of_year_str[4:])\n",
+ "\n",
+ " start_date = datetime(year, 1, 1)\n",
+ "\n",
+ " target_date = start_date + timedelta(days=day_of_year - 1)\n",
+ "\n",
+ " return target_date.strftime('%m-%d-%Y')\n",
+ "\n",
+ "def resample_using_cdl(data, data_meta, cdl_meta_reprojected):\n",
+ "\n",
+ " resampled_data = np.zeros(\n",
+ "\n",
+ " (data.shape[0], cdl_meta_reprojected['height'],\n",
+ " cdl_meta_reprojected['width']),\n",
+ "\n",
+ " dtype=data.dtype \n",
+ " )\n",
+ "\n",
+ " for i in range(data.shape[0]):\n",
+ "\n",
+ " reproject(\n",
+ "\n",
+ " source=data[i],\n",
+ "\n",
+ " destination=resampled_data[i],\n",
+ "\n",
+ " src_transform=data_meta['transform'],\n",
+ "\n",
+ " src_crs=data_meta['crs'],\n",
+ "\n",
+ " dst_transform=cdl_meta_reprojected['transform'],\n",
+ "\n",
+ " dst_crs=cdl_meta_reprojected['crs'],\n",
+ "\n",
+ " resampling=Resampling.bilinear\n",
+ "\n",
+ " )\n",
+ " updated_meta = data_meta.copy()\n",
+ "\n",
+ " updated_meta.update({\n",
+ "\n",
+ " 'transform': cdl_meta_reprojected['transform'],\n",
+ "\n",
+ " 'width': cdl_meta_reprojected['width'],\n",
+ "\n",
+ " 'height': cdl_meta_reprojected['height'],\n",
+ "\n",
+ " 'crs': cdl_meta_reprojected['crs']\n",
+ "\n",
+ " })\n",
+ "\n",
+ " return resampled_data, updated_meta"
+ ],
+ "id": "20eef41b",
+ "execution_count": null,
+ "outputs": []
+ },
+ {
+ "cell_type": "code",
+ "metadata": {},
+ "source": [
+ "county_data = load_county_data()\n",
+ "\n",
+ "county_reprojected = county_data.to_crs('EPSG:4326')\n",
+ "\n",
+ "us = county_reprojected.geometry.iloc[0]\n",
+ "\n",
+ "minx, miny, maxx, maxy = us.bounds"
+ ],
+ "id": "d8e89928",
+ "execution_count": null,
+ "outputs": []
+ },
+ {
+ "cell_type": "code",
+ "metadata": {},
+ "source": [
+ "bbox = [[maxx, miny],\n",
+ "\n",
+ " [maxx, maxy],\n",
+ "\n",
+ " [minx, maxy],\n",
+ "\n",
+ " [minx, maxy],\n",
+ "\n",
+ " [minx, miny]]"
+ ],
+ "id": "605536f8",
+ "execution_count": null,
+ "outputs": []
+ },
+ {
+ "cell_type": "code",
+ "metadata": {},
+ "source": [
+ "start = f'06-01-{year}'\n",
+ "\n",
+ "end = f'09-30-{year}'"
+ ],
+ "id": "6b9805e9",
+ "execution_count": null,
+ "outputs": []
+ },
+ {
+ "cell_type": "code",
+ "metadata": {},
+ "source": [
+ "products = [{'layer': '_250m_16_days_NDVI', 'product':'MOD13Q1.061'}]\n"
+ ],
+ "id": "8425b37d",
+ "execution_count": null,
+ "outputs": []
+ },
+ {
+ "cell_type": "code",
+ "metadata": {},
+ "source": [
+ "USERNAME = 'USERNAME'\n",
+ "PASSWORD = 'PASSWORD'"
+ ],
+ "id": "4ce04aef",
+ "execution_count": null,
+ "outputs": []
+ },
+ {
+ "cell_type": "code",
+ "metadata": {},
+ "source": [
+ "#| eval: true\n",
+ "#| echo: false\n",
+ "\n",
+ "USERNAME = 'jfmartinez4'\n",
+ "PASSWORD = 'Potatoaway1!'"
+ ],
+ "id": "184f042c",
+ "execution_count": null,
+ "outputs": []
+ },
+ {
+ "cell_type": "code",
+ "metadata": {},
+ "source": [
+ "token = get_appeears_token(f'{USERNAME}', f'{PASSWORD}')"
+ ],
+ "id": "6d289dad",
+ "execution_count": null,
+ "outputs": []
+ },
+ {
+ "cell_type": "code",
+ "metadata": {},
+ "source": [
+ "task = {\n",
+ " 'task_type': 'area',\n",
+ " \n",
+ " 'task_name': 'MODIS',\n",
+ " \n",
+ " 'params': {\n",
+ " \"geo\": {\n",
+ " \"type\": \"FeatureCollection\", \n",
+ " \"features\": [{\n",
+ " \"type\": \"Feature\",\n",
+ " \"geometry\": {\n",
+ " \"type\": \"MultiPolygon\",\n",
+ " \"coordinates\": [[bbox]]},\n",
+ " \"properties\": {}}],\n",
+ " \"fileName\": \"Polygon\"},\n",
+ " \n",
+ " 'dates': [{\n",
+ " 'startDate': start, 'endDate': end}],\n",
+ " \n",
+ " 'layers': products,\n",
+ " \n",
+ " 'coordinates': bbox,\n",
+ " \n",
+ " \"output\": {\n",
+ " \"format\": {\n",
+ " \"type\": \"geotiff\"},\n",
+ " \"projection\": \"native\"}\n",
+ " \n",
+ " }}"
+ ],
+ "id": "34130174",
+ "execution_count": null,
+ "outputs": []
+ },
+ {
+ "cell_type": "code",
+ "metadata": {},
+ "source": [
+ "task_id = submit_appeears_task(token, task)"
+ ],
+ "id": "7e70c15e",
+ "execution_count": null,
+ "outputs": []
+ },
+ {
+ "cell_type": "code",
+ "metadata": {},
+ "source": [
+ "#| eval: true\n",
+ "#| echo: false\n",
+ "\n",
+ "task_id = 'd09681c4-2da0-4513-a721-a041285438bf'"
+ ],
+ "id": "48fb424d",
+ "execution_count": null,
+ "outputs": []
+ },
+ {
+ "cell_type": "code",
+ "metadata": {},
+ "source": [
+ "bundle = get_appeears_bundle(token, task_id)\n",
+ "bundle"
+ ],
+ "id": "713e982e",
+ "execution_count": null,
+ "outputs": []
+ },
+ {
+ "cell_type": "code",
+ "metadata": {},
+ "source": [
+ "data = {}\n",
+ "\n",
+ "for file in bundle['files']:\n",
+ "\n",
+ " file_id = file['file_id']\n",
+ "\n",
+ " if 'NDVI' in file['file_name']:\n",
+ "\n",
+ " if '_doy' in file['file_name']:\n",
+ "\n",
+ " datesy = file['file_name'].split('_doy')[1][:7]\n",
+ "\n",
+ " doy = day_of_year_to_date2(datesy)\n",
+ "\n",
+ " print(file['file_name'])\n",
+ "\n",
+ " else:\n",
+ "\n",
+ " continue\n",
+ "\n",
+ " file_download = requests.get(\n",
+ "\n",
+ " 'https://appeears.earthdatacloud.nasa.gov/api/bundle/{0}/{1}'.format(task_id, file_id),\n",
+ "\n",
+ " headers={'Authorization': 'Bearer {0}'.format(token)},\n",
+ "\n",
+ " allow_redirects=True,\n",
+ "\n",
+ " stream=True)\n",
+ "\n",
+ " file_download.raise_for_status()\n",
+ "\n",
+ " if not file_download.content:\n",
+ "\n",
+ " print(f\"Warning: Empty file downloaded for {file['file_name']}\")\n",
+ "\n",
+ " continue\n",
+ "\n",
+ " file_content = BytesIO(file_download.content)\n",
+ "\n",
+ " with rasterio.open(file_content) as src_initial:\n",
+ "\n",
+ " src = src_initial.read(1, masked=True)\n",
+ "\n",
+ " src_meta = src_initial.meta\n",
+ "\n",
+ " dst_crs = cdl_data.crs\n",
+ "\n",
+ " transform, width, height = calculate_default_transform(\n",
+ "\n",
+ " src_meta['crs'], dst_crs, src.shape[1], src.shape[0], *src_initial.bounds)\n",
+ "\n",
+ " kwargs = src_meta.copy()\n",
+ "\n",
+ " kwargs.update({\n",
+ "\n",
+ " 'crs': dst_crs,\n",
+ "\n",
+ " 'transform': transform,\n",
+ "\n",
+ " 'width': width,\n",
+ "\n",
+ " 'height': height\n",
+ "\n",
+ " })\n",
+ "\n",
+ " dst = np.zeros((height, width), dtype=src_meta['dtype'])\n",
+ "\n",
+ " reproject(\n",
+ "\n",
+ " source=src,\n",
+ "\n",
+ " destination=dst,\n",
+ "\n",
+ " src_transform=src_meta['transform'],\n",
+ "\n",
+ " src_crs=src_meta['crs'],\n",
+ "\n",
+ " dst_transform=transform,\n",
+ "\n",
+ " dst_crs=dst_crs,\n",
+ "\n",
+ " resampling=Resampling.nearest)\n",
+ "\n",
+ " kwargs.update({'nodata': 0})\n",
+ "\n",
+ " key = 'MODIS_NDVI'\n",
+ "\n",
+ " if key not in data:\n",
+ "\n",
+ " data[key] = {'data': [], 'meta': kwargs, 'doy': []}\n",
+ "\n",
+ " data[key]['data'].append(dst)\n",
+ "\n",
+ " data[key]['doy'].append(doy)\n",
+ "\n",
+ " else:\n",
+ "\n",
+ " continue"
+ ],
+ "id": "ea247046",
+ "execution_count": null,
+ "outputs": []
+ },
+ {
+ "cell_type": "code",
+ "metadata": {},
+ "source": [
+ "ndvi = np.stack(data['MODIS_NDVI']['data']) / 100\n",
+ "\n",
+ "day = np.stack(data['MODIS_NDVI']['doy'])\n",
+ "\n",
+ "meta = data['MODIS_NDVI']['meta']"
+ ],
+ "id": "e6152067",
+ "execution_count": null,
+ "outputs": []
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "## Combine and compare the results\n"
+ ],
+ "id": "1589255e"
+ },
+ {
+ "cell_type": "code",
+ "metadata": {},
+ "source": [
+ "ndvi_resampled, ndvi_resampled_meta = resample_using_cdl(ndvi, meta, cdl_meta_reprojected)\n",
+ "\n",
+ "slr_resampled, slr_resampled_meta = resample_using_cdl(slr_reprojected, slr_meta_reprojected, cdl_meta_reprojected)\n",
+ "\n",
+ "print(\"CDL shape:\", cdl_reprojected.shape)\n",
+ "\n",
+ "print(\"Original NDVI shape:\", ndvi.shape)\n",
+ "\n",
+ "print(\"Resampled NDVI shape:\", ndvi_resampled.shape)\n",
+ "\n",
+ "print(\"Original SLR shape:\", slr_reprojected.shape)\n",
+ "\n",
+ "print(\"Resampled SLR shape:\", slr_resampled.shape)"
+ ],
+ "id": "67f450f1",
+ "execution_count": null,
+ "outputs": []
+ },
+ {
+ "cell_type": "code",
+ "metadata": {},
+ "source": [
+ "def mask_rasters(cdl_raster, slr_raster, ndvi_rasters):\n",
+ "\n",
+ " corn_mask = cdl_raster == 1\n",
+ "\n",
+ " ABOVE_water_mask = slr_raster <= 1\n",
+ "\n",
+ " corn_above_water_mask = corn_mask & ~ABOVE_water_mask\n",
+ "\n",
+ " corn_below_water_mask = corn_mask & ABOVE_water_mask\n",
+ "\n",
+ " corn_above_water_mask = corn_above_water_mask[0,:,:]\n",
+ "\n",
+ " corn_below_water_mask = corn_below_water_mask[0,:,:]\n",
+ "\n",
+ " masked_ndvi_below_water = ndvi_rasters[:, corn_below_water_mask]\n",
+ "\n",
+ " masked_ndvi_above_water = ndvi_rasters[:, corn_above_water_mask]\n",
+ "\n",
+ " return {\n",
+ "\n",
+ " 'corn_below_water': {\n",
+ "\n",
+ " 'mask': corn_below_water_mask,\n",
+ "\n",
+ " 'ndvi': masked_ndvi_below_water,\n",
+ "\n",
+ " 'mean': np.mean(masked_ndvi_below_water, axis=1) if masked_ndvi_below_water.size > 0 else None,\n",
+ "\n",
+ " 'min': np.min(masked_ndvi_below_water, axis=1) if masked_ndvi_below_water.size > 0 else None,\n",
+ "\n",
+ " 'max': np.max(masked_ndvi_below_water, axis=1) if masked_ndvi_below_water.size > 0 else None\n",
+ "\n",
+ " },\n",
+ "\n",
+ " 'corn_above_water': {\n",
+ "\n",
+ " 'mask': corn_above_water_mask,\n",
+ "\n",
+ " 'ndvi': masked_ndvi_above_water,\n",
+ "\n",
+ " 'mean': np.mean(masked_ndvi_above_water, axis=1) if masked_ndvi_above_water.size > 0 else None,\n",
+ "\n",
+ " 'min': np.min(masked_ndvi_above_water, axis=1) if masked_ndvi_above_water.size > 0 else None,\n",
+ "\n",
+ " 'max': np.max(masked_ndvi_above_water, axis=1) if masked_ndvi_above_water.size > 0 else None\n",
+ "\n",
+ " }\n",
+ "\n",
+ " }"
+ ],
+ "id": "34d86791",
+ "execution_count": null,
+ "outputs": []
+ },
+ {
+ "cell_type": "code",
+ "metadata": {},
+ "source": [
+ "def visualize_ndvi_analysis(ndvi_analysis, days):\n",
+ "\n",
+ " plt.figure(figsize=(12, 6))\n",
+ "\n",
+ " plt.subplot(1, 2, 1)\n",
+ "\n",
+ " plt.plot(days, ndvi_analysis['corn_above_water']['mean'], label='Corn Above Water', color='green')\n",
+ "\n",
+ " plt.plot(days, ndvi_analysis['corn_below_water']['mean'], label='Corn Below Water', color='blue')\n",
+ "\n",
+ " plt.title('Mean NDVI')\n",
+ "\n",
+ " plt.xlabel('Days')\n",
+ "\n",
+ " plt.xticks(rotation=90)\n",
+ "\n",
+ " plt.ylabel('NDVI')\n",
+ "\n",
+ " plt.legend()\n",
+ "\n",
+ " plt.subplot(1, 2, 2)\n",
+ "\n",
+ " plt.plot(days, ndvi_analysis['corn_above_water']['max'], label='Corn Above Water', color='green')\n",
+ "\n",
+ " plt.plot(days, ndvi_analysis['corn_below_water']['max'], label='Corn Below Water', color='blue')\n",
+ "\n",
+ " plt.title('Maximum NDVI')\n",
+ "\n",
+ " plt.xlabel('Days')\n",
+ "\n",
+ " plt.xticks(rotation=90)\n",
+ "\n",
+ " plt.ylabel('NDVI')\n",
+ "\n",
+ " plt.legend()\n",
+ "\n",
+ " plt.tight_layout()\n",
+ "\n",
+ " plt.show()"
+ ],
+ "id": "9639f1af",
+ "execution_count": null,
+ "outputs": []
+ },
+ {
+ "cell_type": "code",
+ "metadata": {},
+ "source": [
+ "ndvi_analysis = mask_rasters(cdl_reprojected, slr_resampled, ndvi_resampled)\n",
+ "\n",
+ "visualize_ndvi_analysis(ndvi_analysis, day)"
+ ],
+ "id": "11e222ee",
+ "execution_count": null,
+ "outputs": []
+ },
+ {
+ "cell_type": "code",
+ "metadata": {},
+ "source": [
+ "for condition in ['corn_below_water', 'corn_above_water']:\n",
+ "\n",
+ " print(f\"\\n{condition.replace('_', ' ').title()} Corn Analysis per MODIS:\")\n",
+ "\n",
+ " for stat in ['mean', 'max']:\n",
+ "\n",
+ " print(f\"{stat.capitalize()} NDVI: {ndvi_analysis[condition][stat]}\")"
+ ],
+ "id": "a9886ec5",
+ "execution_count": null,
+ "outputs": []
+ },
+ {
+ "cell_type": "code",
+ "metadata": {},
+ "source": [
+ "def visualize_ndvi_histogram(ndvi_analysis):\n",
+ "\n",
+ " plt.figure(figsize=(12, 6))\n",
+ "\n",
+ " plt.hist(ndvi_analysis['corn_below_water']['ndvi'].flatten(), bins=30, alpha=0.5, label='Corn Below Water', color='blue')\n",
+ "\n",
+ " plt.hist(ndvi_analysis['corn_above_water']['ndvi'].flatten(), bins=30, alpha=0.5, label='Corn Above Water', color='green')\n",
+ "\n",
+ " plt.xlabel(\"NDVI values\")\n",
+ "\n",
+ " plt.ylabel(\"Frequency\")\n",
+ "\n",
+ " plt.title(\"Histogram of NDVI values for corn\")\n",
+ "\n",
+ " plt.legend(loc='upper right')\n",
+ "\n",
+ "visualize_ndvi_histogram(ndvi_analysis)"
+ ],
+ "id": "c1cae7d1",
+ "execution_count": null,
+ "outputs": []
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "What are the insights that we can derive from this module?\n",
+ "\n",
+ "What are the impacts of land being close to the water currently?\n",
+ "\n",
+ "What is the potential loss of productive land?\n",
+ "\n",
+ "Where could there be errors given the sampling methods?\n",
+ "\n",
+ "# What to do next\n",
+ "\n",
+ "1. Clip the MODIS data to the CDL data of only the agricultural values of interest\n",
+ "2. Find the average NDVI of these values\n",
+ "3. Clip this to the SLR area\n",
+ "4. Find the average NDVI of these values\n",
+ "5. Compare the changes for the year of 2018\n",
+ "\n",
+ "This will tell us the impacted area that could be effected by the floods\n",
+ "\n",
+ "As we look into the future we should see this change become more drastic\n",
+ "\n",
+ "# ABOUT NASA ACRES AND THE CLIMATE RESILIENCE NETWORK\n",
+ "\n",
+ "## NASA Acres\n",
+ "\n",
+ "NASA Acres is NASA's U.S.-focused agriculture and food security Consortium. NASA Acres is commissioned under NASA's Applied Sciences Program and led by the University of Maryland, with over 30 partner institutions and growing. The consortium approach brings together public and private sectors, and allows for flexible partnerships and rapid action in delivering the benefits of NASA data and tools to U.S. agriculture decision makers. We work with people across the agriculture sector to develop Earth observation (EO)-based data and tools that strengthen U.S. agriculture.\n",
+ "\n",
+ "## Climate Resilience Network\n",
+ "\n",
+ "Climate Resilience Network has the goal of unifying the University of Maryland's groundbreaking research with the immediate requirements of Maryland's decision-makers to address the complexities of climate change as it affects our home. GEOG's Grand Challenges project is aimed at promoting climate-resilient agriculture in Maryland, creating interactive tools designed to support research based decision-making.italicized text\n",
+ "\n",
+ "COUNTY LEVEL SHAPEFILE "
+ ],
+ "id": "a52e4a4d"
+ }
+ ],
+ "metadata": {
+ "kernelspec": {
+ "name": "python3",
+ "language": "python",
+ "display_name": "Python 3 (ipykernel)"
+ }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 5
+}
\ No newline at end of file
diff --git a/topst_final.py b/topst_final.py
new file mode 100644
index 0000000..0f9bdeb
--- /dev/null
+++ b/topst_final.py
@@ -0,0 +1,425 @@
+# -*- coding: utf-8 -*-
+"""TOPST_FINAL
+
+Automatically generated by Colab.
+
+Original file is located at
+ https://colab.research.google.com/drive/1JKqL4e1enqrS31UHJoKQ3Godj5q_aWJk
+
+CLIMATE SMART AGRICULTURE IN MARYLAND
+A regional view of the past with applications to the present. Python spatial and tabular analysis using multiple variables to look at a single time slice. This code can be a combination of cloud and local based analysis. A purpose for this exercise is to use use public research data for derived conclusions view the impacts climate change has on agriculture a look at what farmers can do and why we need to support them.
+
+NASA Acres has defined the Essential Agricultural Variables or EAVs. These were designed by NASA Acres to define the capabilities of the top satellite data scientists and practitioners who make up NASA Acres Research, Development, and Extension partners, and the needs of decision-making-collaborators already in our network, to identify an initial set of focus.
+1. Cropland and Crop Type Mapping
+2. Crop and Crop Type Area Estimation
+
+
+Derived results will be the county level statistics of impacts of harvestable acreage. To achieve this, a basic workflow we can
+1. Use api to call in the CDL dataset to map crop types.
+
+Then
+1. Access the Sea level rise (elevation dataset) to identify new areas of potential areas that are at risk of future flooding. Clip to county of interest
+2. Compare the CDL with the SLR mask and with out to identify the NOAA estimated loss of land.
+
+Then using MODIS create a time series for NDVI of the masked crop land, derive insights about trends in NDVI
+1. Use NDVI to create a time series looking backwards in time at areas that have experienced flooding to visualize the moving from productive farms to moderate quality.
+
+
+
+
+---
+
+
+This all together would allow us to make a predictive analysis for Maryland in the future under the projections of sea level rise.
+Given the current conditions, subtracting the sea level rise inundated areas.
+Then given the remaining areas, making predictions about the trends of NDVI given the trends in historical areas that are near impacted areas. We cannot say that
+sea level rise is driving the decrease because there are many other factor and decisions that farmers make around how much to plant, how much to harvest, chemicals applied,
+ growing degree days, soil and soil moisture conditions. But this module should provide the ability to make insights about potential causes and impacts.
+
+
+
+The data story that we have derived is about sea level rise in Maryland and the impact that it has on the production levels with in the state. From this module,
+ we can provide students with the ability to draw from multiple data sources, as well as derive insights using historical and future viewing data sets.
+
+
+
+The U.S. Mid-Atlantic has seen higher rates of sea level rise, marshes may be especially vulnerable. In the Chesapeake Bay, sea level rise has already
+ contributed to the degradation of over 80,000 ha (70%) of tidal marsh (Taylor, Curson, Verutes, et al). This view can help us to understand the impacts
+ that small levels of sea level rise have on land.
+
+
+
+
+We can prompt the user to think about future impacts outside the direct sea level rise projections, allowing them to include a more full picture,
+and finally using that picture to identify economic impacts that action or inaction causes. This begs the question, what can the public do to enact
+changes, rather than putting the pressure on farmers to change?
+
+---
+
+
+SOURCES
+
+Taylor, L., Curson, D., Verutes, G.M. et al. Mapping sea level rise impacts to identify climate change adaptation opportunities in the Chesapeake and Delaware Bays,
+ USA. Wetlands Ecol Manage 28, 527541 (2020). https://doi.org/10.1007/s11273-020-09729-w
+
+Mondal, P., Walter, M., Miller, J. et al. The spread and cost of saltwater intrusion in the US Mid-Atlantic. Nat Sustain 6, 13521362 (2023).
+ https://doi.org/10.1038/s41893-023-01186-6
+
+UMD + NOAA short web course
+https://www.mdsg.umd.edu/coastal-climate-resilience/farming-flooding-salt-land-loss#:~:text=Farmers%20and%20woodlot%20owners%20in%20Maryland%20and%20Virginia%20are%20facing,under%20saltier%20or%20wetter%20conditions.
+"""
+
+USERNAME = 'NASA USERNAME'
+PASSWORD = 'NASA PASSWORD'
+
+# !pip install rasterio
+
+import geopandas as gpd
+import rasterio
+import numpy as np
+import pandas as pd
+import matplotlib.pyplot as plt
+from shapely.geometry import box
+import numpy as np
+import requests
+import xml.etree.ElementTree as ET
+import matplotlib.pyplot as plt
+from io import BytesIO
+import os
+from rasterio.warp import calculate_default_transform, reproject, Resampling
+from datetime import datetime, timedelta
+
+def reproject_raster(src_dataset, src_crs, src_transform, dst_crs='EPSG:4326'):
+ dst_crs = rasterio.crs.CRS.from_string(dst_crs)
+ dst_transform, dst_width, dst_height = calculate_default_transform(
+ src_crs, dst_crs,
+ src_dataset.width, src_dataset.height,
+ *src_dataset.bounds)
+ dst_kwargs = src_dataset.meta.copy()
+ dst_kwargs.update({
+ 'crs': dst_crs,
+ 'transform': dst_transform,
+ 'width': dst_width,
+ 'height': dst_height})
+ dst_data = np.zeros((src_dataset.count, dst_height, dst_width), dtype=src_dataset.dtypes[0])
+ for i in range(1, src_dataset.count + 1):
+ reproject(
+ source=rasterio.band(src_dataset, i),
+ destination=dst_data[i-1],
+ src_transform=src_transform,
+ src_crs=src_crs,
+ dst_transform=dst_transform,
+ dst_crs=dst_crs,
+ resampling=Resampling.nearest)
+ return dst_data, dst_kwargs
+
+def format_cdl_url(fips,year):
+ base_url = "https://nassgeodata.gmu.edu/axis2/services/CDLService/GetCDLFile"
+ url = f"{base_url}?year={year}&fips={fips}"
+ return url
+
+def get_srl_raster():
+ slr_path = '/content/MD_East_slr_depth_3_5ft.tif'
+ with rasterio.open(slr_path) as slr_raster:
+ slr_meta = slr_raster.meta
+ slr_meta['nodata'] = 0
+ slr_reprojected, slr_meta_reprojected = reproject_raster(src_dataset=slr_raster, src_crs=slr_meta['crs'], src_transform=slr_meta['transform'])
+ return slr_reprojected, slr_meta_reprojected
+
+"""
+
+# Gather sea level rise (SLR) raster, crop land database (CDL) raster, and County geodatabase
+
+To understand the nature of Maryland agriculture we can begin by using the CDL. This raster dataset comes at a 30 meter spatial resolution. We can access this data through their api. The current configuration allows data to be pulled at the county level.
+
+COUNTY LEVEL CDL FOR FIPS
+https://nassgeodata.gmu.edu/axis2/services/CDLService/GetCDLFile
+
+This view provides insight into the common commodities that are grown in Maryland in the desired year.
+"""
+
+crop_mapping = {0: "Background", 1: "Corn", 2: "Cotton", 3: "Rice", 4: "Sorghum", 5: "Soybeans", 6: "Sunflower", 10: "Peanuts",
+11: "Tobacco", 12: "Sweet Corn", 13: "Pop or Orn Corn", 14: "Mint", 21: "Barley", 22: "Durum Wheat",23: "Spring Wheat", 24: "Winter Wheat", 25: "Other Small Grains", 26: "Dbl Crop WinWht/Soybeans", 27: "Rye",
+28: "Oats", 29: "Millet", 30: "Speltz", 31: "Canola", 32: "Flaxseed", 33: "Safflower", 34: "Rape Seed",35: "Mustard", 36: "Alfalfa", 37: "Other Hay/Non Alfalfa", 38: "Camelina", 39: "Buckwheat", 41: "Sugarbeets",
+42: "Dry Beans", 43: "Potatoes", 44: "Other Crops", 45: "Sugarcane", 46: "Sweet Potatoes", 47: "Misc Vegs & Fruits",48: "Watermelons", 49: "Onions", 50: "Cucumbers", 51: "Chick Peas", 52: "Lentils", 53: "Peas", 54: "Tomatoes",
+55: "Caneberries", 56: "Hops", 57: "Herbs", 58: "Clover/Wildflowers", 59: "Sod/Grass Seed", 60: "Switchgrass",61: "Fallow/Idle Cropland", 63: "Forest", 64: "Shrubland", 65: "Barren", 66: "Cherries", 67: "Peaches",
+68: "Apples", 69: "Grapes", 70: "Christmas Trees", 71: "Other Tree Crops", 72: "Citrus", 74: "Pecans",75: "Almonds", 76: "Walnuts", 77: "Pears", 81: "Clouds/No Data", 82: "Developed", 83: "Water", 87: "Wetlands",
+88: "Nonag/Undefined", 92: "Aquaculture", 111: "Open Water", 112: "Perennial Ice/Snow", 121: "Developed/Open Space",122: "Developed/Low Intensity", 123: "Developed/Med Intensity", 124: "Developed/High Intensity", 131: "Barren",
+141: "Deciduous Forest", 142: "Evergreen Forest", 143: "Mixed Forest", 152: "Shrubland", 176: "Grass/Pasture",190: "Woody Wetlands", 195: "Herbaceous Wetlands", 204: "Pistachios", 205: "Triticale", 206: "Carrots",
+207: "Asparagus", 208: "Garlic", 209: "Cantaloupes", 210: "Prunes", 211: "Olives", 212: "Oranges",213: "Honeydew Melons", 214: "Broccoli", 215: "Avocados", 216: "Peppers", 217: "Pomegranates", 218: "Nectarines",
+219: "Greens", 220: "Plums", 221: "Strawberries", 222: "Squash", 223: "Apricots", 224: "Vetch",225: "Dbl Crop WinWht/Corn", 226: "Dbl Crop Oats/Corn", 227: "Lettuce", 228: "Dbl Crop Triticale/Corn",
+229: "Pumpkins", 230: "Dbl Crop Lettuce/Durum Wht", 231: "Dbl Crop Lettuce/Cantaloupe", 232: "Dbl Crop Lettuce/Cotton",233: "Dbl Crop Lettuce/Barley", 234: "Dbl Crop Durum Wht/Sorghum", 235: "Dbl Crop Barley/Sorghum",
+236: "Dbl Crop WinWht/Sorghum", 237: "Dbl Crop Barley/Corn", 238: "Dbl Crop WinWht/Cotton",239: "Dbl Crop Soybeans/Cotton", 240: "Dbl Crop Soybeans/Oats", 241: "Dbl Crop Corn/Soybeans", 242: "Blueberries",
+243: "Cabbage", 244: "Cauliflower", 245: "Celery", 246: "Radishes", 247: "Turnips", 248: "Eggplants",249: "Gourds", 250: "Cranberries", 254: "Dbl Crop Barley/Soybeans", 255: ""}
+
+county = 'Talbot'
+code = "24041"
+year = 2018
+
+
+service_url = format_cdl_url(code, year)
+response = requests.get(service_url)
+root = ET.fromstring(response.content)
+tif_url = root.find('.//returnURL').text
+cdl_data = rasterio.open(tif_url)
+cdl_meta = cdl_data.meta
+cdl_reprojected, cdl_meta_reprojected = reproject_raster(src_dataset=cdl_data, src_crs=cdl_meta['crs'], src_transform=cdl_meta['transform'])
+
+"""The eastern seaboard has seen changes in the sea level along with increase flooding and the decrease in predictability of water flow from rivers. This daa comes from NOAA (learn more here https://coast.noaa.gov/data/digitalcoast/pdf/slr-inundation-methods.pdf)
+
+This raster is the projection of inudated areas by the year of 2050. Data access can be found here https://coast.noaa.gov/slrdata/Depth_Rasters/MD/index.html
+"""
+
+slr_reprojected, slr_meta_reprojected = get_srl_raster()
+
+"""# Gather NDVI values"""
+
+def get_appeears_token(username, password):
+ response = requests.post('https://appeears.earthdatacloud.nasa.gov/api/login',auth=(username, password))
+ return response.json()['token']
+
+def submit_appeears_task(token, task):
+ response = requests.post('https://appeears.earthdatacloud.nasa.gov/api/task',json=task,headers={'Authorization': f'Bearer {token}'})
+ return response.json()
+
+def get_appeears_bundle(token, task_id):
+ response = requests.get(f'https://appeears.earthdatacloud.nasa.gov/api/bundle/{task_id}',headers={'Authorization': f'Bearer {token}'})
+ return response.json()
+
+def load_county_data():
+ sea_level_path = '/content/Talbot_County.geojson'
+ sea_level_gdf = gpd.read_file(sea_level_path)
+ return sea_level_gdf
+
+def day_of_year_to_date2(day_of_year_str):
+ year = int(day_of_year_str[:4])
+ day_of_year = int(day_of_year_str[4:])
+ start_date = datetime(year, 1, 1)
+ target_date = start_date + timedelta(days=day_of_year - 1)
+ return target_date.strftime('%m-%d-%Y')
+
+def resample_using_cdl(data, data_meta, cdl_meta_reprojected):
+ resampled_data = np.zeros(
+ (data.shape[0], cdl_meta_reprojected['height'], cdl_meta_reprojected['width']),
+ dtype=data.dtype
+ )
+
+ for i in range(data.shape[0]):
+ reproject(
+ source=data[i],
+ destination=resampled_data[i],
+ src_transform=data_meta['transform'],
+ src_crs=data_meta['crs'],
+ dst_transform=cdl_meta_reprojected['transform'],
+ dst_crs=cdl_meta_reprojected['crs'],
+ resampling=Resampling.bilinear
+ )
+ updated_meta = data_meta.copy()
+ updated_meta.update({
+ 'transform': cdl_meta_reprojected['transform'],
+ 'width': cdl_meta_reprojected['width'],
+ 'height': cdl_meta_reprojected['height'],
+ 'crs': cdl_meta_reprojected['crs']
+ })
+
+ return resampled_data, updated_meta
+
+county_data = load_county_data()
+county_reprojected = county_data.to_crs('EPSG:4326')
+us = county_reprojected.geometry.iloc[0]
+minx, miny, maxx, maxy = us.bounds
+
+
+bbox = [[maxx, miny],
+ [maxx, maxy],
+ [minx, maxy],
+ [minx, maxy],
+ [minx, miny]]
+
+
+start = f'06-01-{year}'
+end = f'09-30-{year}'
+
+
+products = [{'layer': '_250m_16_days_NDVI', 'product':'MOD13Q1.061'}]
+
+token = get_appeears_token(f'{USERNAME}', f'{PASSWORD}')
+
+task = {'task_type': 'area',
+ 'task_name': 'MODIS',
+ 'params': {"geo": {"type": "FeatureCollection",
+ "features": [{"type": "Feature",
+ "geometry": {"type": "MultiPolygon",
+ "coordinates": [[bbox]]},
+ "properties": {}}],
+ "fileName": "Polygon"},
+ 'dates': [{'startDate': start, 'endDate': end}],
+ 'layers': products,
+ 'coordinates': bbox,
+ "output": {"format": {"type": "geotiff"},
+ "projection": "native"}}}
+task_id = submit_appeears_task(token, task)
+
+bundle = get_appeears_bundle(token, task_id['task_id'])
+
+data = {}
+
+for file in bundle['files']:
+ file_id = file['file_id']
+ if 'NDVI' in file['file_name']:
+ if '_doy' in file['file_name']:
+ datesy = file['file_name'].split('_doy')[1][:7]
+ doy = day_of_year_to_date2(datesy)
+
+ print(file['file_name'])
+ else:
+ continue
+ file_download = requests.get(
+ 'https://appeears.earthdatacloud.nasa.gov/api/bundle/{0}/{1}'.format(task_id['task_id'], file_id),
+ headers={'Authorization': 'Bearer {0}'.format(token)},
+ allow_redirects=True,
+ stream=True)
+ file_download.raise_for_status()
+ if not file_download.content:
+ print(f"Warning: Empty file downloaded for {file['file_name']}")
+ continue
+ file_content = BytesIO(file_download.content)
+ with rasterio.open(file_content) as src_initial:
+ src = src_initial.read(1, masked=True)
+ src_meta = src_initial.meta
+ dst_crs = cdl_data.crs
+ transform, width, height = calculate_default_transform(
+ src_meta['crs'], dst_crs, src.shape[1], src.shape[0], *src_initial.bounds)
+ kwargs = src_meta.copy()
+ kwargs.update({
+ 'crs': dst_crs,
+ 'transform': transform,
+ 'width': width,
+ 'height': height
+ })
+
+ dst = np.zeros((height, width), dtype=src_meta['dtype'])
+ reproject(
+ source=src,
+ destination=dst,
+ src_transform=src_meta['transform'],
+ src_crs=src_meta['crs'],
+ dst_transform=transform,
+ dst_crs=dst_crs,
+ resampling=Resampling.nearest)
+
+
+ kwargs.update({'nodata': 0})
+ key = 'MODIS_NDVI'
+ if key not in data:
+ data[key] = {'data': [], 'meta': kwargs, 'doy': []}
+ data[key]['data'].append(dst)
+ data[key]['doy'].append(doy)
+ else:
+ continue
+
+ndvi = np.stack(data['MODIS_NDVI']['data']) / 100
+day = np.stack(data['MODIS_NDVI']['doy'])
+meta = data['MODIS_NDVI']['meta']
+
+"""# Combine and compare the results
+
+"""
+
+ndvi_resampled, ndvi_resampled_meta = resample_using_cdl(ndvi, meta, cdl_meta_reprojected)
+slr_resampled, slr_resampled_meta = resample_using_cdl(slr_reprojected, slr_meta_reprojected, cdl_meta_reprojected)
+
+
+print("CDL shape:", cdl_reprojected.shape)
+
+print("Original NDVI shape:", ndvi.shape)
+print("Resampled NDVI shape:", ndvi_resampled.shape)
+
+print("Original SLR shape:", slr_reprojected.shape)
+print("Resampled SLR shape:", slr_resampled.shape)
+
+def mask_rasters(cdl_raster, slr_raster, ndvi_rasters):
+ corn_mask = cdl_raster == 1
+ ABOVE_water_mask = slr_raster <= 1
+ corn_above_water_mask = corn_mask & ~ABOVE_water_mask
+ corn_below_water_mask = corn_mask & ABOVE_water_mask
+
+ corn_above_water_mask = corn_above_water_mask[0,:,:]
+ corn_below_water_mask = corn_below_water_mask[0,:,:]
+ masked_ndvi_below_water = ndvi_rasters[:, corn_below_water_mask]
+ masked_ndvi_above_water = ndvi_rasters[:, corn_above_water_mask]
+
+
+ return {
+ 'corn_below_water': {
+ 'mask': corn_below_water_mask,
+ 'ndvi': masked_ndvi_below_water,
+ 'mean': np.mean(masked_ndvi_below_water, axis=1) if masked_ndvi_below_water.size > 0 else None,
+ 'min': np.min(masked_ndvi_below_water, axis=1) if masked_ndvi_below_water.size > 0 else None,
+ 'max': np.max(masked_ndvi_below_water, axis=1) if masked_ndvi_below_water.size > 0 else None
+ },
+ 'corn_above_water': {
+ 'mask': corn_above_water_mask,
+ 'ndvi': masked_ndvi_above_water,
+ 'mean': np.mean(masked_ndvi_above_water, axis=1) if masked_ndvi_above_water.size > 0 else None,
+ 'min': np.min(masked_ndvi_above_water, axis=1) if masked_ndvi_above_water.size > 0 else None,
+ 'max': np.max(masked_ndvi_above_water, axis=1) if masked_ndvi_above_water.size > 0 else None
+ }
+ }
+
+def visualize_ndvi_analysis(ndvi_analysis, days):
+ plt.figure(figsize=(12, 6))
+
+ plt.subplot(1, 2, 1)
+ plt.plot(days, ndvi_analysis['corn_above_water']['mean'], label='Corn Above Water', color='green')
+ plt.plot(days, ndvi_analysis['corn_below_water']['mean'], label='Corn Below Water', color='blue')
+ plt.title('Mean NDVI')
+ plt.xlabel('Days')
+ plt.xticks(rotation=90)
+ plt.ylabel('NDVI')
+ plt.legend()
+
+
+ plt.subplot(1, 2, 2)
+ plt.plot(days, ndvi_analysis['corn_above_water']['max'], label='Corn Above Water', color='green')
+ plt.plot(days, ndvi_analysis['corn_below_water']['max'], label='Corn Below Water', color='blue')
+ plt.title('Maximum NDVI')
+ plt.xlabel('Days')
+ plt.xticks(rotation=90)
+ plt.ylabel('NDVI')
+ plt.legend()
+
+ plt.tight_layout()
+ plt.show()
+
+ndvi_analysis = mask_rasters(cdl_reprojected, slr_resampled, ndvi_resampled)
+
+visualize_ndvi_analysis(ndvi_analysis, day)
+
+for condition in ['corn_below_water', 'corn_above_water']:
+ print(f"\n{condition.replace('_', ' ').title()} Corn Analysis per MODIS:")
+ for stat in ['mean', 'max']:
+ print(f"{stat.capitalize()} NDVI: {ndvi_analysis[condition][stat]}")
+
+def visualize_ndvi_histogram(ndvi_analysis):
+ plt.figure(figsize=(12, 6))
+
+ plt.hist(ndvi_analysis['corn_below_water']['ndvi'].flatten(), bins=30, alpha=0.5, label='Corn Below Water', color='blue')
+ plt.hist(ndvi_analysis['corn_above_water']['ndvi'].flatten(), bins=30, alpha=0.5, label='Corn Above Water', color='green')
+ plt.xlabel("NDVI values")
+ plt.ylabel("Frequency")
+ plt.title("Histogram of NDVI values for corn")
+ plt.legend(loc='upper right')
+
+
+visualize_ndvi_histogram(ndvi_analysis)
+
+"""What are the insights that we can derive from this module?
+
+What are the impacts of land being close to the water currently?
+
+What is the potential loss of productive land?
+
+Where could there be errors given the sampling methods?
+
+"""