# Load your spatial data (update paths accordingly)

flood_zones <- st_read("NYC-Flood_Vulnerability/NYC Flood Vulnerability_SS.shp")                                              
## Reading layer `NYC Flood Vulnerability_SS' from data source 
##   `C:\Users\default.LAPTOP-BEI6IEOC\Documents\Data Analysis and Viz\Term Project\Project\NYC-Flood_Vulnerability\NYC Flood Vulnerability_SS.shp' 
##   using driver `ESRI Shapefile'
## Simple feature collection with 311 features and 8 fields
## Geometry type: MULTIPOLYGON
## Dimension:     XY
## Bounding box:  xmin: -74.16154 ymin: 40.5276 xmax: -73.72479 ymax: 40.88969
## Geodetic CRS:  WGS84(DD)
schools_all <- st_read("NYC MIddle Shool/NYC Middle Schools.shp")
## Reading layer `NYC Middle Schools' from data source 
##   `C:\Users\default.LAPTOP-BEI6IEOC\Documents\Data Analysis and Viz\Term Project\Project\NYC MIddle Shool\NYC Middle Schools.shp' 
##   using driver `ESRI Shapefile'
## Simple feature collection with 599 features and 18 fields
## Geometry type: POINT
## Dimension:     XY
## Bounding box:  xmin: -8264724 ymin: 4940054 xmax: -8205724 ymax: 4995854
## Projected CRS: WGS 84 / Pseudo-Mercator
# Create two separate layers
schools_in_flood <- st_read("Schools_In_Flood_Zone/Schools_In_Flood_Zone.shp")
## Reading layer `Schools_In_Flood_Zone' from data source 
##   `C:\Users\default.LAPTOP-BEI6IEOC\Documents\Data Analysis and Viz\Term Project\Project\Schools_In_Flood_Zone\Schools_In_Flood_Zone.shp' 
##   using driver `ESRI Shapefile'
## Simple feature collection with 94 features and 17 fields
## Geometry type: POINT
## Dimension:     XY
## Bounding box:  xmin: -8248840 ymin: 4950169 xmax: -8210832 ymax: 4993941
## Projected CRS: WGS 84 / Pseudo-Mercator
schools_not_in_flood <- st_read("Schools_Not_In_Flood_Zone/Schools_Not_In_Flood_Zone.shp")
## Reading layer `Schools_Not_In_Flood_Zone' from data source 
##   `C:\Users\default.LAPTOP-BEI6IEOC\Documents\Data Analysis and Viz\Term Project\Project\Schools_Not_In_Flood_Zone\Schools_Not_In_Flood_Zone.shp' 
##   using driver `ESRI Shapefile'
## Simple feature collection with 599 features and 17 fields
## Geometry type: POINT
## Dimension:     XY
## Bounding box:  xmin: -8264724 ymin: 4940054 xmax: -8205724 ymax: 4995854
## Projected CRS: WGS 84 / Pseudo-Mercator
library(leaflet)
library(sf)

# Reproject all layers to WGS84 (EPSG:4326)
flood_zones <- st_transform(flood_zones, crs = 4326)
schools_in_flood <- st_transform(schools_in_flood, crs = 4326)
schools_not_in_flood <- st_transform(schools_not_in_flood, crs = 4326)

# Build the map
leaflet() %>%
  addProviderTiles("CartoDB.Positron") %>%
  
  addPolygons(data = flood_zones,
              color = "blue",
              weight = 1,
              fillOpacity = 0.4,
              group = "Flood Zones") %>%
  
  addCircleMarkers(data = schools_in_flood,
                   color = "red",
                   radius = 4,
                   stroke = FALSE,
                   fillOpacity = 0.8,
                   label = ~ATS_CODE,  # Change this to a valid column name
                   group = "Schools In Flood Zone") %>%
  
  addCircleMarkers(data = schools_not_in_flood,
                   color = "orange",
                   radius = 4,
                   stroke = FALSE,
                   fillOpacity = 0.8,
                   label = ~ATS_CODE,  # Same here
                   group = "Schools Not In Flood Zone") %>%
  
  addLayersControl(
    overlayGroups = c("Flood Zones", "Schools In Flood Zone", "Schools Not In Flood Zone"),
    options = layersControlOptions(collapsed = FALSE)
  )