Great_Salt_Lake_GEE_Time_Series_Product_Code.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import geemap\n",
    "import ee\n",
    "from datetime import date\n",
    "import datetime\n",
    "from ipyfilechooser import FileChooser\n",
    "import pandas as pd\n",
    "ee.Initialize()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "NA_export_region = ee.Geometry.Polygon(\n",
    "        [[[-113.19195907269143, 41.768481881452075],\n",
    "          [-113.19195907269143, 41.21713179714162],\n",
    "          [-112.44626205120706, 41.21713179714162],\n",
    "          [-112.44626205120706, 41.768481881452075]]])\n",
    "\n",
    "SA_export_region = ee.Geometry.Polygon(\n",
    "        [[[-112.8838305102725, 41.2295267918564],\n",
    "          [-112.8838305102725, 40.638172391988064],\n",
    "          [-111.9390062915225, 40.638172391988064],\n",
    "          [-111.9390062915225, 41.2295267918564]]])\n",
    "\n",
    "def NA_Export_Clip(image):\n",
    "    return image.clip(NA_export_region).copyProperties(image)\n",
    "def SA_Export_Clip(image):\n",
    "    return image.clip(SA_export_region).copyProperties(image)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "### Shapefile Asset intake ###\n",
    "GSL_SA_boundary = ee.FeatureCollection(\"projects/ee-radwinski/assets/GSL_SouthArm_1984_Boundary\")\n",
    "GSL_NA_boundary = ee.FeatureCollection(\"projects/ee-radwinski/assets/GSL_NorthArm_Boundary_Final\")\n",
    "GSL_NA_boundary_old = ee.FeatureCollection(\"projects/ee-radwinski/assets/GSL_NorthArm_1984_Boundary_Expanded\")\n",
    "GSL_total_boundary = GSL_NA_boundary.merge(GSL_SA_boundary).union()\n",
    "BSF_boundary = ee.FeatureCollection(\"projects/ee-radwinski/assets/BSF_Boundary_1984\")\n",
    "BSF_late_boundary = ee.FeatureCollection(\"projects/ee-radwinski/assets/BSF_Extent_1986_Smaller\")\n",
    "NFB_boundary = ee.FeatureCollection(\"projects/ee-radwinski/assets/NFB_Boundary_1984\")\n",
    "NA_SaltPool = ee.FeatureCollection(\"projects/ee-radwinski/assets/SaltPoolBoundaryUpdated\")\n",
    "NA_SaltPool_old = ee.FeatureCollection(\"projects/ee-radwinski/assets/GSL_NA_SaltPool_toMask\")\n",
    "BB_boundary = ee.FeatureCollection(\"projects/ee-radwinski/assets/BonnevilleBasinShape\")\n",
    "BB_cropped_boundary = ee.FeatureCollection(\"projects/ee-radwinski/assets/Bonneville_Basin_Cropped_Boundary\")\n",
    "BB_and_GSL_boundary = ee.FeatureCollection(\"projects/ee-radwinski/assets/All_BB_Seds_and_Shoreline_Boundary\")\n",
    "BRB_boundary = ee.FeatureCollection(\"projects/ee-radwinski/assets/BRB_Boundary\")\n",
    "BRB_SaltPool = ee.FeatureCollection(\"projects/ee-radwinski/assets/BRB_SaltPools\")\n",
    "GSL_export_region = ee.Geometry.Polygon([[[-113.21733207756196, 41.80229778681049], [-113.21733207756196, 40.63715498285718], [-111.92643852287446, 40.63715498285718], [-111.92643852287446, 41.80229778681049]]])\n",
    "BSF_export_region = ee.Geometry.Polygon([[[-114.10722465568696, 40.9324608917365], [-114.10722465568696, 40.55373327911588], [-113.59910697990571, 40.55373327911588], [-113.59910697990571, 40.9324608917365]]])\n",
    "USGS_water_station_region = ee.Geometry.Point([-112.23, 41.066944]).buffer(500)\n",
    "landsat_tile_overflow_area = ee.Geometry.Polygon([[[-111.67547275476988, 40.56441023109152], [-111.70293857508238, 41.55020817425386], [-112.08141905793187, 41.54589403254853], [-112.0778046333898, 41.32558038715473], [-112.07547209577153, 41.17039827691482], [-112.07337928992675, 41.063050532354744], [-112.07193971329473, 40.96368190627725], [-112.06982741402413, 40.83717571457464], [-112.06538696289063, 40.557904349873844]]])\n",
    "farmington_bay_salt_mask_area = ee.Geometry.Polygon(\n",
    "        [[[-112.08273604182838, 41.21786944789319],\n",
    "          [-112.21594527034401, 41.23026430280494],\n",
    "          [-112.29559614925026, 41.212704231634646],\n",
    "          [-112.26401045589088, 41.16206353011265],\n",
    "          [-112.25851729182838, 41.12172967661335],\n",
    "          [-112.23242476253151, 41.09689651511007],\n",
    "          [-112.22418501643776, 41.083441297796334],\n",
    "          [-112.23517134456276, 41.05755813551418],\n",
    "          [-112.22967818050026, 41.05030902417058],\n",
    "          [-112.2046156194651, 41.043835928336655],\n",
    "          [-112.1936292913401, 41.04176440318438],\n",
    "          [-112.18195631770729, 41.032700713996185],\n",
    "          [-112.17268660335182, 40.99980199072419],\n",
    "          [-112.15826704768776, 40.985549413879326],\n",
    "          [-112.14651167265878, 40.96831781679227],\n",
    "          [-112.13930189482674, 40.95016897806219],\n",
    "          [-112.13071882597909, 40.93486822434222],\n",
    "          [-112.1097761379908, 40.90996450684567],\n",
    "          [-112.09226667754159, 40.902439663403406],\n",
    "          [-112.07819044463143, 40.899585188456534],\n",
    "          [-112.0658308254908, 40.89802815019883],\n",
    "          [-112.06720411650643, 40.966503157478996]]])\n",
    "SA_Salt_Mine_Area = ee.Geometry.Polygon(\n",
    "        [[[-112.74155702451525, 40.93565853875015],\n",
    "          [-112.74696435788927, 40.93079532702177],\n",
    "          [-112.74335946897325, 40.91419286803555],\n",
    "          [-112.73443307737169, 40.91023619808089],\n",
    "          [-112.7036198602086, 40.9169818537473],\n",
    "          [-112.69872751096544, 40.91840873116959],\n",
    "          [-112.70379152158556, 40.925542656271915],\n",
    "          [-112.70370569089708, 40.93669597758743]]])\n",
    "### Functions ###\n",
    "def image_dater(image):\n",
    "    date = ee.Number(image.date().format('YYYY-MM-dd'))\n",
    "    return image.set({'Date_Filter': date})\n",
    "\n",
    "def image_grab(img_col, img_selector):\n",
    "    new_col = img_col.filter(ee.Filter.eq('Date_Filter', img_selector)) \n",
    "    new_col_list = new_col.toList(new_col.size())\n",
    "    return ee.Image(new_col_list.get(0))\n",
    "\n",
    "def maskL8clouds(image):\n",
    "    cloudBitMask = ee.Number(2).pow(3).int()\n",
    "    CirrusBitMask = ee.Number(2).pow(2).int()\n",
    "    qa = image.select('QA_PIXEL')\n",
    "    cloud_mask = qa.bitwiseAnd(cloudBitMask).eq(0)\n",
    "    cirrus_mask = qa.bitwiseAnd(CirrusBitMask).eq(0)\n",
    "    return image.updateMask(cloud_mask).updateMask(cirrus_mask)\n",
    "\n",
    "def MaskCloudsS2(image):\n",
    "  SCL = image.select('SCL')\n",
    "  CloudMask = SCL.neq(9)\n",
    "  return image.updateMask(CloudMask)\n",
    "\n",
    "def MaskWaterS2(image):\n",
    "  SCL = image.select('SCL')\n",
    "  CloudMask = SCL.neq(6)\n",
    "  return image.updateMask(CloudMask)\n",
    "\n",
    "def landsat_ndwi_fn(image):\n",
    "    ndwi_calc = image.normalizedDifference(['SR_B3', 'SR_B5']) #green-NIR / green+NIR -- full NDWI image\n",
    "    standing_water = ndwi_calc.rename('ndwi').copyProperties(image) \n",
    "    return standing_water\n",
    "\n",
    "def sentinel_ndwi_fn(image):\n",
    "    ndwi_calc = image.normalizedDifference(['B3', 'B8']) #green-NIR / green+NIR -- full NDWI image\n",
    "    standing_water = ndwi_calc.rename('ndwi').copyProperties(image) \n",
    "    return standing_water\n",
    "\n",
    "landsat_ndvi_threshold = 0.105 \n",
    "\n",
    "def landsat_ndvi_fn(image):\n",
    "    ndvi_calc = image.normalizedDifference(['SR_B5', 'SR_B4']) #NIR-RED/NIR+RED -- full NDVI image\n",
    "    vegetation = ndvi_calc.updateMask(ndvi_calc.gte(landsat_ndvi_threshold)).rename('ndvi').copyProperties(image) # subsets the image to just water pixels\n",
    "    return vegetation\n",
    "\n",
    "def landsat_BRB_ndvi_fn(image):\n",
    "    ndvi_calc = image.normalizedDifference(['SR_B5', 'SR_B4']) #NIR-RED/NIR+RED -- full NDVI image\n",
    "    vegetation = ndvi_calc.updateMask(ndvi_calc.gte(0.09)).rename('ndvi').copyProperties(image) # subsets the image to just water pixels\n",
    "    return vegetation\n",
    "\n",
    "sentinel_ndvi_threshold = 0.185 \n",
    "\n",
    "def sentinel_ndvi_fn(image):\n",
    "    ndvi_calc = image.normalizedDifference(['B8', 'B4']) #NIR-RED/NIR+RED -- full NDVI image\n",
    "    vegetation = ndvi_calc.updateMask(ndvi_calc.gte(sentinel_ndvi_threshold)).rename('ndvi').copyProperties(image) # subsets the image to just water pixels\n",
    "    return vegetation\n",
    "\n",
    "def landsat_halite_fn(image):\n",
    "    halite_index = image.normalizedDifference(['SR_B4', 'SR_B6']) # red-swir1 / red+swir1\n",
    "    halite = halite_index.updateMask(halite_index.gte(0.1)).rename('halite').copyProperties(image)\n",
    "    return halite \n",
    "\n",
    "landsat_halite_threshold = 0.345 \n",
    "\n",
    "def landsat_halite_fn_conservative(image):\n",
    "    halite_index = image.normalizedDifference(['SR_B4', 'SR_B6']) # red-swir1 / red+swir1\n",
    "    halite = halite_index.updateMask(halite_index.gte(landsat_halite_threshold)).rename('halite').copyProperties(image)  \n",
    "    return halite \n",
    "\n",
    "landsat_gypsum_threshold = 0.153 \n",
    "\n",
    "def landsat_gypsum_fn(image):\n",
    "    gypsum_index = image.normalizedDifference(['SR_B6', 'SR_B7'])\n",
    "    gypsum = gypsum_index.updateMask(gypsum_index.gte(landsat_gypsum_threshold)).rename('gypsum').copyProperties(image)\n",
    "    return gypsum\n",
    "\n",
    "sentinel_halite_threshold = 0.58 \n",
    "\n",
    "def sentinel_halite_fn(image):\n",
    "    halite_index = image.normalizedDifference(['B4', 'B11']) # red-swir1 / red+swir1\n",
    "    halite = halite_index.updateMask(halite_index.gte(sentinel_halite_threshold)).rename('halite').copyProperties(image) \n",
    "    return halite \n",
    "\n",
    "sentinel_gypsum_threshold = 0.3 \n",
    "\n",
    "def sentinel_gypsum_fn(image):\n",
    "    gypsum_index = image.normalizedDifference(['B11', 'B12']) # red-swir1 / red+swir1\n",
    "    gypsum = gypsum_index.updateMask(gypsum_index.gte(sentinel_gypsum_threshold)).rename('gypsum').copyProperties(image) \n",
    "    return gypsum\n",
    "\n",
    "def MaskToWaterL8(image):\n",
    "  WaterBitMask = ee.Number(2).pow(7).int()\n",
    "  qa = image.select('QA_PIXEL')\n",
    "  water_extract = qa.bitwiseAnd(WaterBitMask).neq(0)\n",
    "  return image.updateMask(water_extract).copyProperties(image)\n",
    "\n",
    "def MaskWaterLandsat(image):\n",
    "    WaterBitMask = ee.Number(2).pow(7).int()\n",
    "    qa = image.select('QA_PIXEL')\n",
    "    water_extract = qa.bitwiseAnd(WaterBitMask).eq(0)\n",
    "    water_pixels = qa.bitwiseAnd(WaterBitMask).neq(0)\n",
    "    water_count = image.select('SR_B5').rename('water_count').updateMask(water_pixels).reduceRegion(\n",
    "        reducer = ee.Reducer.count(),\n",
    "        geometry = BSF_late_boundary,\n",
    "        scale=30)\n",
    "    masked_image = image.updateMask(water_extract).copyProperties(image).set(water_count)\n",
    "    return masked_image\n",
    "\n",
    "def landsat5bandrename(img):\n",
    "    return img.select('SR_B1', 'SR_B2', 'SR_B3', 'SR_B4', 'SR_B5', 'SR_B7', 'QA_PIXEL').rename('SR_B2', 'SR_B3', 'SR_B4', 'SR_B5', 'SR_B6', 'SR_B7', 'QA_PIXEL')\n",
    "\n",
    "def NorthArmClip(image):\n",
    "    return image.clip(GSL_NA_boundary).copyProperties(image)\n",
    "\n",
    "def SouthArmClip(image):\n",
    "    return image.clip(GSL_SA_boundary).copyProperties(image)\n",
    "\n",
    "def TotalGSLClip(image):\n",
    "    return image.clip(GSL_total_boundary).copyProperties(image)\n",
    "\n",
    "def BRBClip(image):\n",
    "    return image.clip(BRB_boundary).copyProperties(image)\n",
    "\n",
    "def MaskSaltPool(image):\n",
    "    mask = ee.Image.constant(1).clip(NA_SaltPool).mask().eq(0)\n",
    "    return image.updateMask(mask).copyProperties(image)\n",
    "\n",
    "def MaskBRBSaltPool(image):\n",
    "    mask = ee.Image.constant(1).clip(BRB_SaltPool).mask().eq(0)\n",
    "    return image.updateMask(mask).copyProperties(image)\n",
    "\n",
    "def MaskLandsatOverflow(image):\n",
    "    mask = ee.Image.constant(1).clip(landsat_tile_overflow_area).mask().eq(0)\n",
    "    return image.updateMask(mask).copyProperties(image)\n",
    "\n",
    "def FarmingtonBaySaltMask(image):\n",
    "    mask = ee.Image.constant(1).clip(farmington_bay_salt_mask_area).mask().eq(0)\n",
    "    mask2 = ee.Image.constant(1).clip(SA_Salt_Mine_Area).mask().eq(0)\n",
    "    return image.updateMask(mask).updateMask(mask2).copyProperties(image)\n",
    "\n",
    "def GSLNA_WaterPixelCount(image):\n",
    "    area_image = ee.Image.pixelArea()\n",
    "    ndwi_calc = image.normalizedDifference(['SR_B3', 'SR_B5']).rename('ndwi')\n",
    "    histogram = ndwi_calc.select('ndwi').reduceRegion(\n",
    "        reducer = ee.Reducer.histogram(255, 2),\n",
    "        geometry = GSL_NA_boundary.geometry().buffer(6000),\n",
    "        scale = 30,\n",
    "        bestEffort= True,)\n",
    "    threshold = OtsuThreshold(histogram.get('ndwi')).add(0.15) \n",
    "    mask = ndwi_calc.select('ndwi').gte(threshold)\n",
    "    final = image.addBands(area_image)\n",
    "    stats = final.select('area').updateMask(mask).rename('ndwi').reduceRegion(\n",
    "        reducer = ee.Reducer.sum(),\n",
    "        geometry= GSL_NA_boundary,\n",
    "        scale=30,\n",
    "        bestEffort = True)\n",
    "    return image.set('ndwi', stats.get('ndwi')) \n",
    "\n",
    "def BRB_WaterPixelCount(image):\n",
    "    area_image = ee.Image.pixelArea()\n",
    "    ndwi_calc = image.normalizedDifference(['SR_B3', 'SR_B5']).rename('ndwi')\n",
    "    histogram = ndwi_calc.select('ndwi').reduceRegion(\n",
    "        reducer = ee.Reducer.histogram(255, 2),\n",
    "        geometry = BRB_boundary.geometry().buffer(6000),\n",
    "        scale = 30,\n",
    "        bestEffort= True,)\n",
    "    threshold = OtsuThreshold(histogram.get('ndwi')).add(0.15) \n",
    "    mask = ndwi_calc.select('ndwi').gte(threshold)\n",
    "    final = image.addBands(area_image)\n",
    "    stats = final.select('area').updateMask(mask).rename('ndwi').reduceRegion(\n",
    "        reducer = ee.Reducer.sum(),\n",
    "        geometry= BRB_boundary,\n",
    "        scale=30,\n",
    "        bestEffort = True)\n",
    "    return image.set('ndwi', stats.get('ndwi'))\n",
    "\n",
    "def GSLNA_WaterPixelCount_Sent(image):\n",
    "    area_image = ee.Image.pixelArea()\n",
    "    mask = image.select('ndwi').gte(sentinel_ndwi_threshold)\n",
    "    final = image.addBands(area_image)\n",
    "    stats = final.select('area').updateMask(mask).rename('ndwi').reduceRegion(\n",
    "        reducer = ee.Reducer.sum(),\n",
    "        geometry= GSL_NA_boundary,\n",
    "        scale=10,\n",
    "        maxPixels = 1e12)\n",
    "    return image.set('ndwi', stats.get('ndwi'))\n",
    "\n",
    "def BRB_WaterPixelCount_Sent(image):\n",
    "    area_image = ee.Image.pixelArea()\n",
    "    mask = image.select('ndwi').gte(sentinel_ndwi_threshold)\n",
    "    final = image.addBands(area_image)\n",
    "    stats = final.select('area').updateMask(mask).rename('ndwi').reduceRegion(\n",
    "        reducer = ee.Reducer.sum(),\n",
    "        geometry= BRB_boundary,\n",
    "        scale=10,\n",
    "        maxPixels = 1e12)\n",
    "    return image.set('ndwi', stats.get('ndwi'))\n",
    "\n",
    "def GSLSA_WaterPixelCount(image):\n",
    "    area_image = ee.Image.pixelArea()\n",
    "    ndwi_calc = image.normalizedDifference(['SR_B3', 'SR_B5']).rename('ndwi')\n",
    "    histogram = ndwi_calc.select('ndwi').reduceRegion(\n",
    "        reducer = ee.Reducer.histogram(255, 2),\n",
    "        geometry = GSL_SA_boundary.geometry().buffer(6000),\n",
    "        scale = 30,\n",
    "        bestEffort= True,)\n",
    "    threshold = OtsuThreshold(histogram.get('ndwi')).add(0.15) \n",
    "    mask = ndwi_calc.select('ndwi').gte(threshold)\n",
    "    final = image.addBands(area_image)\n",
    "    stats = final.select('area').updateMask(mask).rename('ndwi').reduceRegion(\n",
    "        reducer = ee.Reducer.sum(),\n",
    "        geometry= GSL_SA_boundary,\n",
    "        scale=30,\n",
    "        bestEffort = True)\n",
    "    return image.set('ndwi', stats.get('ndwi'))\n",
    "\n",
    "def GSLSA_WaterPixelCount_Sent(image):\n",
    "    area_image = ee.Image.pixelArea()\n",
    "    mask = image.select('ndwi').gte(sentinel_ndwi_threshold)\n",
    "    final = image.addBands(area_image)\n",
    "    stats = final.select('area').updateMask(mask).rename('ndwi').reduceRegion(\n",
    "        reducer = ee.Reducer.sum(),\n",
    "        geometry= GSL_SA_boundary,\n",
    "        scale=10,\n",
    "        maxPixels = 1e12)\n",
    "    return image.set('ndwi', stats.get('ndwi'))\n",
    "\n",
    "def GSLNA_VegPixelCount(image):\n",
    "    area_image = ee.Image.pixelArea()\n",
    "    mask = image.select('ndvi').gte(landsat_ndvi_threshold)\n",
    "    final = image.addBands(area_image)\n",
    "    stats = final.select('area').updateMask(mask).rename('ndvi').reduceRegion(\n",
    "        reducer = ee.Reducer.sum(),\n",
    "        geometry= GSL_NA_boundary,\n",
    "        scale=30,\n",
    "        maxPixels = 1e12)\n",
    "    return image.set('ndvi', stats.get('ndvi'))\n",
    "\n",
    "def BRB_VegPixelCount(image):\n",
    "    area_image = ee.Image.pixelArea()\n",
    "    mask = image.select('ndvi').gte(landsat_ndvi_threshold)\n",
    "    final = image.addBands(area_image)\n",
    "    stats = final.select('area').updateMask(mask).rename('ndvi').reduceRegion(\n",
    "        reducer = ee.Reducer.sum(),\n",
    "        geometry= BRB_boundary,\n",
    "        scale=30,\n",
    "        maxPixels = 1e12)\n",
    "    return image.set('ndvi', stats.get('ndvi'))\n",
    "\n",
    "def GSLNA_VegPixelCount_Sent(image):\n",
    "    area_image = ee.Image.pixelArea()\n",
    "    mask = image.select('ndvi').gte(sentinel_ndvi_threshold)\n",
    "    final = image.addBands(area_image)\n",
    "    stats = final.select('area').updateMask(mask).rename('ndvi').reduceRegion(\n",
    "        reducer = ee.Reducer.sum(),\n",
    "        geometry= GSL_NA_boundary,\n",
    "        scale=10,\n",
    "        maxPixels = 1e12)\n",
    "    return image.set('ndvi', stats.get('ndvi'))\n",
    "\n",
    "def BRB_VegPixelCount_Sent(image):\n",
    "    area_image = ee.Image.pixelArea()\n",
    "    mask = image.select('ndvi').gte(sentinel_ndvi_threshold)\n",
    "    final = image.addBands(area_image)\n",
    "    stats = final.select('area').updateMask(mask).rename('ndvi').reduceRegion(\n",
    "        reducer = ee.Reducer.sum(),\n",
    "        geometry= BRB_boundary,\n",
    "        scale=10,\n",
    "        maxPixels = 1e12)\n",
    "    return image.set('ndvi', stats.get('ndvi'))\n",
    "\n",
    "def GSLSA_VegPixelCount(image):\n",
    "    area_image = ee.Image.pixelArea()\n",
    "    mask = image.select('ndvi').gte(landsat_ndvi_threshold)\n",
    "    final = image.addBands(area_image)\n",
    "    stats = final.select('area').updateMask(mask).rename('ndvi').reduceRegion(\n",
    "        reducer = ee.Reducer.sum(),\n",
    "        geometry= GSL_SA_boundary,\n",
    "        scale=30,\n",
    "        maxPixels = 1e12)\n",
    "    return image.set('ndvi', stats.get('ndvi'))\n",
    "\n",
    "def GSLSA_VegPixelCount_Sent(image):\n",
    "    area_image = ee.Image.pixelArea()\n",
    "    mask = image.select('ndvi').gte(sentinel_ndvi_threshold)\n",
    "    final = image.addBands(area_image)\n",
    "    stats = final.select('area').updateMask(mask).rename('ndvi').reduceRegion(\n",
    "        reducer = ee.Reducer.sum(),\n",
    "        geometry= GSL_SA_boundary,\n",
    "        scale=10,\n",
    "        maxPixels = 1e12)\n",
    "    return image.set('ndvi', stats.get('ndvi'))\n",
    "\n",
    "def halite_pixel_count_landsat(image):\n",
    "    stats = image.gt(landsat_halite_threshold).select('halite').pixelArea().rename('halite').reduceRegion(\n",
    "        reducer = ee.Reducer.sum(),\n",
    "        geometry= BSF_boundary,\n",
    "        scale=30)       \n",
    "    return image.set('halite', stats.get('halite'))\n",
    "\n",
    "def halite_pixel_count_landsat_GSLNA(image):\n",
    "    area_image = ee.Image.pixelArea()\n",
    "    mask = image.select('halite').gte(landsat_halite_threshold)\n",
    "    final = image.addBands(area_image)\n",
    "    stats = final.select('area').updateMask(mask).rename('halite').reduceRegion(\n",
    "        reducer = ee.Reducer.sum(),\n",
    "        geometry= GSL_NA_boundary,\n",
    "        scale=30,\n",
    "        maxPixels = 1e12)      \n",
    "    return image.set('halite', stats.get('halite'))\n",
    "\n",
    "def halite_pixel_count_landsat_BRB(image):\n",
    "    area_image = ee.Image.pixelArea()\n",
    "    mask = image.select('halite').gte(landsat_halite_threshold)\n",
    "    final = image.addBands(area_image)\n",
    "    stats = final.select('area').updateMask(mask).rename('halite').reduceRegion(\n",
    "        reducer = ee.Reducer.sum(),\n",
    "        geometry= BRB_boundary,\n",
    "        scale=30,\n",
    "        maxPixels = 1e12)      \n",
    "    return image.set('halite', stats.get('halite'))\n",
    "\n",
    "def halite_pixel_count_landsat_GSLSA(image):\n",
    "    area_image = ee.Image.pixelArea()\n",
    "    mask = image.select('halite').gte(landsat_halite_threshold)\n",
    "    final = image.addBands(area_image)\n",
    "    stats = final.select('area').updateMask(mask).rename('halite').reduceRegion(\n",
    "        reducer = ee.Reducer.sum(),\n",
    "        geometry= GSL_SA_boundary,\n",
    "        scale=30,\n",
    "        maxPixels = 1e12)       \n",
    "    return image.set('halite', stats.get('halite'))\n",
    "\n",
    "def gypsum_pixel_count_landsat_GSLNA(image):\n",
    "    area_image = ee.Image.pixelArea()\n",
    "    mask = image.select('gypsum').gte(landsat_gypsum_threshold)\n",
    "    final = image.addBands(area_image)\n",
    "    stats = final.select('area').updateMask(mask).rename('gypsum').reduceRegion(\n",
    "        reducer = ee.Reducer.sum(),\n",
    "        geometry= GSL_NA_boundary,\n",
    "        scale=30,\n",
    "        maxPixels = 1e12)       \n",
    "    return image.set('gypsum', stats.get('gypsum'))\n",
    "\n",
    "def gypsum_pixel_count_landsat_BRB(image):\n",
    "    area_image = ee.Image.pixelArea()\n",
    "    mask = image.select('gypsum').gte(landsat_gypsum_threshold)\n",
    "    final = image.addBands(area_image)\n",
    "    stats = final.select('area').updateMask(mask).rename('gypsum').reduceRegion(\n",
    "        reducer = ee.Reducer.sum(),\n",
    "        geometry= BRB_boundary,\n",
    "        scale=30,\n",
    "        maxPixels = 1e12)       \n",
    "    return image.set('gypsum', stats.get('gypsum'))\n",
    "\n",
    "def gypsum_pixel_count_landsat_GSLSA(image):\n",
    "    area_image = ee.Image.pixelArea()\n",
    "    mask = image.select('gypsum').gte(landsat_gypsum_threshold)\n",
    "    final = image.addBands(area_image)\n",
    "    stats = final.select('area').updateMask(mask).rename('gypsum').reduceRegion(\n",
    "        reducer = ee.Reducer.sum(),\n",
    "        geometry= GSL_SA_boundary,\n",
    "        scale=30,\n",
    "        maxPixels = 1e12)        \n",
    "    return image.set('gypsum', stats.get('gypsum'))\n",
    "\n",
    "def carbonate_pixel_count_landsat_GSLNA(image):\n",
    "    area_image = ee.Image.pixelArea()\n",
    "    mask = image.select('carb_muds').lt(landsat_gypsum_threshold)\n",
    "    final = image.addBands(area_image)\n",
    "    stats = final.select('area').updateMask(mask).rename('carb_muds').reduceRegion(\n",
    "        reducer = ee.Reducer.sum(),\n",
    "        geometry= GSL_NA_boundary,\n",
    "        scale=30,\n",
    "        maxPixels = 1e12)       \n",
    "    return image.set('carb_muds', stats.get('carb_muds'))\n",
    "\n",
    "def carbonate_pixel_count_landsat_BRB(image):\n",
    "    area_image = ee.Image.pixelArea()\n",
    "    mask = image.select('carb_muds').lt(landsat_gypsum_threshold)\n",
    "    final = image.addBands(area_image)\n",
    "    stats = final.select('area').updateMask(mask).rename('carb_muds').reduceRegion(\n",
    "        reducer = ee.Reducer.sum(),\n",
    "        geometry= BRB_boundary,\n",
    "        scale=30,\n",
    "        maxPixels = 1e12)       \n",
    "    return image.set('carb_muds', stats.get('carb_muds'))\n",
    "\n",
    "def carbonate_pixel_count_landsat_GSLSA(image):\n",
    "    area_image = ee.Image.pixelArea()\n",
    "    mask = image.select('carb_muds').lt(landsat_gypsum_threshold)\n",
    "    final = image.addBands(area_image)\n",
    "    stats = final.select('area').updateMask(mask).rename('carb_muds').reduceRegion(\n",
    "        reducer = ee.Reducer.sum(),\n",
    "        geometry= GSL_SA_boundary,\n",
    "        scale=30,\n",
    "        maxPixels = 1e12)        \n",
    "    return image.set('carb_muds', stats.get('carb_muds'))\n",
    "\n",
    "def FLHV_pixel_count_landsat_GSLNA(image):\n",
    "    stats = image.select('FLHV').reduceRegion(\n",
    "        reducer = ee.Reducer.count(),\n",
    "        geometry= GSL_NA_boundary,\n",
    "        scale=30,\n",
    "        maxPixels=90000000)       \n",
    "    return image.set(stats)\n",
    "\n",
    "def FLHV_pixel_count_landsat_GSLSA(image):\n",
    "    stats = image.select('FLHV').reduceRegion(\n",
    "        reducer = ee.Reducer.count(),\n",
    "        geometry= GSL_SA_boundary,\n",
    "        scale=30,\n",
    "        maxPixels=90000000)       \n",
    "    return image.set(stats)\n",
    "\n",
    "def FLHV_avg_landsat_GSLNA(image):\n",
    "    stats = image.select('FLHV').reduceRegion(\n",
    "        reducer = ee.Reducer.mean(),\n",
    "        geometry= GSL_NA_boundary,\n",
    "        scale=30,\n",
    "        maxPixels=90000000)       \n",
    "    return image.set(stats)\n",
    "\n",
    "def FLHV_avg_landsat_GSLSA(image):\n",
    "    stats = image.select('FLHV').reduceRegion(\n",
    "        reducer = ee.Reducer.mean(),\n",
    "        geometry= GSL_SA_boundary,\n",
    "        scale=30,\n",
    "        maxPixels=90000000)       \n",
    "    return image.set(stats)\n",
    "\n",
    "def total_pixel_count_landsat_GSLNA(image):\n",
    "    area_image = ee.Image.pixelArea()\n",
    "    mask = ee.Image.constant(1).clip(NA_SaltPool).mask().eq(0)\n",
    "    final = image.addBands(area_image)\n",
    "    stats = final.select('area').updateMask(mask).rename('SumArea').reduceRegion(\n",
    "        reducer = ee.Reducer.sum(),\n",
    "        geometry= GSL_NA_boundary,\n",
    "        scale=30,\n",
    "        maxPixels = 1e12)       \n",
    "    return image.set('SumArea', stats.get('SumArea'))\n",
    "\n",
    "def total_pixel_count_landsat_BRB(image):\n",
    "    area_image = ee.Image.pixelArea()\n",
    "    mask = ee.Image.constant(1).mask().eq(0)\n",
    "    final = image.addBands(area_image)\n",
    "    stats = final.select('area').rename('SumArea').reduceRegion(\n",
    "        reducer = ee.Reducer.sum(),\n",
    "        geometry= BRB_boundary,\n",
    "        scale=30,\n",
    "        maxPixels = 1e12)       \n",
    "    return image.set('SumArea', stats.get('SumArea'))\n",
    "\n",
    "def total_pixel_count_landsat_GSLSA(image):\n",
    "    area_image = ee.Image.pixelArea()\n",
    "    final = image.addBands(area_image)\n",
    "    stats = final.select('area').rename('SumArea').reduceRegion(\n",
    "        reducer = ee.Reducer.sum(),\n",
    "        geometry= GSL_SA_boundary,\n",
    "        scale=30,\n",
    "        maxPixels = 1e12)         \n",
    "    return image.set('SumArea', stats.get('SumArea'))\n",
    "\n",
    "def total_pixel_count_sentinel_GSLNA(image):\n",
    "    area_image = ee.Image.pixelArea()\n",
    "    mask = ee.Image.constant(1).clip(NA_SaltPool).mask().eq(0)\n",
    "    final = image.addBands(area_image)\n",
    "    stats = final.select('area').updateMask(mask).rename('SumArea').reduceRegion(\n",
    "        reducer = ee.Reducer.sum(),\n",
    "        geometry= GSL_NA_boundary,\n",
    "        scale=30,\n",
    "        maxPixels = 1e12)       \n",
    "    return image.set('SumArea', stats.get('SumArea'))\n",
    "\n",
    "def total_pixel_count_sentinel_BRB(image):\n",
    "    area_image = ee.Image.pixelArea()\n",
    "    mask = ee.Image.constant(1).mask().eq(0)\n",
    "    final = image.addBands(area_image)\n",
    "    stats = final.select('area').updateMask(mask).rename('SumArea').reduceRegion(\n",
    "        reducer = ee.Reducer.sum(),\n",
    "        geometry= BRB_boundary,\n",
    "        scale=30,\n",
    "        maxPixels = 1e12)       \n",
    "    return image.set('SumArea', stats.get('SumArea'))\n",
    "\n",
    "def total_pixel_count_sentinel_GSLSA(image):\n",
    "    area_image = ee.Image.pixelArea()\n",
    "    final = image.addBands(area_image)\n",
    "    stats = final.select('area').updateMask(image.select('B3').gt(0)).rename('SumArea').reduceRegion(\n",
    "        reducer = ee.Reducer.sum(),\n",
    "        geometry= GSL_SA_boundary,\n",
    "        scale=30,\n",
    "        maxPixels = 1e12)         \n",
    "    return image.set('SumArea', stats.get('SumArea'))\n",
    "\n",
    "def halite_pixel_count_sent_SA(image):\n",
    "    area_image = ee.Image.pixelArea()\n",
    "    mask = image.select('halite').gte(sentinel_halite_threshold)\n",
    "    final = image.addBands(area_image)\n",
    "    stats = final.select('area').updateMask(mask).rename('halite').reduceRegion(\n",
    "        reducer = ee.Reducer.sum(),\n",
    "        geometry= GSL_SA_boundary,\n",
    "        scale=10,\n",
    "        maxPixels = 1e12)      \n",
    "    return image.set('halite', stats.get('halite'))\n",
    "\n",
    "def halite_pixel_count_sent_NA(image):\n",
    "    area_image = ee.Image.pixelArea()\n",
    "    mask = image.select('halite').gte(sentinel_halite_threshold)\n",
    "    final = image.addBands(area_image)\n",
    "    stats = final.select('area').updateMask(mask).rename('halite').reduceRegion(\n",
    "        reducer = ee.Reducer.sum(),\n",
    "        geometry= GSL_NA_boundary,\n",
    "        scale=10,\n",
    "        maxPixels = 1e12)       \n",
    "    return image.set('halite', stats.get('halite'))\n",
    "\n",
    "def halite_pixel_count_sent_BRB(image):\n",
    "    area_image = ee.Image.pixelArea()\n",
    "    mask = image.select('halite').gte(sentinel_halite_threshold)\n",
    "    final = image.addBands(area_image)\n",
    "    stats = final.select('area').updateMask(mask).rename('halite').reduceRegion(\n",
    "        reducer = ee.Reducer.sum(),\n",
    "        geometry= BRB_boundary,\n",
    "        scale=10,\n",
    "        maxPixels = 1e12)       \n",
    "    return image.set('halite', stats.get('halite'))\n",
    "\n",
    "def carbonate_pixel_count_sent_SA(image):\n",
    "    area_image = ee.Image.pixelArea()\n",
    "    mask = image.select('carb_muds').lt(sentinel_gypsum_threshold)\n",
    "    final = image.addBands(area_image)\n",
    "    stats = final.select('area').updateMask(mask).rename('carb_muds').reduceRegion(\n",
    "        reducer = ee.Reducer.sum(),\n",
    "        geometry= GSL_SA_boundary,\n",
    "        scale=10,\n",
    "        maxPixels = 1e12)        \n",
    "    return image.set('carb_muds', stats.get('carb_muds'))\n",
    "\n",
    "def carbonate_pixel_count_sent_NA(image):\n",
    "    area_image = ee.Image.pixelArea()\n",
    "    mask = image.select('carb_muds').lt(sentinel_gypsum_threshold)\n",
    "    final = image.addBands(area_image)\n",
    "    stats = final.select('area').updateMask(mask).rename('carb_muds').reduceRegion(\n",
    "        reducer = ee.Reducer.sum(),\n",
    "        geometry= GSL_NA_boundary,\n",
    "        scale=10,\n",
    "        maxPixels = 1e12)      \n",
    "    return image.set('carb_muds', stats.get('carb_muds'))\n",
    "\n",
    "def carbonate_pixel_count_sent_BRB(image):\n",
    "    area_image = ee.Image.pixelArea()\n",
    "    mask = image.select('carb_muds').lt(sentinel_gypsum_threshold)\n",
    "    final = image.addBands(area_image)\n",
    "    stats = final.select('area').updateMask(mask).rename('carb_muds').reduceRegion(\n",
    "        reducer = ee.Reducer.sum(),\n",
    "        geometry= BRB_boundary,\n",
    "        scale=10,\n",
    "        maxPixels = 1e12)      \n",
    "    return image.set('carb_muds', stats.get('carb_muds'))\n",
    "\n",
    "def gypsum_pixel_count_sent_SA(image):\n",
    "    area_image = ee.Image.pixelArea()\n",
    "    mask = image.select('gypsum').gte(sentinel_gypsum_threshold)\n",
    "    final = image.addBands(area_image)\n",
    "    stats = final.select('area').updateMask(mask).rename('gypsum').reduceRegion(\n",
    "        reducer = ee.Reducer.sum(),\n",
    "        geometry= GSL_SA_boundary,\n",
    "        scale=10,\n",
    "        maxPixels = 1e12)      \n",
    "    return image.set('gypsum', stats.get('gypsum'))\n",
    "\n",
    "def gypsum_pixel_count_sent_NA(image):\n",
    "    area_image = ee.Image.pixelArea()\n",
    "    mask = image.select('gypsum').gte(sentinel_gypsum_threshold)\n",
    "    final = image.addBands(area_image)\n",
    "    stats = final.select('area').updateMask(mask).rename('gypsum').reduceRegion(\n",
    "        reducer = ee.Reducer.sum(),\n",
    "        geometry= GSL_NA_boundary,\n",
    "        scale=10,\n",
    "        maxPixels = 1e12)       \n",
    "    return image.set('gypsum', stats.get('gypsum'))\n",
    "\n",
    "def gypsum_pixel_count_sent_BRB(image):\n",
    "    area_image = ee.Image.pixelArea()\n",
    "    mask = image.select('gypsum').gte(sentinel_gypsum_threshold)\n",
    "    final = image.addBands(area_image)\n",
    "    stats = final.select('area').updateMask(mask).rename('gypsum').reduceRegion(\n",
    "        reducer = ee.Reducer.sum(),\n",
    "        geometry= BRB_boundary,\n",
    "        scale=10,\n",
    "        maxPixels = 1e12)       \n",
    "    return image.set('gypsum', stats.get('gypsum'))\n",
    "\n",
    "def BDA_pixel_count_sent_SA(image):\n",
    "    stats = image.select('2BDA').reduceRegion(\n",
    "        reducer = ee.Reducer.count(),\n",
    "        geometry= GSL_SA_boundary,\n",
    "        scale=20,\n",
    "        maxPixels = 1e12)       \n",
    "    return image.set(stats)\n",
    "\n",
    "def BDA_pixel_count_sent_NA(image):\n",
    "    stats = image.select('2BDA').reduceRegion(\n",
    "        reducer = ee.Reducer.count(),\n",
    "        geometry= GSL_NA_boundary,\n",
    "        scale=20,\n",
    "        maxPixels = 1e12)       \n",
    "    return image.set(stats)\n",
    "\n",
    "def BDA_average_sent_NA(image):\n",
    "    stats = image.select('2BDA').reduceRegion(\n",
    "        reducer = ee.Reducer.mean(),\n",
    "        geometry= GSL_NA_boundary,\n",
    "        scale=20,\n",
    "        maxPixels = 1e12)       \n",
    "    return image.set(stats)\n",
    "\n",
    "def BDA_average_sent_SA(image):\n",
    "    stats = image.select('2BDA').reduceRegion(\n",
    "        reducer = ee.Reducer.mean(),\n",
    "        geometry= GSL_SA_boundary,\n",
    "        scale=20,\n",
    "        maxPixels = 1e12)       \n",
    "    return image.set(stats)\n",
    "\n",
    "def OtsuThreshold(histogram):\n",
    "    counts = ee.Array(ee.Dictionary(histogram).get('histogram'))\n",
    "    means = ee.Array(ee.Dictionary(histogram).get('bucketMeans'))\n",
    "    size = means.length().get([0])\n",
    "    total = counts.reduce(ee.Reducer.sum(), [0]).get([0])\n",
    "    sum = means.multiply(counts).reduce(ee.Reducer.sum(), [0]).get([0])\n",
    "    mean = sum.divide(total)\n",
    "    indices = ee.List.sequence(1, size)\n",
    "\n",
    "    def func_xxx(i):\n",
    "        aCounts = counts.slice(0, 0, i)\n",
    "        aCount = aCounts.reduce(ee.Reducer.sum(), [0]).get([0])\n",
    "        aMeans = means.slice(0, 0, i)\n",
    "        aMean = (\n",
    "            aMeans.multiply(aCounts)\n",
    "            .reduce(ee.Reducer.sum(), [0])\n",
    "            .get([0])\n",
    "            .divide(aCount)\n",
    "        )\n",
    "        bCount = total.subtract(aCount)\n",
    "        bMean = sum.subtract(aCount.multiply(aMean)).divide(bCount)\n",
    "        return aCount.multiply(aMean.subtract(mean).pow(2)).add(\n",
    "            bCount.multiply(bMean.subtract(mean).pow(2)))\n",
    "\n",
    "    bss = indices.map(func_xxx)\n",
    "    return means.sort(bss).get([-1])\n",
    "\n",
    "def extract_water_NA_ls(image):\n",
    "    histogram = image.select('ndwi').reduceRegion(\n",
    "        reducer = ee.Reducer.histogram(255, 2),\n",
    "        geometry = GSL_NA_boundary.geometry().buffer(6000),\n",
    "        scale = 30,\n",
    "        bestEffort= True,)\n",
    "    threshold = OtsuThreshold(histogram.get('ndwi')).add(0.15) \n",
    "    water = image.select('ndwi').updateMask(image.select('ndwi').gte(threshold)).rename('ndwi').copyProperties(image)\n",
    "    return water.set({\"threshold\": threshold}).copyProperties(image)\n",
    "\n",
    "def extract_water_SA_ls(image):\n",
    "    histogram = image.select('ndwi').reduceRegion(\n",
    "        reducer = ee.Reducer.histogram(255, 2),\n",
    "        geometry = GSL_SA_boundary.geometry().buffer(6000),\n",
    "        scale = 30,\n",
    "        bestEffort= True,)\n",
    "    threshold = OtsuThreshold(histogram.get('ndwi')).add(0.175)\n",
    "    water = image.select('ndwi').updateMask(image.select('ndwi').gte(threshold)).rename('ndwi').copyProperties(image)\n",
    "    return water.set({\"threshold\": threshold})\n",
    "\n",
    "def extract_water_BRB_ls(image):\n",
    "    histogram = image.select('ndwi').reduceRegion(\n",
    "        reducer = ee.Reducer.histogram(255, 2),\n",
    "        geometry = BRB_boundary.geometry().buffer(6000),\n",
    "        scale = 30,\n",
    "        bestEffort= True,)\n",
    "    threshold = OtsuThreshold(histogram.get('ndwi')).add(0.175) \n",
    "    water = image.select('ndwi').updateMask(image.select('ndwi').gte(threshold)).rename('ndwi').copyProperties(image)\n",
    "    return water.set({\"threshold\": threshold})\n",
    "\n",
    "def GSL_NA_halite_watermask(image):\n",
    "    histogram = image.select('ndwi').reduceRegion(\n",
    "        reducer = ee.Reducer.histogram(255, 2),\n",
    "        geometry = GSL_NA_boundary.geometry().buffer(6000),\n",
    "        scale = 30,\n",
    "        bestEffort= True,)\n",
    "    threshold = OtsuThreshold(histogram.get('ndwi')).add(0.15) \n",
    "    mask = image.select('ndwi').lte(threshold)\n",
    "    return image.updateMask(mask)\n",
    "\n",
    "def BRB_halite_watermask(image):\n",
    "    histogram = image.select('ndwi').reduceRegion(\n",
    "        reducer = ee.Reducer.histogram(255, 2),\n",
    "        geometry = BRB_boundary.geometry().buffer(6000),\n",
    "        scale = 30,\n",
    "        bestEffort= True,)\n",
    "    threshold = OtsuThreshold(histogram.get('ndwi')).add(0.175) \n",
    "    mask = image.select('ndwi').lte(threshold)\n",
    "    return image.updateMask(mask)\n",
    "\n",
    "def GSL_SA_halite_watermask(image):\n",
    "    histogram = image.select('ndwi').reduceRegion(\n",
    "        reducer = ee.Reducer.histogram(255, 2),\n",
    "        geometry = GSL_SA_boundary.geometry().buffer(6000),\n",
    "        scale = 30,\n",
    "        bestEffort= True,)\n",
    "    threshold = OtsuThreshold(histogram.get('ndwi')).add(0.175)\n",
    "    mask = image.select('ndwi').lte(threshold)\n",
    "    return image.updateMask(mask)\n",
    "\n",
    "sentinel_ndwi_threshold = 0.06\n",
    "\n",
    "def GSL_S2_watermask(image):\n",
    "    mask = image.select('ndwi').lte(sentinel_ndwi_threshold) \n",
    "    return image.updateMask(mask)    \n",
    "\n",
    "def GSL_SA_S2_watermask(image):\n",
    "    mask = image.select('ndwi').lte(0.04) \n",
    "    return image.updateMask(mask) \n",
    "\n",
    "def BRB_S2_watermask(image):\n",
    "    mask = image.select('ndwi').lte(0.04) \n",
    "    return image.updateMask(mask)\n",
    "\n",
    "def extract_water_NA_S2(image):\n",
    " \n",
    "    water = image.select('ndwi').updateMask(image.select('ndwi').gte(sentinel_ndwi_threshold)).rename('ndwi').copyProperties(image) \n",
    "    return water\n",
    "\n",
    "def extract_water_SA_S2(image):\n",
    "\n",
    "    water = image.select('ndwi').updateMask(image.select('ndwi').gte(sentinel_ndwi_threshold)).rename('ndwi').copyProperties(image) \n",
    "    return water\n",
    "\n",
    "def extract_water_BRB_S2(image):\n",
    "\n",
    "    water = image.select('ndwi').updateMask(image.select('ndwi').gte(sentinel_ndwi_threshold)).rename('ndwi').copyProperties(image) \n",
    "    return water\n",
    "\n",
    "def CollectionStitch(col1, col2, dates_list): #this function mosaics north and south images only if their dates match. Ignores scenes without a partner.\n",
    "    image_list = []\n",
    "    for date in dates_list:\n",
    "        if date in col1.aggregate_array('Date_Filter').getInfo() and date in col2.aggregate_array('Date_Filter').getInfo(): #conditional statement to ensure this will only run if there are two images with matching dates. Avoids errors\n",
    "            filtered_col1 = image_grab(col1, date)\n",
    "            filtered_col2 = image_grab(col2, date)\n",
    "            merged_col = ee.ImageCollection.fromImages([filtered_col1, filtered_col2])\n",
    "            mosaic = merged_col.mosaic().copyProperties(filtered_col1) #new collection images contain all image properties of the northern landsat image\n",
    "            image_list.append(mosaic)\n",
    "        else: None #If the condition isnt met, do nothing and keep going through the list\n",
    "    new_col = ee.ImageCollection.fromImages(image_list)\n",
    "    return new_col\n",
    "\n",
    "def temperature_bands(img):\n",
    "    scale1 = ['ST_ATRAN', 'ST_EMIS']\n",
    "    scale2 = ['ST_DRAD', 'ST_TRAD', 'ST_URAD']\n",
    "    scale1_names = ['transmittance', 'emissivity']\n",
    "    scale2_names = ['downwelling', 'B10_radiance', 'upwelling']\n",
    "    scale1_bands = img.select(scale1).multiply(0.0001).rename(scale1_names) #Scaled to new L8 collection\n",
    "    scale2_bands = img.select(scale2).multiply(0.001).rename(scale2_names) #Scaled to new L8 collection\n",
    "    return img.addBands(scale1_bands).addBands(scale2_bands).copyProperties(img)\n",
    "\n",
    "def landsat_LST(image):\n",
    "    # Based on Sekertekin, A., & Bonafoni, S. (2020) https://doi.org/10.3390/rs12020294\n",
    "    \n",
    "    k1 = 774.89\n",
    "    k2 = 1321.08\n",
    "    LST = image.expression(\n",
    "        '(k2/log((k1/((B10_rad - upwelling - transmittance*(1 - emissivity)*downwelling)/(transmittance*emissivity)))+1)) - 273.15',\n",
    "        {'k1': k1,\n",
    "        'k2': k2,\n",
    "        'B10_rad': image.select('B10_radiance'),\n",
    "        'upwelling': image.select('upwelling'),\n",
    "        'transmittance': image.select('transmittance'),\n",
    "        'emissivity': image.select('emissivity'),\n",
    "        'downwelling': image.select('downwelling')}).rename('LST')\n",
    "    return image.addBands(LST).copyProperties(image) \n",
    "\n",
    "def average_GSL_NA_temperature(image):\n",
    "    histogram = image.select('ndwi').reduceRegion(\n",
    "        reducer = ee.Reducer.histogram(255, 2),\n",
    "        geometry = GSL_NA_boundary.geometry().buffer(6000),\n",
    "        scale = 30,\n",
    "        bestEffort= True,)\n",
    "    threshold = OtsuThreshold(histogram.get('ndwi')).add(0.15) \n",
    "    water = image.select('ndwi').updateMask(image.select('ndwi').gte(threshold))\n",
    "    stats = image.select('LST').updateMask(image.select('ndwi').gte(threshold)).reduceRegion(\n",
    "    reducer = ee.Reducer.mean(),\n",
    "    geometry = GSL_NA_boundary,\n",
    "    scale=30)  \n",
    "    return image.updateMask(image.select('ndwi').gte(threshold)).copyProperties(image).set(stats) \n",
    "\n",
    "def average_GSL_SA_temperature(image):\n",
    "    histogram = image.select('ndwi').reduceRegion(\n",
    "        reducer = ee.Reducer.histogram(255, 2),\n",
    "        geometry = GSL_SA_boundary.geometry().buffer(6000),\n",
    "        scale = 30,\n",
    "        bestEffort= True,)\n",
    "    threshold = OtsuThreshold(histogram.get('ndwi')).add(0.175) \n",
    "    water = image.select('ndwi').updateMask(image.select('ndwi').gte(threshold)) \n",
    "    stats = image.select('LST').updateMask(image.select('ndwi').gte(threshold)).reduceRegion(\n",
    "    reducer = ee.Reducer.mean(),\n",
    "    geometry = GSL_SA_boundary,\n",
    "    scale=30)  \n",
    "    return image.updateMask(image.select('ndwi').gte(threshold)).copyProperties(image).set(stats)\n",
    "def average_GSL_waterstation_temperature(image):\n",
    "    histogram = image.select('ndwi').reduceRegion(\n",
    "        reducer = ee.Reducer.histogram(255, 2),\n",
    "        geometry = GSL_SA_boundary.geometry().buffer(6000),\n",
    "        scale = 30,\n",
    "        bestEffort= True,)\n",
    "    threshold = OtsuThreshold(histogram.get('ndwi')).add(0.175)\n",
    "    water = image.select('ndwi').updateMask(image.select('ndwi').gte(threshold)) \n",
    "    stats = image.select('LST').updateMask(image.select('ndwi').gte(threshold)).reduceRegion(\n",
    "    reducer = ee.Reducer.mean(),\n",
    "    geometry = USGS_water_station_region,\n",
    "    scale=30)  \n",
    "    return image.set(stats).updateMask(image.select('ndwi').gte(threshold))\n",
    "\n",
    "def MaskToCarbonatesL8(image):\n",
    "    gypsum = image.select('gypsum')\n",
    "    halite = image.select('halite')\n",
    "    halite_index = image.normalizedDifference(['SR_B4', 'SR_B6'])\n",
    "    gypsum_index = image.normalizedDifference(['SR_B6', 'SR_B7'])\n",
    "    mask1 = gypsum_index.updateMask(halite_index.lt(landsat_halite_threshold)).updateMask(gypsum_index.lt(landsat_gypsum_threshold)).rename('carb_muds').copyProperties(image)\n",
    "    return mask1\n",
    "\n",
    "def MaskToCarbonatesS2(image):\n",
    "    halite_index = image.normalizedDifference(['B4', 'B11'])\n",
    "    gypsum_index = image.normalizedDifference(['B11', 'B12'])\n",
    "    mask1 = gypsum_index.updateMask(halite_index.lt(sentinel_halite_threshold)).updateMask(gypsum_index.lt(sentinel_gypsum_threshold)).rename('carb_muds').copyProperties(image)\n",
    "    return mask1\n",
    "\n",
    "############# CHLOROPHYLL SPECTRAL INDICES #############\n",
    "\n",
    "### Sentinel 2 Band Algorithm adapted from Dall’Olmo, G., & Gitelson, A. A. (2006). Effect of bio-optical parameter variability and uncertainties \n",
    "#in reflectance measurements on the remote estimation of chlorophyll-a concentration in turbid productive waters: \n",
    "# Modeling results. Applied Optics, 45(15), 3577. https://doi.org/10.1364/AO.45.003577\n",
    "def landsat_MaskToWater_NA(image):\n",
    "    histogram = image.select('ndwi').reduceRegion(\n",
    "    reducer = ee.Reducer.histogram(255, 2),\n",
    "    geometry = GSL_NA_boundary.geometry().buffer(6000),\n",
    "    scale = 30,\n",
    "    bestEffort= True,)\n",
    "    threshold = OtsuThreshold(histogram.get('ndwi')).add(0.15) \n",
    "    mask = image.select('ndwi').gte(threshold)\n",
    "    return image.updateMask(mask).copyProperties(image)\n",
    "\n",
    "def landsat_MaskToWater_SA(image):\n",
    "    histogram = image.select('ndwi').reduceRegion(\n",
    "    reducer = ee.Reducer.histogram(255, 2),\n",
    "    geometry = GSL_SA_boundary.geometry().buffer(6000),\n",
    "    scale = 30,\n",
    "    bestEffort= True,)\n",
    "    threshold = OtsuThreshold(histogram.get('ndwi')).add(0.175) \n",
    "    mask = image.select('ndwi').gte(threshold)\n",
    "    return image.updateMask(mask).copyProperties(image)\n",
    "\n",
    "def landsat_MaskToWater_BRB(image):\n",
    "    histogram = image.select('ndwi').reduceRegion(\n",
    "    reducer = ee.Reducer.histogram(255, 2),\n",
    "    geometry = BRB_boundary.geometry().buffer(6000),\n",
    "    scale = 30,\n",
    "    bestEffort= True,)\n",
    "    threshold = OtsuThreshold(histogram.get('ndwi')).add(0.175) \n",
    "    mask = image.select('ndwi').gte(threshold)\n",
    "    return image.updateMask(mask).copyProperties(image)\n",
    "\n",
    "def sentinel_MaskToWater(image):\n",
    "    water = image.updateMask(image.select('ndwi').gte(sentinel_ndwi_threshold)).copyProperties(image) \n",
    "    return water\n",
    "\n",
    "def landsat_FLHV(image):\n",
    "    FLHV = image.expression(\n",
    "    '(GREEN - (RED + (CA-RED)))', {\n",
    "        'CA': image.select('SR_B1'),\n",
    "        'GREEN': image.select('SR_B3'),\n",
    "        'RED': image.select('SR_B4')})\n",
    "    mask = FLHV.updateMask(FLHV.gte(3500)).rename('FLHV').copyProperties(image)\n",
    "    return mask\n",
    "\n",
    "def landsat_KIVU(image):\n",
    "    KIVU = image.expression(\n",
    "    '(BLUE - RED) / GREEN', {\n",
    "        'BLUE': image.select('SR_B2'),\n",
    "        'GREEN': image.select('SR_B3'),\n",
    "        'RED': image.select('SR_B4')})\n",
    "    mask = KIVU.updateMask(KIVU.gte(-1)).rename('FLHV').copyProperties(image) \n",
    "    return mask\n",
    "\n",
    "def sentinel_2BDA(image):\n",
    "    S2BDA = image.expression(\n",
    "    '(RED_EDGE / RED)', {\n",
    "        'RED': image.select('B4'),\n",
    "        'RED_EDGE': image.select('B5')})\n",
    "    mask = S2BDA.updateMask(S2BDA.gte(-1)).rename('2BDA').copyProperties(image)\n",
    "    return mask\n",
    "\n",
    "def sentinel_veg_mask(image):\n",
    "    ndvi_calc = image.normalizedDifference(['B8', 'B4']) #NIR-RED/NIR+RED -- full NDVI image\n",
    "    mask = image.updateMask(ndvi_calc.lt(sentinel_ndvi_threshold)).copyProperties(image) \n",
    "    return mask\n",
    "\n",
    "def landsat_veg_mask(image):\n",
    "    ndvi_calc = image.normalizedDifference(['SR_B5', 'SR_B4']) #NIR-RED/NIR+RED -- full NDVI image\n",
    "    mask = image.updateMask(ndvi_calc.lt(landsat_ndvi_threshold)).copyProperties(image) \n",
    "    return mask\n",
    "\n",
    "def landsat_halite_mask(image):\n",
    "    halite_index = image.normalizedDifference(['SR_B4', 'SR_B6']) # red-swir1 / red+swir1\n",
    "    mask = image.updateMask(halite_index.lt(landsat_halite_threshold)).copyProperties(image) \n",
    "    return mask \n",
    "\n",
    "def sentinel_halite_mask(image):\n",
    "    halite_index = image.normalizedDifference(['B4', 'B11']) # red-swir1 / red+swir1\n",
    "    mask = image.updateMask(halite_index.lt(sentinel_halite_threshold)).copyProperties(image) \n",
    "    return mask "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "### Tile definition ###\n",
    "#Great Salt Lake\n",
    "SLC_S_tile_N = '12TUM' #sentinel\n",
    "SLC_S_tile_S = '12TUL' #sentinel\n",
    "SLC_L_N_row = 31 #landsat\n",
    "SLC_L_S_row = 32 #landsat\n",
    "SLC_L_path = 38 #landsat\n",
    "# Bear River Bay\n",
    "BRB_tile_row = 31\n",
    "BRB_tile_path = 38\n",
    "#Bonneville basin\n",
    "BSF_L_path = 39 #landsat\n",
    "BSF_L_row_N = 31 #landsat\n",
    "BSF_L_row_S = 32 #landsat\n",
    "BSF_S_tile_N = '12TTL' #sentinel\n",
    "BSF_S_tile_S = '12TTK' #sentinel\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [],
   "source": [
    "### Landsat dataset definition and filtering ###\n",
    "#Landsat 5 definition and band renaming\n",
    "landsat5 = ee.ImageCollection(\"LANDSAT/LT05/C02/T1_L2\").map(landsat5bandrename)\n",
    "\n",
    "# Note: used a sort function to sort all images to be in chronological order. Otherwise defaults to L8 images first, then L9 images out of sort.\n",
    "start_date = ee.Date(str('2022-11-20'))\n",
    "end_date = ee.Date(str('2023-05-30')) #date.today()))\n",
    "### Great Salt Lake ###\n",
    "#Landsat 5, 8, & 9 -- spanning from 1984-3-16 to 2012-05-05 and 2013-03-18 to present\n",
    "GSL_landsat_N = ee.ImageCollection(\"LANDSAT/LC08/C02/T1_L2\").merge(ee.ImageCollection(\"LANDSAT/LC09/C02/T1_L2\")).merge(landsat5).filterDate(start_date \\\n",
    "    .format('YYYY-MM-dd'), end_date.format('YYYY-MM-dd')).filter(ee.Filter.And(ee.Filter.eq('WRS_PATH', SLC_L_path), \\\n",
    "        ee.Filter.eq('WRS_ROW', SLC_L_N_row))).filter(ee.Filter.lte('CLOUD_COVER', 10)).map(image_dater).sort('Date_Filter') #.map(maskL8clouds)\n",
    "\n",
    "GSL_landsat_S = ee.ImageCollection(\"LANDSAT/LC08/C02/T1_L2\").merge(ee.ImageCollection(\"LANDSAT/LC09/C02/T1_L2\")).merge(landsat5).filterDate(start_date \\\n",
    "    .format('YYYY-MM-dd'), end_date.format('YYYY-MM-dd')).filter(ee.Filter.And(ee.Filter.eq('WRS_PATH', SLC_L_path), \\\n",
    "        ee.Filter.eq('WRS_ROW', SLC_L_S_row))).filter(ee.Filter.lte('CLOUD_COVER', 10)).map(image_dater).sort('Date_Filter') #.map(maskL8clouds)\n",
    "\n",
    "GSL_landsat_dates_N = GSL_landsat_N.aggregate_array('Date_Filter').getInfo()\n",
    "GSL_landsat_dates_S = GSL_landsat_S.aggregate_array('Date_Filter').getInfo()\n",
    "GSL_landsat_MR = GSL_landsat_N.aggregate_array('Date_Filter').getInfo()[-1] #This is the date of image we want to grab\n",
    "\n",
    "\n",
    "# bear river bay\n",
    "BRB_landsat = ee.ImageCollection(\"LANDSAT/LC08/C02/T1_L2\").merge(ee.ImageCollection(\"LANDSAT/LC09/C02/T1_L2\")).merge(landsat5).filterDate(start_date \\\n",
    "    .format('YYYY-MM-dd'), end_date.format('YYYY-MM-dd')).filter(ee.Filter.And(ee.Filter.eq('WRS_PATH', BRB_tile_path), \\\n",
    "        ee.Filter.eq('WRS_ROW', BRB_tile_row))).filter(ee.Filter.lte('CLOUD_COVER', 10)).map(image_dater).sort('Date_Filter') #.map(maskL8clouds)\n",
    "\n",
    "BRB_landsat_dates = BRB_landsat.aggregate_array('Date_Filter').getInfo()\n",
    "BRB_landsat_MR = BRB_landsat.aggregate_array('Date_Filter').getInfo()[-1] #This is the date of image we want to grab\n",
    "\n",
    "#Standard cloud cover has been 10%\n",
    "\n",
    "\n",
    "### Bonneville Salt Flats ###\n",
    "BSF_landsat_N = ee.ImageCollection(\"LANDSAT/LC08/C02/T1_L2\").merge(ee.ImageCollection(\"LANDSAT/LC09/C02/T1_L2\")).merge(landsat5).filterDate(start_date \\\n",
    "    .format('YYYY-MM-dd'), end_date.format('YYYY-MM-dd')).filter(ee.Filter.And(ee.Filter.eq('WRS_PATH', BSF_L_path), \\\n",
    "        ee.Filter.eq('WRS_ROW', BSF_L_row_N))).filter(ee.Filter.lte('CLOUD_COVER', 10)).map(image_dater).sort('Date_Filter')\n",
    "\n",
    "BSF_landsat_S = ee.ImageCollection(\"LANDSAT/LC08/C02/T1_L2\").merge(ee.ImageCollection(\"LANDSAT/LC09/C02/T1_L2\")).merge(landsat5).filterDate(start_date \\\n",
    "    .format('YYYY-MM-dd'), end_date.format('YYYY-MM-dd')).filter(ee.Filter.And(ee.Filter.eq('WRS_PATH', BSF_L_path), \\\n",
    "        ee.Filter.eq('WRS_ROW', BSF_L_row_S))).filter(ee.Filter.lte('CLOUD_COVER', 10)).map(image_dater).sort('Date_Filter')\n",
    "\n",
    "BSF_landsat_dates_N = BSF_landsat_N.aggregate_array('Date_Filter').getInfo()\n",
    "BSF_landsat_dates_S = BSF_landsat_S.aggregate_array('Date_Filter').getInfo()\n",
    "BSF_landsat_MR = BSF_landsat_N.aggregate_array('Date_Filter').getInfo()[-1]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [],
   "source": [
    "landsat_GSL_mosaic_col = CollectionStitch(GSL_landsat_N, GSL_landsat_S, GSL_landsat_dates_S)\n",
    "landsat_GSL_NA_mosaic_col = CollectionStitch(BSF_landsat_N, BSF_landsat_S, BSF_landsat_dates_N)\n",
    "\n",
    "landsat_GSL_mosaic_col_dates = landsat_GSL_mosaic_col.aggregate_array('Date_Filter').getInfo() #dates to use for visualization and to pair with each observation\n",
    "landsat_GSL_NA_mosaic_col_dates = landsat_GSL_NA_mosaic_col.aggregate_array('Date_Filter').getInfo()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "### total GSL area collections for reference ###\n",
    "GSL_NA_area = landsat_GSL_NA_mosaic_col.map(MaskSaltPool).map(NorthArmClip).map(total_pixel_count_landsat_GSLNA) #.map(MaskSaltPool)\n",
    "GSL_SA_area = landsat_GSL_mosaic_col.map(SouthArmClip).map(total_pixel_count_landsat_GSLSA)\n",
    "BRB_area = BRB_landsat.map(BRBClip).map(total_pixel_count_landsat_BRB)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "metadata": {},
   "outputs": [],
   "source": [
    "### Sentinel dataset definition and filtering ###\n",
    "start_date_s = ee.Date(str('2022-11-22'))\n",
    "end_date_s = ee.Date(str('2023-05-30')) \n",
    "#Great Salt Lake#\n",
    "#No data pixels of 10 and clouds of 50 set to parameters used for data collection. \n",
    "GSL_sentinel_N = ee.ImageCollection(\"COPERNICUS/S2_SR\").filter(ee.Filter.inList('MGRS_TILE', [SLC_S_tile_N])).filterDate(start_date_s.format('YYYY-MM-dd'), end_date_s.format('YYYY-MM-dd')).map(image_dater) \\\n",
    "    .select(['B12','B11', 'B8', 'B5', 'B4', 'B3', 'B2', 'QA60', 'SCL']).filter(ee.Filter.lte('NODATA_PIXEL_PERCENTAGE', 10)).filter(ee.Filter.lte('CLOUDY_PIXEL_PERCENTAGE', 10)) #.filter(ee.Filter.lte('NODATA_PIXEL_PERCENTAGE', 10)).filter(ee.Filter.lte('CLOUDY_PIXEL_PERCENTAGE', 20))\n",
    "\n",
    "GSL_sentinel_S = ee.ImageCollection(\"COPERNICUS/S2_SR\").filter(ee.Filter.inList('MGRS_TILE', [SLC_S_tile_S])).filterDate(start_date_s.format('YYYY-MM-dd'), end_date_s.format('YYYY-MM-dd')).map(image_dater) \\\n",
    "    .select(['B12', 'B11', 'B8', 'B5', 'B4', 'B3', 'B2', 'QA60', 'SCL']).filter(ee.Filter.lte('NODATA_PIXEL_PERCENTAGE', 10)) #.filter(ee.Filter.lte('CLOUDY_PIXEL_PERCENTAGE', 10)) #cloud filter provides too little images at lte to 5, so setting to 10\n",
    "GSL_sent_dates_N = GSL_sentinel_N.aggregate_array('Date_Filter').getInfo()\n",
    "GSL_sent_dates_S = GSL_sentinel_S.aggregate_array('Date_Filter').getInfo()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 28,
   "metadata": {},
   "outputs": [],
   "source": [
    "### Stitched sentinel collections ###\n",
    "sentinel_GSL_mosaic_col = CollectionStitch(GSL_sentinel_N, GSL_sentinel_S, GSL_sent_dates_N)\n",
    "\n",
    "sentinel_GSL_mosaic_col_dates = sentinel_GSL_mosaic_col.aggregate_array('Date_Filter').getInfo()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "### total GSL area collections for reference ###\n",
    "GSL_NA_area_S2 = sentinel_GSL_mosaic_col.map(NorthArmClip).map(total_pixel_count_sentinel_GSLNA)\n",
    "GSL_NA_area_S2_mask = sentinel_GSL_mosaic_col.map(NorthArmClip).map(MaskSaltPool).map(total_pixel_count_sentinel_GSLNA)\n",
    "GSL_SA_area_S2 = sentinel_GSL_mosaic_col.map(SouthArmClip).map(total_pixel_count_sentinel_GSLSA)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [],
   "source": [
    "### NDWI and LST collections ###\n",
    "### Landsat\n",
    "GSL_SA_NDWI = landsat_GSL_mosaic_col.map(GSLSA_WaterPixelCount).map(landsat_ndwi_fn).map(extract_water_SA_ls).map(SouthArmClip)\n",
    "GSL_NA_NDWI_pre93 = landsat_GSL_NA_mosaic_col.filter(ee.Filter.lt('Date_Filter', '1994-01-01')).map(GSLNA_WaterPixelCount).map(landsat_ndwi_fn).map(extract_water_NA_ls).map(NorthArmClip) #.map(MaskSaltPool).map(GSLNA_WaterPixelCount)\n",
    "GSL_NA_NDWI_post93 = landsat_GSL_NA_mosaic_col.filter(ee.Filter.gt('Date_Filter', '1994-01-01')).map(MaskSaltPool).map(GSLNA_WaterPixelCount).map(landsat_ndwi_fn).map(extract_water_NA_ls).map(NorthArmClip) #.map(MaskSaltPool).map(GSLNA_WaterPixelCount)\n",
    "GSL_NA_NDWI = GSL_NA_NDWI_pre93.merge(GSL_NA_NDWI_post93).sort('Date_Filter')\n",
    "BRB_NDWI = BRB_landsat.map(BRB_WaterPixelCount).map(landsat_ndwi_fn).map(extract_water_BRB_ls).map(BRBClip)\n",
    "\n",
    "unclipped_NA_NDWI = landsat_GSL_NA_mosaic_col.map(landsat_ndwi_fn).map(MaskSaltPool)\n",
    "unclipped_SA_NDWI = landsat_GSL_mosaic_col.map(landsat_ndwi_fn)\n",
    "unclipped_BRB_NDWI = BRB_landsat.map(landsat_ndwi_fn)\n",
    "\n",
    "toyFilter = ee.Filter.equals(leftField='Date_Filter', rightField='Date_Filter')\n",
    "innerJoin = ee.Join.inner()\n",
    "def joinedCol(image):\n",
    "    return ee.Image.cat(image.get('primary'), image.get('secondary'))\n",
    "\n",
    "combined_NA_NDWI_col = ee.ImageCollection(innerJoin.apply(GSL_NA_NDWI, landsat_GSL_NA_mosaic_col, toyFilter).map(joinedCol))\n",
    "combined_NA_NDWI_co_unclipped = ee.ImageCollection(innerJoin.apply(unclipped_NA_NDWI, landsat_GSL_NA_mosaic_col, toyFilter).map(joinedCol))\n",
    "\n",
    "GSL_NA_LST = combined_NA_NDWI_co_unclipped.map(temperature_bands).map(MaskSaltPool).map(landsat_LST).map(average_GSL_NA_temperature).map(NorthArmClip)\n",
    "\n",
    "#merging NDWI col with normal RGB col to use for clipping image when making LST map\n",
    "combined_SA_NDWI_col = ee.ImageCollection(innerJoin.apply(GSL_SA_NDWI, landsat_GSL_mosaic_col, toyFilter).map(joinedCol))\n",
    "combined_SA_NDWI_co_unclipped = ee.ImageCollection(innerJoin.apply(unclipped_SA_NDWI, landsat_GSL_mosaic_col, toyFilter).map(joinedCol))\n",
    "\n",
    "GSL_SA_LST = combined_SA_NDWI_co_unclipped.map(temperature_bands).map(landsat_LST).map(average_GSL_SA_temperature).map(SouthArmClip)\n",
    "GSL_SA_USGS_LST = combined_SA_NDWI_col.map(temperature_bands).map(landsat_LST).map(average_GSL_waterstation_temperature)\n",
    "\n",
    "combined_BRB_NDWI_col = ee.ImageCollection(innerJoin.apply(BRB_NDWI, BRB_landsat, toyFilter).map(joinedCol))\n",
    "combined_BRB_NDWI_col_unclipped = ee.ImageCollection(innerJoin.apply(unclipped_BRB_NDWI, BRB_landsat, toyFilter).map(joinedCol))\n",
    "\n",
    "### NDWI stats dataframes ###\n",
    "# these dataframes contain the water area and associated date and are to be exported #\n",
    "\n",
    "landsat_GSLNA_ndwi_df = pd.merge(pd.DataFrame({'NA_Water_Area': GSL_NA_NDWI.aggregate_array('ndwi').getInfo()}),\n",
    "                                                pd.DataFrame({'Date': landsat_GSL_NA_mosaic_col_dates}), right_index=True, left_index=True, how='outer') #new\n",
    "landsat_GSLSA_ndwi_df = pd.merge(pd.DataFrame({'SA_Water_Area': GSL_SA_NDWI.aggregate_array('ndwi').getInfo()}),\n",
    "                                                pd.DataFrame({'Date': landsat_GSL_mosaic_col_dates}), right_index=True, left_index=True, how='outer')\n",
    "\n",
    "landsat_BRB_ndwi_df = pd.merge(pd.DataFrame({'BRB_Water_Area': BRB_NDWI.aggregate_array('ndwi').getInfo()}),\n",
    "                                                pd.DataFrame({'Date': BRB_landsat_dates}), right_index=True, left_index=True, how='outer') #new\n",
    "\n",
    "\n",
    "landsat_GSLNA_LST_df = pd.merge(pd.DataFrame({'Temperature (C)': GSL_NA_LST.aggregate_array('LST').getInfo()}),\n",
    "                                                pd.DataFrame({'Date': landsat_GSL_NA_mosaic_col_dates}), right_index=True, left_index=True, how='outer')\n",
    "landsat_GSLSA_LST_df = pd.merge(pd.DataFrame({'Temperature (C)': GSL_SA_LST.aggregate_array('LST').getInfo()}),\n",
    "                                                pd.DataFrame({'Date': landsat_GSL_mosaic_col_dates}), right_index=True, left_index=True, how='outer')\n",
    "landsat_GSL_USGS_LST_df = pd.merge(pd.DataFrame({'Temperature (C)': GSL_SA_USGS_LST.aggregate_array('LST').getInfo()}),\n",
    "                                                pd.DataFrame({'Date': landsat_GSL_mosaic_col_dates}), right_index=True, left_index=True, how='outer')\n",
    "\n",
    "##Other useful dataframes!!\n",
    "landsat_GSL_NA_NDWI_threshold_df = pd.merge(pd.DataFrame({'Threshold': GSL_NA_NDWI.aggregate_array('threshold').getInfo()}),\n",
    "                                                pd.DataFrame({'Date': landsat_GSL_NA_mosaic_col_dates}), right_index=True, left_index=True, how='outer')\n",
    "landsat_GSL_NA_cloud_cover_df = pd.merge(pd.DataFrame({'Threshold': GSL_NA_NDWI.aggregate_array('CLOUD_COVER').getInfo()}),\n",
    "                                                pd.DataFrame({'Date': landsat_GSL_NA_mosaic_col_dates}), right_index=True, left_index=True, how='outer')\n",
    "landsat_GSL_SA_NDWI_threshold_df = pd.merge(pd.DataFrame({'Threshold': GSL_SA_NDWI.aggregate_array('threshold').getInfo()}),\n",
    "                                                pd.DataFrame({'Date': landsat_GSL_mosaic_col_dates}), right_index=True, left_index=True, how='outer')\n",
    "landsat_GSL_SA_cloud_cover_df = pd.merge(pd.DataFrame({'Threshold': GSL_SA_NDWI.aggregate_array('CLOUD_COVER').getInfo()}),\n",
    "                                                pd.DataFrame({'Date': landsat_GSL_mosaic_col_dates}), right_index=True, left_index=True, how='outer')\n",
    "\n",
    "\n",
    "landsat_BRB_NDWI_threshold_df = pd.merge(pd.DataFrame({'Threshold': BRB_NDWI.aggregate_array('threshold').getInfo()}),\n",
    "                                                pd.DataFrame({'Date': BRB_landsat_dates}), right_index=True, left_index=True, how='outer')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 29,
   "metadata": {},
   "outputs": [],
   "source": [
    "### NDWI collections ###\n",
    "### Sentinel\n",
    "sent_GSL_SA_NDWI = sentinel_GSL_mosaic_col.map(sentinel_ndwi_fn).map(extract_water_SA_S2).map(SouthArmClip).map(GSLSA_WaterPixelCount_Sent) \n",
    "sent_GSL_NA_NDWI = sentinel_GSL_mosaic_col.map(sentinel_ndwi_fn).map(extract_water_NA_S2).map(MaskSaltPool).map(NorthArmClip).map(GSLNA_WaterPixelCount_Sent) \n",
    "\n",
    "sent_GSL_NDWI_unclipped = sentinel_GSL_mosaic_col.combine(sentinel_GSL_mosaic_col.map(sentinel_ndwi_fn))\n",
    "\n",
    "### NDWI stats dataframes ###\n",
    "# these dataframes contain the water area and associated date and are to be exported #\n",
    "\n",
    "sentinel_GSLNA_ndwi_df = pd.merge(pd.DataFrame({'NA_Water_Area': sent_GSL_NA_NDWI.aggregate_array('ndwi').getInfo()}),\n",
    "                                                pd.DataFrame({'Date': sent_GSL_NA_NDWI.aggregate_array('Date_Filter').getInfo()}), right_index=True, left_index=True, how='outer') \n",
    "sentinel_GSLSA_ndwi_df = pd.merge(pd.DataFrame({'SA_Water_Area': sent_GSL_SA_NDWI.aggregate_array('ndwi').getInfo()}),\n",
    "                                                pd.DataFrame({'Date': sent_GSL_SA_NDWI.aggregate_array('Date_Filter').getInfo()}), right_index=True, left_index=True, how='outer') \n",
    "\n",
    "sentinel_GSL_NA_cloud_cover_df = pd.merge(pd.DataFrame({'%_Cloudy_Pixels': GSL_sentinel_N.aggregate_array('CLOUDY_PIXEL_PERCENTAGE').getInfo()}),\n",
    "                                                pd.DataFrame({'Date': sentinel_GSL_mosaic_col_dates}), right_index=True, left_index=True, how='outer')\n",
    "sentinel_GSL_SA_cloud_cover_df = pd.merge(pd.DataFrame({'%_Cloudy_Pixels': GSL_sentinel_S.aggregate_array('CLOUDY_PIXEL_PERCENTAGE').getInfo()}),\n",
    "                                                pd.DataFrame({'Date': GSL_sent_dates_S}), right_index=True, left_index=True, how='outer')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {},
   "outputs": [],
   "source": [
    "### NDVI collections ###\n",
    "### Landsat\n",
    "GSL_SA_NDVI = landsat_GSL_mosaic_col.map(landsat_ndvi_fn).map(SouthArmClip).map(GSLSA_VegPixelCount)\n",
    "GSL_SA_NDVI_S2Area = landsat_GSL_mosaic_col.map(landsat_ndvi_fn).map(MaskLandsatOverflow).map(SouthArmClip).map(GSLSA_VegPixelCount)\n",
    "GSL_NA_NDVI = landsat_GSL_NA_mosaic_col.map(MaskSaltPool).map(landsat_ndvi_fn).map(NorthArmClip).map(GSLNA_VegPixelCount)\n",
    "\n",
    "BRB_NDVI = BRB_landsat.map(landsat_ndvi_fn).map(BRBClip).map(BRB_VegPixelCount)\n",
    "\n",
    "### NDVI stats dataframes ###\n",
    "# these dataframes contain the vegetation area and associated date and are to be exported #\n",
    "#Landsat\n",
    "\n",
    "landsat_GSLNA_ndvi_df = pd.merge(pd.DataFrame({'NA_Veg_Area': GSL_NA_NDVI.aggregate_array('ndvi').getInfo()}),\n",
    "                                                pd.DataFrame({'Date': landsat_GSL_NA_mosaic_col_dates}), right_index=True, left_index=True, how='outer') #Take the counted number of pixels from ndwi image and multiply by 900 square meters, which is the area per pixel. Join with date column dataframe\n",
    "landsat_GSLSA_ndvi_df = pd.merge(pd.DataFrame({'SA_Veg_Area': GSL_SA_NDVI.aggregate_array('ndvi').getInfo()}),\n",
    "                                                pd.DataFrame({'Date': landsat_GSL_mosaic_col_dates}), right_index=True, left_index=True, how='outer') #Take the counted number of pixels from ndwi image and multiply by 900 square meters, which is the area per pixel. Then join with date column dataframe\n",
    "landsat_GSLSA_ndvi_S2Area_df = pd.merge(pd.DataFrame({'SA_Veg_Area': GSL_SA_NDVI_S2Area.aggregate_array('ndvi').getInfo()}),\n",
    "                                                pd.DataFrame({'Date': landsat_GSL_mosaic_col_dates}), right_index=True, left_index=True, how='outer')\n",
    "\n",
    "# Landsat\n",
    "landsat_BRB_ndvi_df = pd.merge(pd.DataFrame({'NA_Veg_Area': BRB_NDVI.aggregate_array('ndvi').getInfo()}),\n",
    "                                                pd.DataFrame({'Date': BRB_landsat_dates}), right_index=True, left_index=True, how='outer') #Take the counted number of pixels from ndwi image and multiply by 900 square meters, which is the area per pixel. Join with date column dataframe"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 30,
   "metadata": {},
   "outputs": [],
   "source": [
    "### NDVI collections ###\n",
    "### Sentinel\n",
    "sent_GSL_SA_NDVI = sentinel_GSL_mosaic_col.map(sentinel_ndvi_fn).map(MaskSaltPool).map(SouthArmClip).map(GSLSA_VegPixelCount_Sent)\n",
    "sent_GSL_NA_NDVI = sentinel_GSL_mosaic_col.map(sentinel_ndvi_fn).map(MaskSaltPool).map(NorthArmClip).map(GSLNA_VegPixelCount_Sent)\n",
    "\n",
    "### NDVI stats dataframes ###\n",
    "# these dataframes contain the vegetation area and associated date and are to be exported #\n",
    "\n",
    "sent_GSLNA_ndvi_df = pd.merge(pd.DataFrame({'NA_Veg_Area': sent_GSL_NA_NDVI.aggregate_array('ndvi').getInfo()}),\n",
    "                                                pd.DataFrame({'Date': sent_GSL_NA_NDVI.aggregate_array('Date_Filter').getInfo()}), right_index=True, left_index=True, how='outer') \n",
    "sent_GSLSA_ndvi_df = pd.merge(pd.DataFrame({'SA_Veg_Area': sent_GSL_SA_NDVI.aggregate_array('ndvi').getInfo()}),\n",
    "                                                pd.DataFrame({'Date': sent_GSL_SA_NDVI.aggregate_array('Date_Filter').getInfo()}), right_index=True, left_index=True, how='outer') "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {},
   "outputs": [],
   "source": [
    "### Halite index collections ###\n",
    "### Landsat\n",
    "#GSl shoreline\n",
    "GSL_NA_halite = combined_NA_NDWI_co_unclipped.map(GSL_NA_halite_watermask).map(MaskSaltPool).map(landsat_veg_mask).map(landsat_halite_fn_conservative).map(NorthArmClip).map(halite_pixel_count_landsat_GSLNA) \n",
    "GSL_SA_halite = combined_SA_NDWI_co_unclipped.map(GSL_SA_halite_watermask).map(landsat_veg_mask).map(FarmingtonBaySaltMask).map(landsat_halite_fn_conservative).map(SouthArmClip).map(halite_pixel_count_landsat_GSLSA)\n",
    "BRB_halite = combined_BRB_NDWI_col_unclipped.map(BRB_halite_watermask).map(landsat_veg_mask).map(MaskBRBSaltPool).map(landsat_halite_fn_conservative).map(BRBClip).map(halite_pixel_count_landsat_BRB)\n",
    "\n",
    "### Halite stats datatframe ###\n",
    "\n",
    "landsat_halite_GSL_NA_df = pd.merge(pd.DataFrame({'Halite_Area': GSL_NA_halite.aggregate_array('halite').getInfo()}),\n",
    "                                                pd.DataFrame({'Date': landsat_GSL_NA_mosaic_col_dates}), right_index=True, left_index=True, how='outer') \n",
    "landsat_halite_GSL_SA_df = pd.merge(pd.DataFrame({'Halite_Area': GSL_SA_halite.aggregate_array('halite').getInfo()}),\n",
    "                                                pd.DataFrame({'Date': landsat_GSL_mosaic_col_dates}), right_index=True, left_index=True, how='outer')   \n",
    "landsat_halite_BRB_df = pd.merge(pd.DataFrame({'Halite_Area': BRB_halite.aggregate_array('halite').getInfo()}),\n",
    "                                                pd.DataFrame({'Date': BRB_landsat_dates}), right_index=True, left_index=True, how='outer') "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 31,
   "metadata": {},
   "outputs": [],
   "source": [
    "### Halite index collections ###\n",
    "### Sentinel\n",
    "GSL_halite_index_sentinel = sent_GSL_NDWI_unclipped.map(GSL_S2_watermask).map(FarmingtonBaySaltMask).map(sentinel_veg_mask).map(MaskSaltPool).map(sentinel_halite_fn).map(TotalGSLClip)\n",
    "GSL_SA_halite_sent = GSL_halite_index_sentinel.map(halite_pixel_count_sent_SA)\n",
    "GSL_NA_halite_sent = GSL_halite_index_sentinel.map(halite_pixel_count_sent_NA)\n",
    "\n",
    "### Halite stats dataframes ###\n",
    "\n",
    "sentinel_SA_halite_df = pd.merge(pd.DataFrame({'Halite_Area': GSL_SA_halite_sent.aggregate_array('halite').getInfo()}),\n",
    "                                                pd.DataFrame({'Date': GSL_SA_halite_sent.aggregate_array('Date_Filter').getInfo()}), right_index=True, left_index=True, how='outer')\n",
    "sentinel_NA_halite_df = pd.merge(pd.DataFrame({'Halite_Area': GSL_NA_halite_sent.aggregate_array('halite').getInfo()}),\n",
    "                                                pd.DataFrame({'Date': GSL_NA_halite_sent.aggregate_array('Date_Filter').getInfo()}), right_index=True, left_index=True, how='outer')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {},
   "outputs": [],
   "source": [
    "### Gypsum index collections ###\n",
    "### Landsat\n",
    "#GSl shoreline\n",
    "GSL_NA_gypsum = combined_NA_NDWI_co_unclipped.map(GSL_NA_halite_watermask).map(MaskSaltPool).map(landsat_halite_mask).map(landsat_veg_mask).map(landsat_gypsum_fn).map(NorthArmClip).map(gypsum_pixel_count_landsat_GSLNA) \n",
    "GSL_SA_gypsum = combined_SA_NDWI_co_unclipped.map(GSL_SA_halite_watermask).map(FarmingtonBaySaltMask).map(landsat_veg_mask).map(landsat_halite_mask).map(landsat_gypsum_fn).map(SouthArmClip).map(gypsum_pixel_count_landsat_GSLSA)\n",
    "BRB_gypsum = combined_BRB_NDWI_col_unclipped.map(BRB_halite_watermask).map(landsat_veg_mask).map(landsat_halite_mask).map(landsat_gypsum_fn).map(BRBClip).map(gypsum_pixel_count_landsat_BRB)\n",
    "\n",
    "### Gypsum stats datatframe ###\n",
    "\n",
    "landsat_gypsum_GSL_NA_df = pd.merge(pd.DataFrame({'Gypsum_Area': GSL_NA_gypsum.aggregate_array('gypsum').getInfo()}),\n",
    "                                                pd.DataFrame({'Date': landsat_GSL_NA_mosaic_col_dates}), right_index=True, left_index=True, how='outer') \n",
    "landsat_gypsum_GSL_SA_df = pd.merge(pd.DataFrame({'Gypsum_Area': GSL_SA_gypsum.aggregate_array('gypsum').getInfo()}),\n",
    "                                                pd.DataFrame({'Date': landsat_GSL_mosaic_col_dates}), right_index=True, left_index=True, how='outer')  \n",
    "landsat_gypsum_BRB_df = pd.merge(pd.DataFrame({'Gypsum_Area': BRB_gypsum.aggregate_array('gypsum').getInfo()}),\n",
    "                                                pd.DataFrame({'Date': BRB_landsat_dates}), right_index=True, left_index=True, how='outer') "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 32,
   "metadata": {},
   "outputs": [],
   "source": [
    "### Gypsum index collections ###\n",
    "### Sentinel\n",
    "#GSl shoreline\n",
    "GSL_gypsum_sent = sent_GSL_NDWI_unclipped.map(GSL_S2_watermask).map(MaskSaltPool).map(FarmingtonBaySaltMask).map(sentinel_halite_mask).map(sentinel_veg_mask).map(sentinel_gypsum_fn).map(MaskSaltPool).map(TotalGSLClip)\n",
    "GSL_NA_gypsum_sent = GSL_gypsum_sent.map(gypsum_pixel_count_sent_NA)\n",
    "GSL_SA_gypsum_sent = GSL_gypsum_sent.map(gypsum_pixel_count_sent_SA)\n",
    "\n",
    "### Gypsum stats datatframe ###\n",
    "\n",
    "gypsum_GSL_NA_S2_df = pd.merge(pd.DataFrame({'Gypsum_Area': GSL_NA_gypsum_sent.aggregate_array('gypsum').getInfo()}),\n",
    "                                                pd.DataFrame({'Date': GSL_NA_gypsum_sent.aggregate_array('Date_Filter').getInfo()}), right_index=True, left_index=True, how='outer') \n",
    "gypsum_GSL_SA_S2_df = pd.merge(pd.DataFrame({'Gypsum_Area': GSL_SA_gypsum_sent.aggregate_array('gypsum').getInfo()}),\n",
    "                                                pd.DataFrame({'Date': GSL_SA_gypsum_sent.aggregate_array('Date_Filter').getInfo()}), right_index=True, left_index=True, how='outer') "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {},
   "outputs": [],
   "source": [
    "### Other sediments isolation ###\n",
    "#Landsat\n",
    "GSL_NA_carbonates = combined_NA_NDWI_co_unclipped.map(GSL_NA_halite_watermask).map(MaskSaltPool).map(landsat_veg_mask).map(MaskToCarbonatesL8).map(NorthArmClip).map(carbonate_pixel_count_landsat_GSLNA) #all seds except halite, gypsum, water, and vegetation, .map(MaskSaltPool)\n",
    "GSL_SA_carbonates = combined_SA_NDWI_co_unclipped.map(GSL_SA_halite_watermask).map(MaskSaltPool).map(landsat_veg_mask).map(MaskToCarbonatesL8).map(SouthArmClip).map(carbonate_pixel_count_landsat_GSLSA)\n",
    "BRB_carbonates = combined_BRB_NDWI_col_unclipped.map(BRB_halite_watermask).map(landsat_veg_mask).map(MaskToCarbonatesL8).map(BRBClip).map(carbonate_pixel_count_landsat_BRB)\n",
    "\n",
    "### Carbs stats datatframe ###\n",
    "\n",
    "landsat_carbonates_GSL_NA_df = pd.merge(pd.DataFrame({'Carb_muds_Area': GSL_NA_carbonates.aggregate_array('carb_muds').getInfo()}),\n",
    "                                                pd.DataFrame({'Date': landsat_GSL_NA_mosaic_col_dates}), right_index=True, left_index=True, how='outer') \n",
    "landsat_carbonates_GSL_SA_df = pd.merge(pd.DataFrame({'Carb_muds_Area': GSL_SA_carbonates.aggregate_array('carb_muds').getInfo()}),\n",
    "                                                pd.DataFrame({'Date': landsat_GSL_mosaic_col_dates}), right_index=True, left_index=True, how='outer')  \n",
    "landsat_carbonates_BRB_df = pd.merge(pd.DataFrame({'Carb_muds_Area': BRB_carbonates.aggregate_array('carb_muds').getInfo()}),\n",
    "                                                pd.DataFrame({'Date': BRB_landsat_dates}), right_index=True, left_index=True, how='outer') "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 33,
   "metadata": {},
   "outputs": [],
   "source": [
    "### Other sediments isolation ###\n",
    "#Sentinel\n",
    "GSL_carbonates_sn = sent_GSL_NDWI_unclipped.map(GSL_S2_watermask).map(sentinel_veg_mask).map(MaskSaltPool).map(MaskToCarbonatesS2).map(TotalGSLClip)\n",
    "GSL_NA_carbonates_sn = GSL_carbonates_sn.map(NorthArmClip).map(carbonate_pixel_count_sent_NA)\n",
    "GSL_SA_carbonates_sn = GSL_carbonates_sn.map(SouthArmClip).map(carbonate_pixel_count_sent_SA)\n",
    "\n",
    "### Carb muds stats datatframe ###\n",
    "\n",
    "carbonates_GSL_NA_S2_df = pd.merge(pd.DataFrame({'Carb_muds_Area': GSL_NA_carbonates_sn.aggregate_array('carb_muds').getInfo()}),\n",
    "                                                pd.DataFrame({'Date': GSL_NA_carbonates_sn.aggregate_array('Date_Filter').getInfo()}), right_index=True, left_index=True, how='outer')\n",
    "carbonates_GSL_SA_S2_df = pd.merge(pd.DataFrame({'Carb_muds_Area': GSL_SA_carbonates_sn.aggregate_array('carb_muds').getInfo()}),\n",
    "                                                pd.DataFrame({'Date': GSL_SA_carbonates_sn.aggregate_array('Date_Filter').getInfo()}), right_index=True, left_index=True, how='outer') "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [],
   "source": [
    "### Chlorophyll / cyanobacteria detection ###\n",
    "### Landsat Fluoresence Line Height adapted from Beck, R., Zhan, S., Liu, H., Tong, S., Yang, B., Xu, M., Ye, Z., Huang, Y., \n",
    "# Shu, S., Wu, Q., Wang, S., Berling, K., Murray, A., Emery, E., Reif, M., Harwood, J., Young, J., Nietch, C., Macke, D., … Su, H. (2016). \n",
    "# Comparison of satellite reflectance algorithms for estimating chlorophyll-a in a temperate reservoir using coincident hyperspectral \n",
    "#aircraft imagery and dense coincident surface observations. Remote Sensing of Environment, 178, 15–30. https://doi.org/10.1016/j.rse.2016.03.002\n",
    "\n",
    "# GSL_NA_FLHV = combined_NA_NDWI_col.map(landsat_MaskToWater_NA).map(MaskSaltPool).map(landsat_KIVU).map(NorthArmClip).map(FLHV_pixel_count_landsat_GSLNA) #for landsat 5\n",
    "# GSL_SA_FLHV = combined_SA_NDWI_col.map(landsat_MaskToWater_SA).map(landsat_KIVU).map(SouthArmClip).map(FLHV_pixel_count_landsat_GSLSA)\n",
    "\n",
    "# GSL_NA_FLHV = combined_NA_NDWI_col.map(landsat_MaskToWater_NA).map(MaskSaltPool).map(landsat_FLHV).map(NorthArmClip).map(FLHV_pixel_count_landsat_GSLNA) #for landsat 8, need to figure out way to normalize results\n",
    "# GSL_SA_FLHV = combined_SA_NDWI_col.map(landsat_MaskToWater_SA).map(landsat_FLHV).map(SouthArmClip).map(FLHV_pixel_count_landsat_GSLSA)\n",
    "\n",
    "GSL_NA_avg_FLHV = combined_NA_NDWI_co_unclipped.map(landsat_MaskToWater_NA).map(MaskSaltPool).map(landsat_KIVU).map(NorthArmClip).map(FLHV_avg_landsat_GSLNA) #for getting average chla value for all images using the KIVU method\n",
    "GSL_SA_avg_FLHV = combined_SA_NDWI_co_unclipped.map(landsat_MaskToWater_SA).map(landsat_KIVU).map(SouthArmClip).map(FLHV_avg_landsat_GSLSA)\n",
    "\n",
    "### Algae stats datatframe ###\n",
    "landsat_avg_algae_GSL_NA_df = pd.merge(pd.DataFrame({'Algae_Area': GSL_NA_avg_FLHV.aggregate_array('FLHV').getInfo()}),\n",
    "                                                pd.DataFrame({'Date': landsat_GSL_NA_mosaic_col_dates}), right_index=True, left_index=True, how='outer') \n",
    "landsat_avg_algae_GSL_SA_df = pd.merge(pd.DataFrame({'Algae_Area': GSL_SA_avg_FLHV.aggregate_array('FLHV').getInfo()}),\n",
    "                                                pd.DataFrame({'Date': landsat_GSL_mosaic_col_dates}), right_index=True, left_index=True, how='outer') \n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [],
   "source": [
    "### Chlorophyll / cyanobacteria detection ###\n",
    "GSL_2BDA = sent_GSL_NDWI_unclipped.map(sentinel_MaskToWater).map(MaskSaltPool).map(sentinel_2BDA).map(TotalGSLClip)\n",
    "GSL_2BDA_NA = GSL_2BDA.map(BDA_pixel_count_sent_NA).map(NorthArmClip)\n",
    "GSL_2BDA_SA = GSL_2BDA.map(BDA_pixel_count_sent_SA).map(SouthArmClip)\n",
    "\n",
    "GSL_2BDA_avg_NA = GSL_2BDA.map(BDA_average_sent_NA).map(NorthArmClip)\n",
    "GSL_2BDA_avg_SA = GSL_2BDA.map(BDA_average_sent_SA).map(SouthArmClip)\n",
    "\n",
    "### 2BDA cyanobacteria stats datatframe ###\n",
    "\n",
    "algae_avg_GSL_NA_S2_df = pd.merge(pd.DataFrame({'Avg_Chl_a_2BDA': GSL_2BDA_avg_NA.aggregate_array('2BDA').getInfo()}),\n",
    "                                                pd.DataFrame({'Date': GSL_2BDA_avg_NA.aggregate_array('Date_Filter').getInfo()}), right_index=True, left_index=True, how='outer') \n",
    "algae_avg_GSL_SA_S2_df = pd.merge(pd.DataFrame({'Avg_Chl_a_2BDA': GSL_2BDA_avg_SA.aggregate_array('2BDA').getInfo()}),\n",
    "                                                pd.DataFrame({'Date': GSL_2BDA_avg_NA.aggregate_array('Date_Filter').getInfo()}), right_index=True, left_index=True, how='outer')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 34,
   "metadata": {},
   "outputs": [],
   "source": [
    "### Save file center ###\n",
    "#Uncomment needed lines and insert a path to the folder and filename desired for file export where '' is shown in \"to_csv('INSERT/PATH/HERE')\"\n",
    "\n",
    "# landsat_GSL_NA_NDWI_threshold_df.to_csv('')\n",
    "# landsat_GSL_NA_cloud_cover_df.to_csv('')\n",
    "# landsat_GSL_SA_NDWI_threshold_df.to_csv('')\n",
    "# landsat_GSL_SA_cloud_cover_df.to_csv('')\n",
    "# sentinel_GSL_NA_cloud_cover_df.to_csv('')\n",
    "# sentinel_GSL_SA_cloud_cover_df.to_csv('')\n",
    "\n",
    "# landsat_halite_GSL_NA_df.to_csv('')\n",
    "# landsat_halite_GSL_SA_df.to_csv('')\n",
    "# landsat_gypsum_GSL_NA_df.to_csv('')\n",
    "# landsat_gypsum_GSL_SA_df.to_csv('')\n",
    "# landsat_carbonates_GSL_NA_df.to_csv('')\n",
    "# landsat_carbonates_GSL_SA_df.to_csv('')\n",
    "# landsat_GSLNA_ndwi_df.to_csv('')\n",
    "# landsat_GSLSA_ndwi_df.to_csv('')\n",
    "# landsat_algae_GSL_NA_df.to_csv('')\n",
    "# landsat_algae_GSL_SA_df.to_csv('')\n",
    "# landsat_avg_algae_GSL_NA_df.to_csv('')\n",
    "# landsat_avg_algae_GSL_SA_df.to_csv('')\n",
    "# landsat_GSLSA_ndvi_df.to_csv('')\n",
    "# landsat_GSLSA_ndvi_S2Area_df.to_csv('')\n",
    "# landsat_GSLNA_ndvi_df.to_csv('')\n",
    "# landsat_GSLNA_LST_df.to_csv('')\n",
    "# landsat_GSLSA_LST_df.to_csv('')\n",
    "# landsat_GSL_USGS_LST_df.to_csv('')\n",
    "\n",
    "# landsat_BRB_NDWI_threshold_df.to_csv('')\n",
    "# landsat_BRB_ndwi_df.to_csv('')\n",
    "# landsat_halite_BRB_df.to_csv('')\n",
    "# landsat_gypsum_BRB_df.to_csv('')\n",
    "# landsat_carbonates_BRB_df.to_csv('')\n",
    "# landsat_BRB_ndvi_df.to_csv('')\n",
    "\n",
    "# sentinel_SA_halite_df.to_csv('')\n",
    "# sentinel_NA_halite_df.to_csv('')\n",
    "# gypsum_GSL_NA_S2_df.to_csv('')\n",
    "# gypsum_GSL_SA_S2_df.to_csv('')\n",
    "# carbonates_GSL_NA_S2_df.to_csv('')\n",
    "# carbonates_GSL_SA_S2_df.to_csv('')\n",
    "# sentinel_GSLNA_ndwi_df.to_csv('')\n",
    "# sentinel_GSLSA_ndwi_df.to_csv('')\n",
    "# algae_GSL_NA_S2_df.to_csv('')\n",
    "# algae_GSL_SA_S2_df.to_csv('')\n",
    "# algae_avg_GSL_NA_S2_df.to_csv('')\n",
    "# algae_avg_GSL_SA_S2_df.to_csv('')\n",
    "# sent_GSLNA_ndvi_df.to_csv('')\n",
    "# sent_GSLSA_ndvi_df.to_csv('')\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "### Export file center ###\n",
    "NA_export_col = landsat_GSL_NA_mosaic_col.select(['SR_B4', 'SR_B3', 'SR_B2']).map(NA_Export_Clip)\n",
    "SA_export_col = landsat_GSL_mosaic_col.select(['SR_B4', 'SR_B3', 'SR_B2']).map(SA_Export_Clip)\n",
    "NA_S2_export_col = sentinel_GSL_mosaic_col.select(['B4', 'B3', 'B2']).map(NA_Export_Clip)\n",
    "SA_S2_export_col = sentinel_GSL_mosaic_col.select(['B4', 'B3', 'B2']).map(SA_Export_Clip)\n",
    "\n",
    "### Export collections one at a time -- separated cells below for indivdual export calls ###\n",
    "\n",
    "# geemap.ee_export_image_collection_to_drive(NA_export_col, descriptions=landsat_GSL_NA_mosaic_col_dates, folder='GSLNA_TrueColor_Images', scale=60, region=NA_export_region, skipEmptyTiles=True)\n",
    "# geemap.ee_export_image_collection_to_drive(SA_export_col, descriptions=landsat_GSL_mosaic_col_dates, folder='GSLSA_TrueColor_Images', scale=60, region=SA_export_region, skipEmptyTiles=True)\n",
    "\n",
    "# geemap.ee_export_image_collection_to_drive(NA_S2_export_col, descriptions=sentinel_GSL_mosaic_col_dates, folder='GSLNA_TrueColor_Images', scale=60, region=NA_export_region, skipEmptyTiles=True)\n",
    "geemap.ee_export_image_collection_to_drive(SA_S2_export_col, descriptions=sentinel_GSL_mosaic_col_dates, folder='GSLSA_TrueColor_Images', scale=60, region=SA_export_region, skipEmptyTiles=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "geemap.ee_export_image_collection_to_drive(GSL_SA_NDWI, descriptions=landsat_GSL_mosaic_col_dates, folder='GSLSA_NDWI_Timeseries_1984_1994', scale=30, region=GSL_export_region, skipEmptyTiles=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "geemap.ee_export_image_collection_to_drive(GSL_SA_NDVI, descriptions=landsat_GSL_mosaic_col_dates, folder='GSLSA_NDVI_Timeseries_1984_1994', scale=30, region=GSL_export_region, skipEmptyTiles=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "geemap.ee_export_image_collection_to_drive(GSL_NA_NDVI, descriptions=landsat_GSL_NA_mosaic_col_dates, folder='GSLNA_NDVI_Timeseries_1984_1994', scale=30, region=GSL_export_region, skipEmptyTiles=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "geemap.ee_export_image_collection_to_drive(GSL_SA_halite, descriptions=landsat_GSL_mosaic_col_dates, folder='GSLSA_halite_Timeseries_1984_1994', scale=30, region=GSL_export_region, skipEmptyTiles=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "geemap.ee_export_image_collection_to_drive(GSL_NA_halite, descriptions=landsat_GSL_NA_mosaic_col_dates, folder='GSLNA_halite_Timeseries_1984_1994_saltpool', scale=30, region=GSL_export_region, skipEmptyTiles=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "geemap.ee_export_image_collection_to_drive(GSL_SA_gypsum, descriptions=landsat_GSL_mosaic_col_dates, folder='GSLSA_gypsum_Timeseries_1984_1994', scale=30, region=GSL_export_region, skipEmptyTiles=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "geemap.ee_export_image_collection_to_drive(GSL_NA_gypsum, descriptions=landsat_GSL_NA_mosaic_col_dates, folder='GSLNA_gypsum_Timeseries_1984_1994_saltpool', scale=30, region=GSL_export_region, skipEmptyTiles=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "geemap.ee_export_image_collection_to_drive(GSL_SA_carbonates, descriptions=landsat_GSL_mosaic_col_dates, folder='GSLSA_carbonates_Timeseries_1984_1994', scale=30, region=GSL_export_region, skipEmptyTiles=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "geemap.ee_export_image_collection_to_drive(GSL_NA_carbonates, descriptions=landsat_GSL_NA_mosaic_col_dates, folder='GSLNA_carbonates_Timeseries_1984_1994_saltpool', scale=30, region=GSL_export_region, skipEmptyTiles=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "geemap.ee_export_image_collection_to_drive(GSL_NA_LST, descriptions=landsat_GSL_NA_mosaic_col_dates, folder='GSLNA_LST_14_23', scale=30, region=GSL_export_region, skipEmptyTiles=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "geemap.ee_export_image_collection_to_drive(GSL_SA_LST, descriptions=landsat_GSL_mosaic_col_dates, folder='GSLSA_LST_14_23', scale=30, region=GSL_export_region, skipEmptyTiles=True)"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "gee-old",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.9.13"
  },
  "orig_nbformat": 4,
  "vscode": {
   "interpreter": {
    "hash": "61c4f21bc5437dca716a52008504e78a1129b65c35f1f049b4da8503bf984a77"
   }
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}