diff --git a/application.py b/application.py
index a62270d3..3796d1b8 100644
--- a/application.py
+++ b/application.py
@@ -6,6 +6,7 @@
 from config import SITE_OFFLINE, geojson_names
 from marshmallow import Schema, fields, validate, ValidationError
 import re
+import threading
 
 from luts import (
     fire_weather_ops,
@@ -28,6 +29,13 @@
 
 app.register_blueprint(routes)
 
+# attach a custom attribute directly to the Flask app instance
+# this dictionary is the in-memory store for uploaded polygons
+# each key will be a unique ID and each value is a dict of metadata
+# see routes/upload_polygon.py for usage
+app.uploaded_polygons = {}
+app.store_lock = threading.Lock()
+
 
 def get_service_categories():
     """
diff --git a/fetch_data.py b/fetch_data.py
index 1c507e36..af7aeca4 100644
--- a/fetch_data.py
+++ b/fetch_data.py
@@ -501,3 +501,18 @@ def cftime_value_to_ymd(time_value, base_date):
     """Convert a time value in days since the base date to a year, month, day tuple."""
     date = base_date + datetime.timedelta(days=time_value)
     return date.year, date.month, date.day
+
+
+def create_gdf_from_geojson(geojson, in_crs, out_crs):
+    """Create a GeoDataFrame from GeoJSON and CRS.
+    Both input and output CRS must be specified.
+    Args:
+        geojson (dict): GeoJSON geometry
+        in_crs (int): EPSG CRS code
+        out_crs (int): EPSG CRS code
+    Returns:
+        poly (GeoDataFrame): GeoDataFrame of the polygon
+    """
+    geojson_dict = json.loads(geojson)
+    poly = gpd.GeoDataFrame.from_features(geojson_dict).set_crs(in_crs).to_crs(out_crs)
+    return poly
diff --git a/routes/__init__.py b/routes/__init__.py
index 4956249a..113372ea 100644
--- a/routes/__init__.py
+++ b/routes/__init__.py
@@ -40,3 +40,4 @@ def enforce_site_offline():
 from .places import *
 from .era5wrf import *
 from .fire_weather import *
+from .upload_polygon import *
diff --git a/routes/cmip6_downscaled.py b/routes/cmip6_downscaled.py
index a47b434d..9152e27a 100644
--- a/routes/cmip6_downscaled.py
+++ b/routes/cmip6_downscaled.py
@@ -1,17 +1,33 @@
 import asyncio
 import logging
-from flask import Blueprint, render_template, request
+import shapely
+import numpy as np
+from flask import Blueprint, render_template, request, current_app, jsonify
 
 # local imports
 from generate_urls import generate_wcs_query_url
-from generate_requests import generate_wcs_getcov_str
-from fetch_data import fetch_data, describe_via_wcps
+from generate_requests import generate_wcs_getcov_str, generate_netcdf_wcs_getcov_str
+from fetch_data import (
+    fetch_data,
+    describe_via_wcps,
+    get_poly,
+    create_gdf_from_geojson,
+    fetch_bbox_netcdf,
+)
 from validate_request import (
     latlon_is_numeric_and_in_geodetic_range,
     construct_latlon_bbox_from_coverage_bounds,
     validate_latlon_in_bboxes,
     project_latlon,
     generate_time_index_from_coverage_metadata,
+    check_for_uploaded_polygon,
+    validate_var_id,
+)
+from zonal_stats import (
+    get_scale_factor,
+    rasterize_polygon,
+    interpolate,
+    calculate_zonal_means_vectorized,
 )
 from postprocessing import postprocess, prune_nulls_with_max_intensity
 from csv_functions import create_csv
@@ -36,6 +52,29 @@ async def get_cmip6_metadata(cov_id):
     return metadata
 
 
+async def fetch_cmip6_downscaled_area_data(cov_id, polygon):
+    """
+    Make an async request for CMIP6 downscaled daily data for provided coverage within a specified polygon
+
+    Args:
+        cov_id (str): coverage ID
+        polygon (shapely.geometry.Polygon): polygon geometry
+
+    Returns:
+        list of data results within the specified polygon
+    """
+
+    wcs_str = generate_netcdf_wcs_getcov_str(polygon.total_bounds, cov_id=cov_id)
+
+    # Generate the URL for the WCS query
+    url = generate_wcs_query_url(wcs_str)
+
+    # Fetch the data
+    area_data_list = await fetch_bbox_netcdf([url])
+
+    return area_data_list
+
+
 async def fetch_cmip6_downscaled_point_data(cov_id, x, y):
     """
     Make an async request for CMIP6 downscaled daily data for provided coverage at a specified point
@@ -86,6 +125,36 @@ def package_cmip6_downscaled_data(metadata, point_data_list):
     return di
 
 
+def package_cmip6_downscaled_data_area(metadata, time_series_means):
+    """
+    Package CMIP6 downscaled daily values into human-readable JSON format
+
+    Args:
+        metadata (dict): coverage metadata
+        time_series_means (list): list of zonal mean values calculated from Rasdaman WCPS query results
+
+    Returns:
+        di (dict): time series dictionary of date/value pairs
+    """
+    di = dict()
+    try:
+        time_series = generate_time_index_from_coverage_metadata(metadata)
+
+        for idx, value in enumerate(time_series_means):
+            time = time_series[idx].date().strftime("%Y-%m-%d")
+
+            # Handle None and NaN values (required for zonal stats values derived from dataset)
+            if value is None or np.isnan(value):
+                di[time] = np.nan
+            else:
+                di[time] = round(float(value))
+    except Exception as exc:
+        logger.error(f"Error packaging area data: {exc}")
+        raise exc
+
+    return di
+
+
 @routes.route("/cmip6_downscaled/")
 def cmip6_downscaled_about():
     return render_template("/documentation/cmip6_downscaled.html")
@@ -254,3 +323,213 @@ def run_fetch_cmip6_downscaled_point_data(lat, lon, varname, model, scenario):
     point_data_list = asyncio.run(fetch_cmip6_downscaled_point_data(cov_id, x, y))
     results = package_cmip6_downscaled_data(metadata, point_data_list)
     return results
+
+
+@routes.route("/cmip6_downscaled/area/<place_id>")
+def cmip6_downscaled_area(place_id):
+    """
+    Fetch CMIP6 downscaled daily data for a specified area, compute zonal stats for each variable/model/scenario,
+    then combine them all into a single nested dictionary.
+
+    Args:
+        place_id (str): place identifier
+
+    Returns:
+        dict: time series data for the specified point for provided variables/models/scenarios
+
+    Notes:
+        example request (all variables): /cmip6_downscaled/area/1908030609
+        example request (specific variable): /cmip6_downscaled/area/1908030609?vars=tasmax
+        example request (specific model): /cmip6_downscaled/area/1908030609models=6ModelAvg
+        example request (specific scenario): /cmip6_downscaled/area/1908030609?scenarios=ssp585
+    """
+
+    ####### BYO-POLYGON IMPLEMENTATION ########
+
+    if check_for_uploaded_polygon(place_id):
+        # if place_id matches uploaded polygon key, create a GeoDataFrame from it, converting to 3338 for zonal stats
+        polygon = create_gdf_from_geojson(
+            current_app.uploaded_polygons[place_id]["geojson"],
+            current_app.uploaded_polygons[place_id]["crs"],
+            3338,
+        )
+
+    else:
+        # otherwise, do the standard process of validating place_id and getting polygon from GeoServer
+        poly_type = validate_var_id(place_id)
+        if type(poly_type) is tuple:
+            return poly_type
+        try:
+            polygon = get_poly(place_id)
+        except:
+            return render_template("422/invalid_area.html"), 422
+
+    # NOTE: uncomment this line if you want to test that your uploaded polygon is being read correctly
+    # return jsonify(polygon.to_json())
+
+    ###### STANDARD ZONAL STATS IMPLEMENTATION ########
+
+    # Split and assign optional HTTP GET parameters.
+    if request.args.get("vars"):
+        vars = request.args.get("vars").split(",")
+        if not all(var in all_cmip6_downscaled_vars for var in vars):
+            return render_template("400/bad_request.html"), 400
+        logger.debug(f"Results limited to vars: {vars}")
+    else:
+        vars = all_cmip6_downscaled_vars
+
+    if request.args.get("models"):
+        models = request.args.get("models").split(",")
+        if not all(model in all_cmip6_downscaled_models for model in models):
+            return render_template("400/bad_request.html"), 400
+        logger.debug(f"Results limited to models: {models}")
+    else:
+        models = all_cmip6_downscaled_models
+
+    if request.args.get("scenarios"):
+        scenarios = request.args.get("scenarios").split(",")
+        if not all(
+            scenario in all_cmip6_downscaled_scenarios for scenario in scenarios
+        ):
+            return render_template("400/bad_request.html"), 400
+        logger.debug(f"Results limited to scenarios: {scenarios}")
+    else:
+        scenarios = all_cmip6_downscaled_scenarios
+
+    try:
+        results = fetch_all_requested_combos_area(polygon, vars, models, scenarios)
+        if isinstance(results, tuple):
+            return results
+        # TODO: implement CSV export for area-based requests, making sure to handle BYO-Polygon name
+        # if request.args.get("format") == "csv":
+        #     place_id = request.args.get("community")
+        #     return create_csv(
+        #         results,
+        #         "cmip6_downscaled",
+        #         place_id, # <<< how to use use uploaded polygon name downstream?
+        #         vars=vars,
+        #     )
+    except ValueError:
+        return render_template("400/bad_request.html"), 400
+    except Exception as exc:
+        if hasattr(exc, "status") and exc.status == 404:
+            return render_template("404/no_data.html"), 404
+        return render_template("500/server_error.html"), 500
+
+    return results
+
+
+def fetch_all_requested_combos_area(polygon, vars, models, scenarios):
+    """
+    Query each variable/model/scenario coverage individually and combine them all into a nested dictionary.
+
+    Args:
+        polygon (GeoDataFrame): requested polygon
+        vars (list): list of variable names
+        models (list): list of model names
+        scenarios (list): list of scenario names
+
+    Returns:
+        dict: combined time series data for the specified area
+    """
+    results = {}
+    for varname in vars:
+        for model in models:
+            for scenario in scenarios:
+                # Return immediately if an exception is encountered for any coverage.
+                # All coverages share the same BBOX and structure, so an exception for one
+                # can be assumed to be an exception for all.
+                result = run_fetch_cmip6_downscaled_area_data(
+                    polygon, varname, model, scenario
+                )
+                if isinstance(result, tuple) and result[1] in [400, 404, 422]:
+                    return result
+                else:
+                    if model not in results:
+                        results[model] = {}
+                    if scenario not in results[model]:
+                        results[model][scenario] = {}
+                    for time, value in result.items():
+                        if time not in results[model][scenario]:
+                            results[model][scenario][time] = {}
+                        results[model][scenario][time][varname] = value
+
+    results = prune_nulls_with_max_intensity(postprocess(results, "cmip6_downscaled"))
+    return results
+
+
+def run_fetch_cmip6_downscaled_area_data(polygon, varname, model, scenario):
+    """
+    Fetch CMIP6 downscaled daily data for a single variable/model/scenario combo over a specified area.
+
+    Args:
+        polygon (GeoDataFrame): requested polygon
+        varname (str): variable name
+        model (str): model name
+        scenario (str): scenario name
+
+    Returns:
+        dict: time series data for the specified area for a single variable/model/scenario combo
+    """
+
+    # If we have made it this far, the model and scenario are valid.
+    # If they are not found for the variable, return empty results.
+    if model not in cmip6_downscaled_options[varname]:
+        return {}
+    if scenario not in cmip6_downscaled_options[varname][model]:
+        return {}
+
+    cov_id = f"cmip6_downscaled_{varname}_{model}_{scenario}_wcs"
+    cov_id = cov_id.replace("-", "_")
+    metadata = asyncio.run(get_cmip6_metadata(cov_id))
+
+    # TODO: validate that polygon is within coverage bounds
+
+    area_dataset = asyncio.run(fetch_cmip6_downscaled_area_data(cov_id, polygon))
+
+    time_series_means = calculate_cmip6_downscaled_zonal_stats(
+        polygon, area_dataset, varname
+    )
+
+    results = package_cmip6_downscaled_data_area(metadata, time_series_means)
+    return results
+
+
+def calculate_cmip6_downscaled_zonal_stats(polygon, area_dataset, varname):
+    """Process zonal statistics for a dataset.
+
+    Args:
+        polygon (GeoDataFrame): Target polygon
+        area_dataset (xarray.Dataset): Input dataset
+        varname (str): Variable name
+
+    Returns:
+        list: daily zonal means (one value per time step)
+    """
+    ds = area_dataset
+
+    # get scale factor once, not per variable or time slice!
+    spatial_resolution = ds.rio.resolution()
+    grid_cell_area_m2 = abs(spatial_resolution[0]) * abs(spatial_resolution[1])
+    polygon_area_m2 = polygon.area
+    scale_factor = get_scale_factor(grid_cell_area_m2, polygon_area_m2)
+
+    # create an initial array for the basis of polygon rasterization
+    # why? polygon rasterization bogs down hard when doing it in the loop
+    da_i = interpolate(
+        ds.isel(ansi=0), varname, "X", "Y", scale_factor, method="nearest"
+    )
+
+    rasterized_polygon_array = rasterize_polygon(da_i, "X", "Y", polygon)
+
+    da_i_3d = interpolate(ds, varname, "X", "Y", scale_factor, method="nearest")
+    # calculate zonal stats for the entire time series
+    # rename the ansi dimensions to "time" so this function will work
+    da_i_3d = da_i_3d.rename({"ansi": "time"})
+    time_series_means = calculate_zonal_means_vectorized(
+        da_i_3d, rasterized_polygon_array, "X", "Y"
+    )
+
+    print(varname, len(time_series_means))
+
+    return time_series_means
diff --git a/routes/upload_polygon.py b/routes/upload_polygon.py
new file mode 100644
index 00000000..c645319f
--- /dev/null
+++ b/routes/upload_polygon.py
@@ -0,0 +1,188 @@
+from flask import Blueprint, request, jsonify, render_template, current_app
+import geopandas as gpd
+import io, zipfile, secrets, time, tempfile, os
+from datetime import datetime, timezone
+
+from . import routes
+
+upload_polygon = Blueprint("upload_polygon", __name__)
+
+
+@routes.route("/upload_polygon", methods=["POST"])
+def upload_polygon():
+    """
+    Accepts an uploaded shapefile (as a ZIP containing .shp, .shx, .dbf, .prj)
+    and stores the geometry in memory with a unique ID. The ID can then be used to
+    reference the polygon in subsequent requests from other endpoints. Optionally,
+    a user-defined name for the polygon can be provided via form data - this would
+    be used for downstream display purposes, for naming output files, etc.
+
+    Upload via curl example:
+        curl -F "file=@my_shapefile.zip" \
+             -F "name=my custom polygon" \
+             http://localhost:5000/upload_polygon 
+
+    #TODO: implement form in HTML template for browser-based upload
+    #TODO: implement additional validation (file size limit, filename checks, upload limit per IP, etc.)
+    #TODO: implement cleanup of old in-memory polygons after expiry time (e.g., 1 hour)
+    #TODO: encourage users to create polygons via https://geojson.io/#map=3.92/63/-154.56 , which
+    # automatically zips the shapefile components for download.
+
+    Example response JSON:
+
+        {
+            "expires_at": "2025-11-04T18:48:02.522820+00:00",
+            "name": "my custom polygon",
+            "polygon_id": "f-EMX4YL",
+            "uploaded_at": "2025-11-04T17:48:02.522820+00:00"
+        }
+
+    Example of current_app.uploaded_polygons content after upload:
+
+        {
+            'f-EMX4YL': {
+                'geojson': ,
+                'expires_at': "2025-11-04T18:48:02.522820+00:00",
+                'uploaded_at': "2025-11-04T17:48:02.522820+00:00",
+                'name': 'my custom polygon',
+                'crs': 4326
+            }
+        }
+
+    Returns:
+        JSON with polygon_id, polygon name, and expiry time, or error message.
+    """
+
+    # check for file in request
+    if "file" not in request.files:
+        return (
+            jsonify({"error": "No file part in request"}),
+            400,
+        )
+    file = request.files["file"]
+    if file.filename == "":
+        return (
+            jsonify({"error": "No selected file"}),
+            400,
+        )
+
+    # optional user-defined polygon name
+    user_name = request.form.get("name", "").strip()
+    if not user_name:
+        user_name = "unnamed_polygon"
+
+    # read the uploaded zip into memory
+    file_bytes = file.read()
+
+    # validate that it's a zip containing a .shp file
+    try:
+        with zipfile.ZipFile(io.BytesIO(file_bytes)) as zf:
+            shapefile_names = [n for n in zf.namelist() if n.lower().endswith(".shp")]
+            if not shapefile_names:
+                return (
+                    jsonify({"error": "No .shp file found in ZIP"}),
+                    400,
+                )
+    except zipfile.BadZipFile:
+        return (
+            jsonify({"error": "Uploaded file is not a valid ZIP archive"}),
+            400,
+        )
+
+    # write zip file contents to a temporary file so GeoPandas/Fiona can read it
+
+    # NOTE: this actually writes to disk because GeoPandas (via Fiona) cannot read
+    # directly from an in-memory object (BytesIO) when dealing with shapefiles inside ZIP archives.
+    # Shapefiles consist of multiple files (.shp, .shx, .dbf, .prj) that need to be accessed together.
+    # Fiona delegates to GDAL, which expects either a folder path, or a “virtual file system” path
+    # like zip://path_to_zip!inner_path.
+    # we write the uploaded ZIP to disk temporarily and give GDAL that real path,
+    # then delete the temp file after reading.
+
+    with tempfile.NamedTemporaryFile(suffix=".zip", delete=False) as tmp:
+        tmp.write(file_bytes)
+        tmp_path = tmp.name
+
+    try:
+        with zipfile.ZipFile(tmp_path) as zf:
+            shapefile_names = [n for n in zf.namelist() if n.lower().endswith(".shp")]
+            if not shapefile_names:
+                return (
+                    jsonify({"error": "No .shp file found in ZIP"}),
+                    400,
+                )
+            shapefile_path = shapefile_names[0]
+
+        # note the exclamation mark syntax:
+        gdf = gpd.read_file(f"zip://{tmp_path}!{shapefile_path}")
+
+    except Exception as e:
+        return (
+            jsonify({"error": f"Invalid shapefile ZIP: {str(e)}"}),
+            400,
+        )
+    finally:
+        # clean up temporary file
+        if os.path.exists(tmp_path):
+            os.remove(tmp_path)
+
+    # check for features
+    if gdf.empty:
+        return (
+            jsonify({"error": "Uploaded shapefile contains no features"}),
+            400,
+        )
+
+    # check for a valid CRS
+    if gdf.crs is None:
+        return (
+            jsonify(
+                {
+                    "error": "Uploaded shapefile has no defined coordinate reference system (CRS)"
+                }
+            ),
+            400,
+        )
+    else:  # get 4 digit EPSG code as int to save in uploaded polygons dict
+        crs_int = gdf.crs.to_epsg()
+
+    # combine all polygons into one geometry
+    # if the user uploads a shapefile with multiple features, they are treated as one feature
+    # we could decide to loop through each feature and process them separately,
+    # but for now we just dissolve them to force a single geometry
+    gdf_dissolved = gdf.dissolve()
+
+    # create valid geojson from this dissolved gdf
+    geojson = gdf_dissolved.to_json()
+
+    # generate a unique polygon ID
+    poly_id = secrets.token_urlsafe(6)
+
+    # store geometry, upload time, and user-defined name in the app's in-memory store
+
+    upload_time = datetime.fromtimestamp(time.time(), tz=timezone.utc).isoformat()
+    # 3600 second = 1 hour expiry
+    expiry_time = datetime.fromtimestamp(
+        time.time() + 3600, tz=timezone.utc
+    ).isoformat()
+
+    with current_app.store_lock:
+        current_app.uploaded_polygons[poly_id] = {
+            "geojson": geojson,
+            "uploaded_at": upload_time,
+            "expires_at": expiry_time,
+            "name": user_name,
+            "crs": crs_int,
+        }
+
+    return (
+        jsonify(
+            {
+                "polygon_id": poly_id,
+                "name": user_name,
+                "expires_at": expiry_time,
+                "uploaded_at": upload_time,
+            }
+        ),
+        200,
+    )
diff --git a/validate_request.py b/validate_request.py
index 44ac7659..cf379494 100644
--- a/validate_request.py
+++ b/validate_request.py
@@ -7,8 +7,9 @@
 import rasterio
 import os.path
 import os
+import shapely
 
-from flask import render_template
+from flask import render_template, current_app, jsonify
 from pyproj import Transformer
 import numpy as np
 import pandas as pd
@@ -565,3 +566,18 @@ def get_coverage_crs_str(coverage_metadata):
         )
 
     return crs.to_string()
+
+
+def check_for_uploaded_polygon(place_id):
+    """Check if a user-defined polygon exists for the given place_id.
+
+    Args:
+        place_id (str): The place identifier to check for an uploaded polygon. This is the a unique token generated when the user uploads a polygon.
+    Returns:
+        True if a user-defined polygon exists, False otherwise.
+    """
+    with current_app.store_lock:
+        polygon = current_app.uploaded_polygons.get(place_id)
+    if polygon:
+        return True
+    return False