diff --git a/src/analysis/forecast_twl.py b/src/analysis/forecast_twl.py
index 4628daa..d0e389f 100644
--- a/src/analysis/forecast_twl.py
+++ b/src/analysis/forecast_twl.py
@@ -14,7 +14,6 @@ logger = setup_logging()
 
 MULTIPROCESS_THREADS = int(os.environ.get("MULTIPROCESS_THREADS", 4))
 
-
 def forecast_twl(
     df_tides,
     df_profiles,
@@ -103,11 +102,12 @@ def mean_slope_for_site_id(
     """
 
     # Get the prestorm beach profile
-    profile = df_profiles.query("site_id =='{}' and profile_type == '{}'".format(site_id, profile_type))
+    profile = df_profiles.loc[(site_id, profile_type)]
     profile_x = profile.index.get_level_values("x").tolist()
     profile_z = profile.z.tolist()
 
-    df_twl_site = df_twl.query("site_id == '{}'".format(site_id))
+    idx = pd.IndexSlice
+    df_twl_site = df_twl.loc[idx[site_id, :], :]
 
     df_beta = df_twl_site.apply(
         lambda row: slope_from_profile(
diff --git a/src/data/csv_to_geojson.py b/src/data/csv_to_geojson.py
index 40516ef..d7481b3 100644
--- a/src/data/csv_to_geojson.py
+++ b/src/data/csv_to_geojson.py
@@ -73,8 +73,8 @@ def R_high_to_geojson(sites_csv, profiles_csv, impacts_csv, output_geojson):
             # Find lat/lon of R_high position
             R_high_z = row["R_high"]
 
-            # Get poststorm profile (or should this be prestorm?)
-            df_profile = df_profiles.query('site_id=="{}" & profile_type=="prestorm"'.format(index))
+            # Get poststorm profile
+            df_profile = df_profiles.loc[(site_id, "prestorm")]
             int_x = crossings(df_profile.index.get_level_values("x").tolist(), df_profile.z.tolist(), R_high_z)
 
             # Take most landward interesection. Continue to next site if there is no intersection