diff --git a/cities.csv b/cities.csv
index 778c300..776e14a 100644
--- a/cities.csv
+++ b/cities.csv
@@ -5,3 +5,4 @@ Portland,city,OR,United States
 Port Angeles,city,WA,United States
 Boston,city,MA,United States
 Boise City,city,ID,United States
+Salt Lake City,city,UT,United States
diff --git a/city_compare.Rmd b/city_compare.Rmd
index e6b0d81..fde1df1 100644
--- a/city_compare.Rmd
+++ b/city_compare.Rmd
@@ -20,6 +20,9 @@ date()
 rm(list=ls())
 library(tidyverse)
 library(tidycensus)
+library(sf)
+library(openmeteo)
+library(maps)
 
 ```
 
@@ -30,12 +33,19 @@ census_api_key(key = substr(read_file(file = "api_keys/census_api_key"), 1, 40))
 
 ```
 
+## Date ranges
+```{r date_range, eval = TRUE, echo = TRUE, results = "show", warning = FALSE, error = TRUE, message = FALSE}
+date_start <- "2010-01-01"
+date_end <- "2024-12-31"
+
+```
 ## Cities to compare
 
 ```{r cities, eval = TRUE, echo = TRUE, results = "show", warning = FALSE, error = TRUE, message = FALSE}
 cities <- read_csv(file = "cities.csv")
 cities <- cities %>%
   mutate(city_name = paste0(City, " ", Type))
+
 ```
 
 # Get data
@@ -49,14 +59,25 @@ for(city in cities$city_name){
     geography = "place",
     variables = "B01003_001",
     state = state,
-    year = 2023
+    year = 2023,
+    geometry = TRUE
   ) %>%
   filter(str_detect(NAME, city))
 }
 
 populations <- bind_rows(populations)
 
-cities <- bind_cols(cities, populations)
+city_center <- populations %>%
+  st_centroid() %>%
+  st_transform(4326) %>%  # Convert to WGS84 (standard lat/lon)
+  mutate(
+    lon = st_coordinates(.)[,1],
+    lat = st_coordinates(.)[,2]
+  ) %>%
+  st_drop_geometry() %>%
+  select(lat, lon)
+
+cities <- bind_cols(cities, populations, city_center)
 
 ggplot(cities) +
   geom_col(aes(x = City,
@@ -64,26 +85,65 @@ ggplot(cities) +
 
 ```
 
-## Weather
+## Map cities
+```{r cities_map, eval = TRUE, echo = TRUE, results = "show", warning = FALSE, error = TRUE, message = FALSE}
+ggplot(data = cities) +
+  geom_polygon(data = map_data(map = "state"),
+               aes(long, lat, group = group),
+               fill = "white", colour = "grey50") +
+  geom_point(aes(x = lon,
+                 y = lat),
+             shape = 21,
+             fill = "lightgreen",
+             color = "black",
+             size = 4)
+
+
+```
+
+## weather
 ```{r weather, eval = TRUE, echo = TRUE, results = "show", warning = FALSE, error = TRUE, message = FALSE}
-populations <- list(NULL)
-for(city in cities$city_name){
-  state <- cities %>% filter(city_name == city) %>% pull(State)
-  populations[[city]] <- get_acs(
-    geography = "place",
-    variables = "B01003_001",
-    state = state,
-    year = 2023
-  ) %>%
-  filter(str_detect(NAME, city))
+weather <- list(NULL)
+for(city in cities$City){
+  city_info <- cities %>% filter(City == city)
+  city_run <- weather_history(
+    location = c(city_info %>% pull(lat), city_info %>% pull(lon)),
+    start = date_start,
+    end = date_end,
+    daily = c("apparent_temperature_max", "apparent_temperature_min"),
+    response_units = list(temperature_unit = "fahrenheit")
+  )
+  city_run$city <- city
+  weather[[city]] <- city_run
 }
 
-populations <- bind_rows(populations)
+weather <- bind_rows(weather)
 
-cities <- bind_cols(cities, populations)
+weather_summary <- weather %>%
+  mutate(year = year(ymd(date)),
+         month = month(ymd(date))) %>%
+  group_by(year, city) %>%
+  summarise(days_above_80 = sum(daily_apparent_temperature_max > 80)) %>%
+  group_by(city) %>%
+  summarise(median_days_above_80 = median(days_above_80))
 
-ggplot(cities) +
-  geom_col(aes(x = City,
-               y = estimate))
+ggplot(data = weather_summary) +
+  geom_col(aes(x = city,
+               y = median_days_above_80)) +
+  labs(title = "Days above 80°F",
+       y = "Median days per year",
+       x = "City",
+       fill = NULL)
+
+ggplot(data = weather %>% pivot_longer(cols = starts_with("daily"), names_to = "max_min", values_to = "temperature") %>% filter(max_min %in% c("daily_apparent_temperature_min", "daily_apparent_temperature_max"))) +
+  geom_violin(aes(x = city,
+                  y = temperature,
+                  fill = max_min)) +
+  scale_fill_manual(labels = c("daily max", "daily min"), 
+                    values = c("firebrick", "dodgerblue")) +
+  labs(title = "Apparent Temperature",
+       y = "°F",
+       x = "City",
+       fill = NULL)
+```
 
-```
\ No newline at end of file