GHG emissions/update accounting within the leveed polygons (rasters)

camilavargasp commented 4 months ago

Email From Karrin Alstad

I am responsible for the GHG emissions/update accounting within the leveed polygons of our delta study area. I am following the directions provided by Steven Deverel from HydroFocus for integrating the larger scale NOAA land cover data (raster data) with more local land cover data sets: BAARI, DARI, and CDFW VegCAMP mapping. I have attached the instructions Steven sent to me and he is someone that can be contacted for further details. I was going to try to recreate the analyses he presents first to make sure I get similar results as he got (included I the document) before I changed the boundary to that used in our specific analyses.

The first difficulty I have run into is just working with raster data. I am getting stronger in my use of R sf package, but I have not yet extracted information from raster data types.

I have downloaded the the CCAP canopy and water .tif files from this site: https://coastalimagery.blob.core.windows.net/ccap-landcover/CCAP_bulk_download/High_Resolution_Land_Cover/Phase_1_Initial_Layers/index.html

From step 4 of the instructions, it seems like I am trying to extract 4 types of water areas from the CCAP data.

I have moved all of the datasets into our project sharepoint folder at:

Documents/NCEAS-DSP working group/2023/Infrastructure Group/Data Sets and Products/GHG

Another question comes up when I look into the CARI and the BAARI datasets which have a file type of .gdb. Apparently this is an ESRII file format and I do not know how to convert this to an R sf file.

Here is a preview of the CARI dataset: https://www.ecoatlas.org/regions/ecoregion/bay-delta

Anyway – lots of challenges in this exercise. Hopefully you are still interested to help me consider a good approach after learning more about it.

The Delta 2002 and Suisun 2003 VegCAMP mapping data have been downloaded in another folder for our project so I will find these and move a copy over to this GHG folder for convenience.

camilavargasp commented 4 months ago

Find the attached file here: scoping_plan_writeup_04252022

angelchen7 commented 4 months ago

Email from Karrin on 2/9:

Hello Angel,

Thank you for your email – and I had not heard about the terra R package. I will look into this right away.

To be honest, no, I have not yet coded steps 1-3,, and so I do not have code to provide yet – except some practice attempts with the raster package. I will try to make some headway this weekend. I will also find and move the last dataset (cdfw) over to the GHG folder so we have a complete set.

Thank you for your thoughts on how to proceed.

Have a good weekend!

Kind regards,

Karrin

angelchen7 commented 4 months ago

Created `terra` demonstration script and invited Karrin to a meeting

I created a script here to show Karrin some of the terra functions that I think will be useful for this upcoming project. Then I invited her to a meeting next week where we can brainstorm some ways for me to help her get started.

angelchen7 commented 3 months ago

Met with Karrin

I forgot to update here but I met with Karrin on Feb 23 and here are some notes from our meeting:

Karrin been getting used to sf
Need to confirm tidal boundary w group (step 2)
Karrin: good to keep going for now, will keep working on the clipping step (step 3)
Ping Chris Kwan about Sharepoint permission issue so we can possibly try to get the Microsoft 365 R package working

I think it's fine to leave Karrin alone for now since she wants to work on it by herself a bit more but I let her know that she can come to me for any issues at any time.

angelchen7 commented 3 months ago

Since it's been a while since I've heard from Karrin, I think it's fine to close the issue for now. I can reopen or open a new issue if she gets back to me.

angelchen7 commented 3 months ago

Reopening issue

I'm reopening this issue because Karrin reached out to me for additional help on working with rasters on March 20 and we've been regularly corresponding via email since.

Been working on the GHGaccounting3.R and GHGaccounting4.R scripts

angelchen7 commented 2 months ago

Summary of progress so far

March 20:

Karrin asking about reprojecting raster data and cropping/masking the raster with a vector boundary
My solution in this script

March 28:

Karrin asking about how to access the raster values to see the different land categories
My solution in this script

March 29:

Karrin asking about how to do step 3 in the instructions (see below)
Then I solved the first part of her problem about getting the midpoint in this script

Hi Angel,

I wonder if I may get your thoughts on step 3 in my instructions:

“We delineated the brackish from freshwater areas as the midpoint between Browns Island and Sherman Lake. We then created two corresponding polygons and clipped all spatial datasets within the polygons.”

I have imported the AllIslands dataset and have calculated the distance between Brown Island and Sherman Lake, but I don’t know how to:

# add the midpoint to the CCAP data / map

# create two corresponding polygons

# clip all spatial datasets within the polygons

Thank you for any suggestions!

# Angel's solution
# get just the Sherman Island
sherman <- islands %>%
  filter(NAME=="SHERMAN ISLAND")

# get just the Browns Island
browns <- islands %>%
  filter(NAME=="BROWNS ISLAND")

# draw a line to represent the least amount of distance between the two islands
line_betw_islands <- sf::st_nearest_points(sherman, browns)

# find the midpoint on that line
midpoint <- sf::st_centroid(line_betw_islands)

# check the plot
ggplot() +
  geom_sf(data = sherman) +
  geom_sf(data = browns) +
  geom_sf(data = line_betw_islands) + 
  geom_sf(data = midpoint)

# reproject the midpoint to use WGS 84 to get the exact lat/lon
midpoint_lat_lon <- sf::st_transform(midpoint, crs = 4326)
midpoint_lat_lon

angelchen7 commented 2 months ago

April 4 update

I touched base with Karrin to see how she was doing in her progress. Then she reached out for assistance (see below).

Hi Angel,

Thanks so much for reaching out!   I have just uploaded GHGaccounting5.R to the Github site.

At the bottom of this program, I am included an example code that uses a point to find the line to the edges of the polygon on the same longitude.

I wonder if there is an easier way to do this, and then I need to create 2 polygons – which I have not done yet.

And then, I need to apply this example approach to my tidal boundary polygon.

Thanks for having a look at this approach if you have time – I wonder if there is a more efficient way to do this?

I appreciate your help~!

After asking for more details on the problem, I came up with this solution:

I've created the separate polygons based on the method in the StackOverflow post that you linked the other day. I'm not sure if you've done so already but here is my own attempt:

# Angel's attempt
# split and create separate polygons from the vertical line
separate_polygons <- lwgeom::st_split(tidal, vertical_line_sf) %>%
  sf::st_collection_extract("POLYGON")

# combine the two smaller polygons together to represent the brackish area
brackish <- sf::st_union(sf::st_geometry(separate_polygons[3,]), sf::st_geometry(separate_polygons[2,]))

# get the remaining freshwater area
freshwater <- separate_polygons[1,]

# check the plot
ggplot() +
  geom_sf(data = brackish, aes(fill = "blue"))+
  geom_sf(data = freshwater, aes(fill = "red"))+
  scale_fill_identity(guide = "legend", labels = c("brackish","freshwater"))

In the end, I got the 2 sf polygons representing the brackish and freshwater areas. Feel free to try it out if it's useful to you!

angelchen7 commented 1 month ago

May 1 Update

Karrin has been reaching out to Steven Deverel in order to ask him for more information on how to replicate his original spatial workflow. Karrin is still working with him to get some of the shapefiles that he used for the Landscape Scenario Planning Tool. She is currently trying to replicate the Appendix A methods in order to provide GHG valuation estimates for the nceas levee failure [RAND redo] model.

I'm unsure when I need to jump in but Karrin already knows that she can reach out to me at any point.

camilavargasp commented 1 month ago

May 7 Update

Email from Karrin

Hi Camila,

I appreciate your offer to consider my analyses.  I will try to be very specific in my questions.  (outlined in GHGaccounting8.R in the R folder of the intrastructure github site.

The first thing I am running into is a problem trying to clip down a large raster file that contains land use types.

I am using the tidal boundary as the extent of the clip. 

Because the raster is so large ,  it was proving difficult to reproject to match the tidal boundary.

Angel's idea was to reproject the vector to the same projection as the CCAP raster data because is was impossible to reproject CCAP (too large)

She advised to first:

# convert tidal into a terra vector object

tidal_tv <- terra::vect(tidal)

However: when i run this re-projection (a  couple of different ways) I get errors:

project(tidal_tv, c2016r, partial = FALSE) # There were 49 warnings (use warnings() to see them)

# reproject tidal_tv with the same CRS as the canopy raster

tidal_tv_reproj <- terra::project(tidal_tv, "EPSG:5070", partial = FALSE) 

Errors:

# There were 49 warnings (use warnings() to see them)

# 1: PROJ: Cannot open https://cdn.proj.org/ca_nrc_NA83SCRS.tif: schannel: CertGetCertificateChain trust error CERT_TRUST_IS_UNTRUSTED_ROOT (GDAL error 1)

I am not sure how to correct these errors in this very first step of my analyses…  :/

mbjones commented 1 month ago

This is a web server / web client configuration error regarding how the web server is set up for SSL:

CertGetCertificateChain trust error CERT_TRUST_IS_UNTRUSTED_ROOT

I just tried the URL from my mac and from included-crab, and it seems fine, so it is likely an issue with the client config. Do you know where this code was run? On whatever host it is on, does the following command work?

curl -s -I https://cdn.proj.org/ca_nrc_NA83SCRS.tif

camilavargasp commented 1 month ago

Thanks @mbjones

It works for me too! I believe Karrin is running her analysis on her work computer. So, there might be an issue there. I'll let Karrin know she can try running her analysis in included-crab.

angelchen7 commented 4 weeks ago

May 30 Update

Back from my vacation and Karrin had a question for me on why her code wasn't working in GHGaccounting10.R. I emailed back saying that her vector data wasn't intersecting with her raster data, and that was why her raster could not be clipped.

NCEAS / learning-hub-organization