{ "cells": [ { "cell_type": "markdown", "id": "63faa58e-08a3-435d-bb8f-79b8849045ca", "metadata": {}, "source": [ "This notebook contains the processing algorithms used for CONTRASTS spectral albedo data. Processing is done in three steps:\n", "1. Conversion from v0 (raw data) to v1: reads in albedo notes, ASD files, cruise event list, and ship track data and converts the data to long format.\n", "2. Conversion from v1 to v2 (automated QC): applies seven automated quality control flags to the data. Saves as H5 files for easy interfacing with pandas and CSV files suitable for export.\n", "3. Conversion from v2 to v3 (final data): updates quality control flags with any manual updates, saves final dataset\n", "\n", "Author: Nora Greenwood, adapted from David-Clemens Sewall and Maddie Smith" ] }, { "cell_type": "code", "execution_count": null, "id": "621afb59-1026-49e9-a7b0-e38eefd1aedd", "metadata": {}, "outputs": [], "source": [ "import os\n", "import glob\n", "import pandas as pd\n", "import numpy as np\n", "from openpyxl import load_workbook\n", "from itertools import compress\n", "import csv\n", "import matplotlib.pyplot as plt\n", "import seaborn as sns\n", "import datetime as dt" ] }, { "cell_type": "markdown", "id": "48b71a65-0dcd-4ba3-baf5-77f5e0996a1e", "metadata": {}, "source": [ "# v0 to v1" ] }, { "cell_type": "code", "execution_count": null, "id": "9f2bc225-0e06-4774-8806-82b7a692f2a2", "metadata": {}, "outputs": [], "source": [ "# paths\n", "v0_path = r\"C:\\Users\\Nora Greenwood\\Desktop\\CONTRASTS\\processing_albedo_data\\v0_CONTRASTS_albedo\"\n", "\n", "# make a list of things in v0_path that 1) are folders and 2) have an asd file in them and 3) have a notes file in them and 4) are a transect\n", "asd_dirs = [subdir for subdir in os.listdir(v0_path)\n", " if os.path.isdir(os.path.join(v0_path, subdir)) and\n", " any(fname.startswith(\"asd\") for fname in os.listdir(os.path.join(v0_path, subdir)))\n", " and any(fname.startswith(\"alb\") for fname in os.listdir(os.path.join(v0_path, subdir)))]\n", "\n", "v1_path = r\"C:\\Users\\Nora Greenwood\\Desktop\\CONTRASTS\\processing_albedo_data\\v1_CONTRASTS_albedo\\spectral\"\n", "notes_path = r\"C:\\Users\\Nora Greenwood\\Desktop\\CONTRASTS\\processing_albedo_data\\v1_CONTRASTS_albedo\\albedonotes\"" ] }, { "cell_type": "code", "execution_count": null, "id": "7743b845-f333-4606-bbcd-ee0bceef59a5", "metadata": {}, "outputs": [], "source": [ "# load in GPS data\n", "gps_path = r\"C:\\Users\\Nora Greenwood\\Desktop\\CONTRASTS\\processing_albedo_data\\PS149_mastertrack.tab\"\n", "# Read in ship track\n", "# data format is list of lists of strings of dateTtime, lat, lon, speed [kn], course [deg]\n", "gps_data = []\n", "with open(gps_path, newline='') as f:\n", " reader = csv.reader(f, delimiter='\\t')\n", " for row in reader:\n", " # Skip empty lines\n", " if not row:\n", " continue\n", " # Check if this line starts the real data\n", " if row[0].startswith('2025'):\n", " # Convert to numbers if appropriate\n", " gps_data.append(row)\n", "\n", "gps_dt = []\n", "lats = []\n", "lons = []\n", "for row in gps_data:\n", " date = row[0]\n", " gps_dt.append(dt.datetime(int(date[0:4]),int(date[5:7]),int(date[8:10]),int(date[11:13]),int(date[14:16])))\n", " lats.append(float(row[1]))\n", " lons.append(float(row[2]))\n", "\n", "gps = {\"datetime\": gps_dt, \"latitude\": lats, \"longitude\": lons}\n", "\n", "# Convert to a dataframe\n", "gps_df = pd.DataFrame(gps)" ] }, { "cell_type": "code", "execution_count": null, "id": "9e781c2a-ec6f-47e7-a0fc-98362a613b1f", "metadata": {}, "outputs": [], "source": [ "# load in event IDs\n", "event_path = r\"C:\\Users\\Nora Greenwood\\Desktop\\CONTRASTS\\processing_albedo_data\\eventIDs.xlsx\"\n", "events = pd.read_excel(event_path,'ASD')" ] }, { "cell_type": "code", "execution_count": null, "id": "45387cc1-e351-46b4-9534-bffd97111aac", "metadata": {}, "outputs": [], "source": [ "def process_asd(site, out_path):\n", " \"\"\"\n", " Reads in albedo notes, ASD files, ship track data, and event list and converts them into long format albedo\n", "\n", " parameters: \n", " -------\n", " site: string\n", " folder names in v0 directory and excerpt of file names. formatted as YYYYMMDD-rFV-sitename\n", " out_path: string\n", " path to output v1 files to\n", " \"\"\"\n", " files = os.listdir(v0_path + '\\\\' + site)\n", "\n", " asd_path = os.path.join(v0_path, site, 'asd-contrasts-'+site+'.xlsx')\n", " note_path = os.path.join(notes_path, 'albedonotes-contrasts-'+site+'-v1.csv')\n", " scan_path = os.path.join(v0_path, site, 'asd-contrasts-'+site) \n", " \n", " ## Parse ASD file\n", " ws = load_workbook(asd_path, read_only=True, data_only=True)['Template']\n", " \n", " # get positions\n", " for value in ws.iter_rows(min_row=2, max_row=2, min_col=48, values_only=True):\n", " positions = list(value)\n", " # duplicate positions because only incident positions are written\n", " for i in range(len(positions)//2):\n", " positions[2*i+1] = positions[2*i]\n", " # Get I/R row\n", " for value in ws.iter_rows(min_row=3, max_row=3, min_col=48, max_col = 300, values_only=True):\n", " ir_row = list(value)\n", " # Get scannumber row\n", " for value in ws.iter_rows(min_row=5, max_row=5, min_col=48, max_col = 300, values_only=True):\n", " sn_row = list(value)\n", " for value in ws.iter_rows(min_row=4, max_row=4, min_col=48, max_col = 300, values_only=True):\n", " rep_row = list(value)\n", " for i in range(len(rep_row)//2):\n", " rep_row[2*i+1] = rep_row[2*i]\n", " # Crop end of lists where scan number row is None\n", " asd_length = max([i for i in range(len(sn_row)) if sn_row[i] != None])+1\n", " notNone = [sn != None for sn in sn_row]\n", " positions = list(compress(positions, notNone))\n", " ir_row = list(compress(ir_row, notNone))\n", " sn_row = list(compress(sn_row, notNone))\n", " rep_row = list(compress(rep_row, notNone))\n", " ## Get data\n", " asd = np.zeros((2151, asd_length))\n", " for i, value in zip(np.arange(2151), \n", " ws.iter_rows(min_row=6, max_row=2156, min_col=48, max_col=48+asd_length-1, values_only=True)):\n", " asd[i,:] = value\n", " # filter out skipped positions\n", " asd = asd[:, ~np.all(np.isnan(asd), axis=0)]\n", " \n", " ## Wavelength index (fixed)\n", " wv_ind = range(350, 2501)\n", "\n", " # pull times from when files were created\n", " sn_times = []\n", " scan_path = os.path.join(v0_path, site, 'asd-'+site)\n", " sn_times = []\n", " for sn in sn_row:\n", " sn_path = os.path.join(scan_path, sn)\n", " sn_times.append(dt.datetime.fromtimestamp(os.path.getmtime(sn_path),dt.UTC))\n", " start_time = min(sn_times)\n", "\n", " # make header\n", " #read in met file again to get additional header information\n", " header = pd.read_csv(note_path,header=None, nrows=2)\n", " header = header.fillna('NaN')\n", " # start_time = header.iloc[1,3]\n", " date = dt.datetime(int(site[0:4]),int(site[4:6]),int(site[6:8])).date()\n", " datetime = dt.datetime.combine(date,start_time.time())\n", " idx = (gps_df['datetime'] - datetime).abs().idxmin()\n", " lat = round(gps_df['latitude'][idx],5)\n", " lon = round(gps_df['longitude'][idx],5)\n", "\n", " event_id = events.loc[events['Site'] == site, 'Event IDs'].squeeze()\n", " \n", " sky_cond = header.iloc[0,3]\n", " other = str(header.iloc[1,5])\n", " other = other.replace(',',';')\n", " line = header.iloc[0,1]\n", " \n", " \n", " ### Construct albedo dataframe\n", " df_alb = pd.DataFrame(data=asd.T, index=sn_row, columns=wv_ind)\n", " df_alb['position'] = positions\n", " df_alb['variable'] = ['incident_solar_counts' if ir=='I' else 'reflected_solar_counts' for ir in ir_row]\n", " df_alb['time_utc'] = sn_times\n", " df_alb['site'] = line\n", " df_alb['repetition'] = rep_row\n", "\n", " # converts to longform by breaking up scans by wavelength\n", " df_alb_lon = df_alb.melt(id_vars=['time_utc', 'site', 'position', 'repetition', 'variable'],\n", " var_name='wavelength')\n", " \n", " # Create albedo values in new df\n", " df_piv = df_alb_lon[['position', 'wavelength', 'variable', 'value','repetition']\n", " ].set_index(['position', 'wavelength','repetition']).pivot(columns='variable', values='value')\n", " s_alb = df_piv['reflected_solar_counts'] / df_piv['incident_solar_counts']\n", " s_alb.name = 'value'\n", " df_alb_temp = df_alb_lon[df_alb_lon['variable']=='reflected_solar_counts'][\n", " ['time_utc', 'site', 'position', 'repetition', 'wavelength']\n", " ].set_index(['position', 'wavelength','repetition']).join(s_alb)\n", " df_alb_temp['variable'] = 'albedo'\n", " df_alb_temp = df_alb_temp.reset_index()\n", " # Create albedo quality flag field\n", " df_alb_qual = df_alb_temp.copy()\n", " df_alb_qual['variable'] = 'asd_quality_flag'\n", " df_alb_qual['value'] = -999.9\n", " df_alb_qual = df_alb_qual.reset_index()\n", " df_alb_qual.drop(columns=['index'], inplace=True)\n", " \n", " ### Construct ancillary dataframes\n", " df_anc = pd.DataFrame(data={'position': positions, 'time_utc': sn_times, 'repetition': rep_row, 'lat': lat, 'lon': lon,\n", " 'variable': ['incident_solar_counts' if ir=='I' else 'reflected_solar_counts' for ir in ir_row]}, index=sn_row)\n", " \n", " ## Create time difference\n", " df_piv = df_anc[['position', 'variable', 'repetition','time_utc']].pivot(columns='variable', index=['position','repetition'], values='time_utc')\n", " sr_td_s = (df_piv['reflected_solar_counts'] - df_piv['incident_solar_counts']).dt.total_seconds()\n", " sr_td_s.name = 'value'\n", " # Join to reflected timestamps\n", " df_td = df_anc[df_anc['variable']=='reflected_solar_counts'][['position', 'time_utc','repetition']].join(sr_td_s, on=['position','repetition'])\n", " \n", " # Add rest of index\n", " df_td['variable'] = 'asd_time_diff_s'\n", " df_td['site'] = line\n", " df_td['wavelength'] = np.nan\n", " ## Create lat lon dataframe\n", " df_ll = df_anc[['time_utc', 'position', 'lat', 'lon','repetition']].melt(\n", " id_vars=['position','repetition','time_utc'], var_name='variable')\n", " df_ll['site'] = line\n", " df_ll['wavelength'] = np.nan\n", " \n", " # Combine data, time difference, and lat-lon into one dataframe\n", " df_out = pd.concat([df_alb_lon, df_alb_temp, df_alb_qual, df_td])\n", "\n", "\n", " #write a new csv with the header and data\n", " # write the header\n", " header = 'Site, ' + line + '\\n Event ID, ' + str(event_id) + '\\n Date, '+ str(date) + '\\n Start time (UTC), '+ str(start_time.time()) + ' \\n Sky conditions, '+ sky_cond+' \\n Other notes, ' +other+' \\n Ship latitude, '+str(lat)+'\\n Ship longitude, '+str(lon)+'\\n'\n", " header = header.replace(', nan\\n',', NaN\\n')\n", " \n", " newname = 'spectral-contrasts-' + site + '-v1.csv'\n", " newfilename = os.path.join(out_path, newname)\n", " \n", " with open(newfilename, 'w') as fp:\n", " fp.write(header)\n", " \n", " # write the rest\n", " df_out.to_csv(newfilename, header=True,index=True, mode='a')\n", "\n", " print('wrote: '+newname)" ] }, { "cell_type": "code", "execution_count": null, "id": "98e20bec-b239-4ec2-bb2a-e80f240b2cfb", "metadata": {}, "outputs": [], "source": [ "# convert data\n", "for i in range(len(asd_dirs)):\n", " try:\n", " process_asd(asd_dirs[i], v1_path)\n", " except Exception as e:\n", " print('')\n", " print(asd_dirs[i])\n", " print(f\"An unexpected error occurred: {e}\")\n", " print('')" ] }, { "cell_type": "markdown", "id": "84170f51-d535-4ac2-b05c-75df68d9e000", "metadata": {}, "source": [ "# v1 to v2" ] }, { "cell_type": "markdown", "id": "6326047a-c0b4-4d49-aa73-a887fed7a238", "metadata": {}, "source": [ "This code should takes v1 data (long format) and adds automated QC flags. Adapted from David Clemen-Sewall's SALVO albedo processing code.\n", "\n", "\n", "First, \n", " 1. Parabolic correction from 750 to 1000 nm for temperature-offset between VNIR and SWIR (procedure from ASD manual).\n", "\n", "Next, the following automated checks for which the entire scan is flagged:\n", "\n", " 2. Visible wavelengths include albedos greater than 1 (Flag 2)\n", " 3. Slope of albedos between 400 and 450 nm is less than 1.6e-4 (Flag 3)\n", " 4. Time between incident and reflected measurement is less than seven seconds (Flag 4)\n", " 5. Slope condition and transition time condition are both met (Flag 5)\n", " \n", "Then, the following automated checks flag portions of the scan that meet certain conditions.\n", "\n", " 6. Variability in the incident by wavelength. Greater than 15% change averaged between the preceeding and following measurement (Flag 6)\n", " 7. Insufficient incident irradiance; below 1250 counts for the SWIR1 detector (1001 - 1800 nm) and below 5500 counts for the SWIR2 detector (1801 - 2500 nm). (Flag 7)\n", "\n", "Finally, albedos where the reflected measurement is negative are set to zero." ] }, { "cell_type": "code", "execution_count": null, "id": "841fee32-801a-4746-a281-ccbdca637345", "metadata": {}, "outputs": [], "source": [ "# paths\n", "v0_path = r\"C:\\Users\\Nora Greenwood\\Desktop\\CONTRASTS\\processing_albedo_data\\v0_CONTRASTS_albedo\"\n", "v1_path = r\"C:\\Users\\Nora Greenwood\\Desktop\\CONTRASTS\\processing_albedo_data\\v1_CONTRASTS_albedo\\spectral\"\n", "v1_notes = r\"C:\\Users\\Nora Greenwood\\Desktop\\CONTRASTS\\processing_albedo_data\\v1_CONTRASTS_albedo\\albedonotes\"\n", "v2_path = r\"C:\\Users\\Nora Greenwood\\Desktop\\CONTRASTS\\processing_albedo_data\\v2_CONTRASTS_albedo\\spectral\"\n", "\n", "v1_dirs = os.listdir(v1_path)" ] }, { "cell_type": "code", "execution_count": null, "id": "5540be88-4eb2-4a9d-ac62-3554d3ad3304", "metadata": {}, "outputs": [], "source": [ "def parabolic_correction(df):\n", " \"\"\"\n", " Apply parabolic correction to rescale albedos between 750 and 1000 nm to smooth\n", " offset from 1000 to 1001 nm\n", " \n", " Parameters:\n", " -----------\n", " df: pandas DataFrame\n", " Contains the long-format spectral albedo data from a single day-site\n", " \n", " Returns:\n", " --------\n", " None. But df is modified in-place.\n", " \"\"\"\n", " \n", " # Get albedo offset from 1000 to 1001 nm\n", " wvs = [1000, 1001]\n", " df_para = df.query('variable == \"albedo\" & wavelength in @wvs').reset_index()\n", " df_para = df.pivot(columns='wavelength', values='value',\n", " index=['time_utc', 'site', 'repetition', 'position', 'variable'])\n", " df_para['factor'] = (df_para[1001] - df_para[1000])/df_para[1000]\n", " df_para.drop(columns=wvs, inplace=True)\n", "\n", " # Create dataframe of multiplicative scales\n", " wvs = np.arange(251).T + 750\n", " arr_scales = ((wvs - 750)**2)/(251**2) # parabolic transformed wavelength coordinate\n", " arr_scales = np.outer(df_para['factor'].values, arr_scales) + 1\n", " df_scales = pd.DataFrame(data=arr_scales, index=df_para.index, columns=wvs)\n", " df_scales = df_scales.melt(var_name='wavelength', ignore_index=False)\n", " df_scales.set_index('wavelength', append=True, inplace=True)\n", " df_scales = df_scales.reorder_levels(\n", " ['time_utc', 'site', 'position', 'repetition', 'wavelength', 'variable']) \n", "\n", " # Update values in df\n", " df.update((df_scales*df).query('variable == \"albedo\"'))" ] }, { "cell_type": "code", "execution_count": null, "id": "9363dadd-f310-49cb-b631-50c0efac3354", "metadata": {}, "outputs": [], "source": [ "def vis_alb_filter(df, wv_max=750, thrsh=1.0):\n", " \"\"\"\n", " Flag entire scans when any visible albedo exceeds thrsh\n", " \n", " Parameters:\n", " -----------\n", " df: pandas DataFrame\n", " Contains the long-format spectral albedo data from a single day-site\n", " wv_max: int, optional\n", " Maximum wavelength considered part of the visible. The default is 750 (nm).\n", " thrh: float, optional\n", " Albedo threshold at and above which to flag scan. The default is 1.0.\n", " \n", " Returns:\n", " --------\n", " None. But df is modified in-place.\n", " \"\"\"\n", " \n", " # Get positions to filter\n", " df_temp = df.query('wavelength <= @wv_max & variable == \"albedo\"')\n", " pos_filt = np.unique(df_temp[df_temp['value']>=1.0].index.get_level_values('position').values)\n", " \n", " # Update flags\n", " df.loc[df.eval('variable == \"asd_quality_flag\" & position in @pos_filt'), 'value'] = 2.0" ] }, { "cell_type": "code", "execution_count": null, "id": "b909add2-2709-4d82-8734-292806f68952", "metadata": {}, "outputs": [], "source": [ "def slope_filter(df, slope_threshold=-.002, wv_min=400, wv_max=450):\n", " \"\"\"\n", " Flags scans in-place slope between wv_min and wv_max is less than slope_threshold\n", " \n", " Parameters:\n", " -----------\n", " df: pandas DataFrame\n", " Contains the long-format spectral albedo data from a single day-site\n", " slope_threshold: float, optional\n", " below this slope, flag entire scan. The default is -1.6e-4\n", " wv_min: int, optional\n", " Lower end of wavelength range. The default is 400 (nm).\n", " wv_max: int, optional\n", " Upper end of wavelength range. The default is 450 (nm).\n", " \n", " Returns:\n", " --------\n", " None. But df is modified in-place.\n", " \"\"\"\n", " \n", " # Extract albedos\n", " df_alb = df.query('variable == \"albedo\" & wavelength >= @wv_min & wavelength <= @wv_max'\n", " ).reset_index().drop(\n", " columns=['time_utc', 'site', 'variable'])\n", " df_alb = df_alb.pivot(columns=['position','repetition'], index='wavelength', values='value')\n", "\n", " # Compute slopes of linear regressions\n", " arr_pos = df_alb.columns.values\n", " arr_slopes = np.zeros(arr_pos.shape)\n", " for i in np.arange(arr_pos.size):\n", " arr_slopes[i] = np.polynomial.Polynomial.fit(\n", " df_alb.index.values, df_alb[arr_pos[i]].values, 1).coef[1]\n", "\n", " # Get positions to filter\n", " pos_filt = arr_pos[arr_slopes < slope_threshold]\n", " \n", " bad_pairs = pd.MultiIndex.from_tuples(pos_filt, names=['position', 'repetition'])\n", " \n", " # Build mask using index.get_level_values\n", " mask = ((df.index.get_level_values('variable') == 'asd_quality_flag') &\n", " (pd.MultiIndex.from_arrays([df.index.get_level_values('position'), df.index.get_level_values('repetition')]).isin(bad_pairs)))\n", " \n", " df.loc[mask, 'value'] = 3" ] }, { "cell_type": "code", "execution_count": null, "id": "96fe9fc2-446f-4887-9454-0dbe9fb8c7f4", "metadata": {}, "outputs": [], "source": [ "def time_diff_filter(df, min_flip_s=7):\n", " \"\"\"\n", " Flags scans in-place if time between incident and reflected measurement is less than min_flip_s\n", " \n", " Parameters:\n", " -----------\n", " df: pandas DataFrame\n", " Contains the long-format spectral albedo data from a single day-site\n", " min_flip_s: float, optional\n", " below this flip time, flag entire scan. The default is 7 (seconds)\n", " \n", " Returns:\n", " --------\n", " None. But df is modified in-place.\n", " \"\"\"\n", " \n", " # Get positions to filter\n", " pos_filt = df.query('variable == \"asd_time_diff_s\" & value <= @min_flip_s').index.get_level_values('position')\n", " \n", " # Update flags\n", " df_a4.loc[df.eval('variable == \"asd_quality_flag\" & position in @pos_filt' +\n", " ' & value != 3.0'), 'value'] = 4.0\n", " df_a4.loc[df.eval('variable == \"asd_quality_flag\" & position in @pos_filt' +\n", " ' & value == 3.0'), 'value'] = 5.0" ] }, { "cell_type": "code", "execution_count": null, "id": "cf9c1f0d-0feb-4c29-8b2f-02177819e827", "metadata": {}, "outputs": [], "source": [ "def inc_var_filter(df, chg_thrsh=0.15, time_diff = 180):\n", " \"\"\"\n", " Flags wavelengths where the incident changes by more than chg_thrsh\n", " \n", " Parameters:\n", " -----------\n", " df: pandas DataFrame\n", " Contains the long-format spectral albedo data from a single day-site\n", " chg_thrsh: float, optional\n", " Above this ratio of change in incident (averaged between scan before and after), flag wavelength. The default is 0.15\n", " time_diff: float, optional\n", " Only consecutive scans within this time window will be taken into consideration for change in incident. The default is 180 seconds.\n", " \n", " Returns:\n", " --------\n", " None. But df is modified in-place.\n", " \"\"\"\n", " \n", " # Get dataframe of incident values\n", " df_inc = df.query('variable == \"incident_solar_counts\"').reset_index()\n", " df_inc = df_inc.pivot(columns='wavelength', values='value',\n", " index=['time_utc', 'site', 'repetition', 'position', 'variable'])\n", " df_inc.sort_index(level='time_utc', inplace=True)\n", "\n", "\n", " # retrieve time stamps for incident scans\n", " time_index = pd.to_datetime(df_inc.index.get_level_values('time_utc'))\n", " # calculate time difference to previous and next scan\n", " dt_prev = time_index.to_series().diff().dt.total_seconds()\n", " dt_next = time_index.to_series().diff(-1).abs().dt.total_seconds()\n", "\n", " # calculate change in incident (as function of wavelength) for prev and next scan\n", " chg_prev = df_inc.diff(periods=1).abs() / df_inc\n", " chg_next = df_inc.diff(periods=-1).abs() / df_inc\n", "\n", " mask = dt_prev.values[:, None] >= time_diff\n", " mask = pd.Series(np.squeeze(mask), index=df_inc.index) \n", " chg_prev[mask] = np.nan\n", " \n", " mask = dt_next.values[:, None] >= time_diff\n", " mask = pd.Series(np.squeeze(mask), index=df_inc.index) \n", " chg_next[mask] = np.nan\n", "\n", " # check how many neighbors to average over (1 or 2)\n", " n_valid = (~chg_prev.isna()).astype(int) + (~chg_next.isna()).astype(int)\n", "\n", " # average prev and next, counting nan as 0 and dividing by non nans\n", " df_chg = (chg_prev.fillna(0) + chg_next.fillna(0)) / n_valid\n", " # remove those with no neighbors considered\n", " df_chg = df_chg.where(n_valid > 0)\n", "\n", " # Get positions and wavelengths exceeding threshold\n", " df_mask = (df_chg > chg_thrsh).melt(ignore_index=False)\n", " df_mask.reset_index(['time_utc', 'variable'], inplace=True)\n", " df_mask['variable'] = 'asd_quality_flag'\n", " df_mask.drop(columns='time_utc', inplace=True)\n", " df_mask.set_index(['wavelength', 'variable'], append=True, inplace=True)\n", " df_mask = df_mask.reorder_levels(\n", " ['site', 'position', 'repetition', 'wavelength', 'variable'])\n", "\n", " \n", "\n", " # Update flags (ugly... indices don't match because incident and reflected have diff timestamps)\n", " df.reset_index('time_utc', inplace=True)\n", " df.loc[df_mask.index[df_mask['value']], 'value'] = 6.0\n", " df.set_index('time_utc', inplace=True, append=True)\n", " df = df.reorder_levels(['time_utc', 'site', 'position', 'repetition', 'wavelength', 'variable'])" ] }, { "cell_type": "code", "execution_count": null, "id": "e5cd3f61-b59a-4abb-acdb-ef6a55e582dd", "metadata": {}, "outputs": [], "source": [ "def varying_light_filter(df, num_thrsh = 425, max_wv = 1201):\n", " \"\"\"\n", " Flags entire scan where more than num_thrsh wavelengths below max_wv recieved flag 6 (varying light) \n", " \n", " Parameters:\n", " -----------\n", " df: pandas DataFrame\n", " Contains the long-format spectral albedo data from a single day-site\n", " num_thrsh: float, optional\n", " Above this number of flagged incident light wavelengths, entire scan is flagged\n", " max_wv: float, optional\n", " Maximum wavelength considered for counting flag 6s\n", " \n", " Returns:\n", " --------\n", " None. But df is modified in-place.\n", " \"\"\"\n", " ## check how many flag 6s there are below max_wv, flag entire scan if more than num_thrsh\n", "\n", " df_qf = df.query('variable == \"asd_quality_flag\" & wavelength <= @max_wv').reset_index().drop(columns=['time_utc', 'site', 'variable'])\n", " df_qf = df_qf.pivot(columns=['position','repetition'], index='wavelength', values='value')\n", "\n", " # Compute slopes of linear regressions\n", " arr_pos = df_qf.columns.values\n", " counts = np.zeros(arr_pos.shape)\n", " for i in np.arange(arr_pos.size):\n", " counts[i] = sum(df_qf[arr_pos[i]].values == 6.0)\n", "\n", " # Get positions to filter\n", " pos_filt = arr_pos[counts > num_thrsh]\n", " \n", " bad_pairs = pd.MultiIndex.from_tuples(pos_filt, names=['position', 'repetition'])\n", " \n", " # Build mask using index.get_level_values\n", " mask = ((df.index.get_level_values('variable') == 'asd_quality_flag') & \n", " (pd.MultiIndex.from_arrays([df.index.get_level_values('position'), df.index.get_level_values('repetition')]).isin(bad_pairs)))\n", " \n", " df.loc[mask, 'value'] = 6" ] }, { "cell_type": "code", "execution_count": null, "id": "d3789348-2c84-43f0-a7de-6be2b1870fc5", "metadata": {}, "outputs": [], "source": [ "def inc_mag_filter(df, swir1_thrsh=400, swir2_thrsh=900):\n", " \"\"\"\n", " Flags wavelengths where the incident is below thresholds\n", " \n", " Parameters:\n", " -----------\n", " df: pandas DataFrame\n", " Contains the long-format spectral albedo data from a single day-site\n", " swir1_thrsh: float, optional\n", " Below this incident value flag swir1 wavelengths. The default is 400 (counts)\n", " swir2_thrsh: float, optional\n", " Below this incident value flag swir2 wavelengths. The default is 900 (counts)\n", " \n", " Returns:\n", " --------\n", " None. But df is modified in-place.\n", " \"\"\"\n", " \n", " df.reset_index('time_utc', inplace=True)\n", " \n", " # swir1\n", " wv_min = 1001\n", " wv_max = 1800\n", " s_mask = df.query('variable == \"incident_solar_counts\" ' +\n", " '& wavelength >= @wv_min & wavelength <= @wv_max')['value'] <= swir1_thrsh\n", " df_mask = s_mask.reset_index('variable')\n", " df_mask['variable'] = 'asd_quality_flag'\n", " df_mask.set_index('variable', append=True, inplace=True)\n", " df.loc[df_mask.index[df_mask['value']], 'value'] = 7.0\n", "\n", " # swir2\n", " wv_min = 1801\n", " s_mask = df.query('variable == \"incident_solar_counts\" ' +\n", " '& wavelength >= @wv_min')['value'] <= swir2_thrsh\n", " df_mask = s_mask.reset_index('variable')\n", " df_mask['variable'] = 'asd_quality_flag'\n", " df_mask.set_index('variable', append=True, inplace=True)\n", " df.loc[df_mask.index[df_mask['value']], 'value'] = 7.0\n", "\n", " df.set_index('time_utc', inplace=True, append=True)\n", " df = df.reorder_levels(\n", " ['time_utc', 'site', 'position', 'repetition', 'wavelength', 'variable'])\n", " " ] }, { "cell_type": "code", "execution_count": null, "id": "cc35241d-28ce-4257-b6d9-95234e076c12", "metadata": {}, "outputs": [], "source": [ "def negative_reflected(df):\n", " \"\"\"\n", " Sets albedo to zero where reflected value is negative (negative incident values are caught with flag 7)\n", " \n", " Parameters:\n", " -----------\n", " df: pandas DataFrame\n", " Contains the long-format spectral albedo data from a single day-site\n", " \n", " Returns:\n", " --------\n", " df\n", " \"\"\"\n", " \n", " df = df.reorder_levels(['time_utc', 'site', 'position', 'repetition', 'wavelength', 'variable'])\n", "\n", " # Get negative reflected values\n", " df_neg = df.query('variable == \"reflected_solar_counts\" & value < 0.0').copy()\n", "\n", " # Update albedo values\n", " df_neg.reset_index('variable', inplace=True)\n", " df_neg['variable'] = 'albedo'\n", " df_neg.set_index('variable', append=True, inplace=True)\n", " df_neg = df_neg.reorder_levels(['time_utc', 'site', 'position', 'repetition', 'wavelength', 'variable'])\n", " df.loc[df_neg.index, 'value'] = 0.0\n", " \n", " return df" ] }, { "cell_type": "code", "execution_count": null, "id": "7eb2d064-5bf2-4f8a-9bba-8548980902a0", "metadata": {}, "outputs": [], "source": [ "def merge_with_notes(fname, df):\n", " \"\"\"\n", " Adds note fields (surface type, depth, notes) to albedo data arrays and flags entire scan if trampling was present (flag 0)\n", " \n", " Parameters:\n", " -----------\n", " fname: string\n", " Contains the v1 file name\n", " df: pandas DataFrame\n", " Contains the long-format spectral albedo data from a single day-site\n", " \n", " Returns:\n", " --------\n", " df_merged\n", " \"\"\"\n", " # now pull relevant stuff from albedo notes\n", " site = fname[19:-7]\n", " notes_path = os.path.join(v1_notes, 'albedonotes-contrasts-'+site+'-v1.csv')\n", " \n", " # read in notes and pull in position data\n", " notes = pd.read_csv(notes_path,header=2)\n", " # filter to just important notes columns\n", " notes_info = notes.loc[:,['Position','Snow or Pond (s/p)','S/P Depth (cm)','Notes on footprint','Trampled','Repetition']]\n", " notes_info = notes_info.rename(columns={'Position': 'position', 'Snow or Pond (s/p)': 'surface type', 'S/P Depth (cm)': 'surface thickness (cm)', 'Notes on footprint': 'notes', 'Trampled': 'trampled', 'Repetition':'repetition'})\n", " notes_info = notes_info.drop_duplicates()\n", " \n", " notes_info = notes_info.set_index(['position','repetition'])\n", " df_merged = df_a4.join(notes_info)\n", " \n", " df_merged['surface thickness (cm)'] = df_merged['surface thickness (cm)'].round(1)\n", "\n", " \n", " trampled = df_a4_merged.query('variable == \"asd_quality_flag\" & value == -999.9 & trampled == 1').copy()\n", " # Update flag values\n", " trampled.reset_index('variable', inplace=True)\n", " trampled['variable'] = 'asd_quality_flag'\n", " trampled.set_index('variable', append=True, inplace=True)\n", " trampled = trampled.reorder_levels(\n", " ['time_utc', 'site', 'position', 'repetition', 'wavelength', 'variable'])\n", " \n", " df_merged = df_merged.reorder_levels(\n", " ['time_utc', 'site', 'position', 'repetition', 'wavelength', 'variable'])\n", " df_merged.loc[trampled.index, 'value'] = 1.0\n", " df_merged = df_merged.drop('trampled',axis=1)\n", " \n", " return df_merged" ] }, { "cell_type": "code", "execution_count": null, "id": "95287406-11eb-486b-9b79-6d34a23e3c98", "metadata": {}, "outputs": [], "source": [ "def v1asd_to_v2(asd_v1, out_path):\n", " \"\"\"\n", " Applies automated quality flags to long-formatted spectral data (v0)\n", " \n", " Parameters:\n", " -----------\n", " asd_v1: string\n", " Contains the v1 file name\n", " out_path: string\n", " path to output v1 files to\n", " \n", " Returns:\n", " --------\n", " Nothing, but creates an h5 file with full dataset (incident and reflected counts) and three CSV files: (1) all data, (2) all corresponding quality flags and (3) just flag 0 data (passed all QC)\n", " \"\"\"\n", " # Read in a1 data\n", " asd_path = os.path.join(v1_path, asd_v1)\n", " df_v1 = pd.read_csv(asd_path, header=8, low_memory=False).set_index(\n", " ['time_utc', 'site', 'position', 'repetition', 'wavelength', 'variable'])\n", " \n", " # Correct VNIR/SWIR1 offset\n", " parabolic_correction(df_v1)\n", " # Apply whole scan filters\n", " vis_alb_filter(df_v1)\n", " slope_filter(df_v1)\n", " time_diff_filter(df_v1)\n", " \n", " # Apply wavelength filters\n", " inc_var_filter(df_v1)\n", " varying_light_filter(df_v1)\n", " inc_mag_filter(df_v1)\n", "\n", " # Set albedos to zero where reflected is negative\n", " df_flagged = negative_reflected(df_v1)\n", "\n", " # merges with notes and adds trampling flag\n", " df_merged = merge_with_notes(asd_v1, df_flagged)\n", " # Replace all flag -999.9 (unprocessed data) with flag 0.0 (data that has passed automated qc)\n", " df_merged.loc[df_merged.eval('variable == \"asd_quality_flag\" & value == -999.9'), 'value'] = 0.0\n", "\n", " # Write df_a4 to disk (useful for manual changes)\n", " v2_name = os.path.basename(asd_v1).replace('-v1.csv', '-fullv2.h5')\n", " df_merged.to_hdf(os.path.join(out_path, v2_name), key='data', mode='w', format='table',\n", " complevel=9, complib='blosc:lz4')\n", " print(\"wrote: \" + v2_name)\n", "\n", " # figure out total position range\n", " if 'Transect' in asd_v1:\n", " # add blank positions as needed\n", " if 'R2_main' in asd_v1:\n", " pos = pd.DataFrame({'position': range(0, 91, 5)})\n", " elif 'R1_mainB' in asd_v1:\n", " pos = pd.DataFrame({'position': range(0, 186, 5)})\n", " elif 'R1_mainA' in asd_v1:\n", " pos = pd.DataFrame({'position': range(0, 201, 5)})\n", " elif 'main' in asd_v1 or 'R1_opticsC' in asd_v1 or 'R2_tTransect' in asd_v1 or 'R3_optics' in asd_v1:\n", " pos = pd.DataFrame({'position': range(0, 101, 5)})\n", " elif 'raftedIce' in asd_v1:\n", " pos = pd.DataFrame({'position': range(15, 101, 5)})\n", " elif 'R2_optics' in asd_v1:\n", " pos = pd.DataFrame({'position': range(0, 91, 5)})\n", " elif 'R1_opticsB' in asd_v1:\n", " pos = pd.DataFrame({'position': range(5, 66, 5)})\n", " \n", " ## make 3 CSV outputs\n", " header = pd.read_csv(asd_path,header=None,on_bad_lines='skip')\n", " start_time = header.iloc[3,1]\n", " date = header.iloc[2,1]\n", " date = dt.datetime(int(date[1:5]),int(date[6:8]),int(date[9:11])).date()\n", " lat = header.iloc[6,1]\n", " lon = header.iloc[7,1]\n", " \n", " sky_cond = header.iloc[4,1]\n", " other = str(header.iloc[5,1])\n", " line = header.iloc[0,1]\n", " event_id = header.iloc[1,1]\n", "\n", " header = 'Site, ' + line + '\\n Event ID, ' + str(event_id) + '\\n Date, '+ str(date) + '\\n Start time (UTC), '+ str(start_time) + ' \\n Sky conditions, '+ sky_cond+' \\n Other notes, ' +other+' \\n Ship latitude, '+str(lat)+'\\n Ship longitude, '+str(lon)+'\\n'\n", " header = header.replace(', nan\\n',', NaN\\n')\n", " \n", " df_temp = df_merged.query('variable in [\"albedo\", \"asd_quality_flag\"] ').reset_index()\n", " df_temp = df_temp.pivot(columns='variable', values='value',\n", " index=['time_utc', 'position', 'repetition', 'wavelength', 'surface type', 'surface thickness (cm)', 'notes'])\n", " \n", " # JUST FLAG 0\n", " df_flag0 = df_temp[df_temp['asd_quality_flag'] == 0.0]['albedo'].reset_index()\n", " \n", " # make sure all blank columns get included for transects\n", " if 'Transect' in asd_v1:\n", " df_flag0 = pos.merge(df_flag0, on='position',how='left')\n", " \n", " df_flag0 = df_flag0.round(3) \n", " df_flag0 = df_flag0.fillna('NaN')\n", " df_flag0 = df_flag0.pivot(columns=['time_utc', 'repetition', 'position', 'surface type', 'surface thickness (cm)', 'notes'],\n", " index='wavelength', values='albedo')\n", " df_flag0 = df_flag0.reindex(range(350, 2501))\n", "\n", " v2_name = os.path.basename(asd_v1).replace('-v1.csv', '-v2-flag0.csv') \n", " newfilename = os.path.join(out_path, v2_name)\n", " with open(newfilename, 'w') as fp:\n", " fp.write(header)\n", " # write the rest\n", " df_flag0.to_csv(newfilename, header=True,index=True, mode='a', na_rep='NaN')\n", "\n", " ## ALL ALBEDO AND ALL FLAGS\n", " df_albedo = df_temp['albedo'].reset_index()\n", " if 'Transect' in asd_v1:\n", " df_albedo = pos.merge(df_albedo, on='position',how='left')\n", " # df_albedo = df_albedo.fillna('NaN')\n", " df_albedo = df_albedo.round(3)\n", " df_albedo = df_albedo.fillna('NaN')\n", " df_albedo = df_albedo.pivot(columns=['time_utc', 'repetition', 'position', 'surface type', 'surface thickness (cm)', 'notes'],\n", " index='wavelength', values='albedo')\n", " df_albedo = df_albedo.reindex(range(350, 2501))\n", "\n", " \n", " df_flags = df_temp['asd_quality_flag'].reset_index()\n", " if 'Transect' in asd_v1:\n", " df_flags = pos.merge(df_flags, on='position', how='left')\n", " df_flags = df_flags.fillna('NaN')\n", " df_flags = df_flags.pivot(columns=['time_utc', 'repetition', 'position', 'surface type', 'surface thickness (cm)', 'notes'],\n", " index='wavelength', values='asd_quality_flag')\n", " df_flags = df_flags.reindex(range(350, 2501))\n", " \n", " alb_name = os.path.basename(asd_v1).replace('-v1.csv', '-v2-albedo.csv')\n", " flag_name = os.path.basename(asd_v1).replace('-v1.csv', '-v2-flags.csv')\n", " \n", " newfilename = os.path.join(out_path, alb_name)\n", " \n", " with open(newfilename, 'w') as fp:\n", " fp.write(header)\n", " \n", " # write the rest\n", " df_albedo.to_csv(newfilename, header=True,index=True, mode='a', na_rep='NaN')\n", "\n", " newfilename = os.path.join(out_path, flag_name)\n", " \n", " with open(newfilename, 'w') as fp:\n", " fp.write(header)\n", " \n", " # write the rest\n", " df_flags.to_csv(newfilename, header=True,index=True, mode='a', na_rep='NaN')\n", " \n", " print(\"wrote: \" + v2_name)\n" ] }, { "cell_type": "code", "execution_count": null, "id": "308cb300-9f42-4a68-b900-24d4f73bc486", "metadata": {}, "outputs": [], "source": [ "for i in range(len(v1_dirs)):\n", " try:\n", " v1asd_to_v2(v1_dirs[i], v2_path)\n", " except Exception as e: \n", " print('')\n", " print(v1_dirs[i])\n", " print(f\"An unexpected error occurred: {e}\")\n", " print('')" ] }, { "cell_type": "markdown", "id": "33d515d9-6103-448d-b489-5d85885c1946", "metadata": {}, "source": [ "# v2 to v3" ] }, { "cell_type": "markdown", "id": "839b22bc-3185-4e61-9ab1-3908e7ed3aeb", "metadata": {}, "source": [ "This code applies manually updates to the automated quality flags created in v2, and also confirms there are no empty rows." ] }, { "cell_type": "code", "execution_count": null, "id": "00b783cd-bc95-4c21-bd17-ed45773d46c7", "metadata": {}, "outputs": [], "source": [ "# paths\n", "\n", "v2_path = r\"C:\\Users\\Nora Greenwood\\Desktop\\CONTRASTS\\processing_albedo_data\\v2_CONTRASTS_albedo\\spectral\"\n", "manual_path = r\"C:\\Users\\Nora Greenwood\\Desktop\\CONTRASTS\\processing_albedo_data\\v2_CONTRASTS_albedo\\manual flag updates\"\n", "v3_path = r\"C:\\Users\\Nora Greenwood\\Desktop\\CONTRASTS\\processing_albedo_data\\v3_CONTRASTS_albedo\\spectral\"\n", "\n", "v2_files = [fname for fname in os.listdir(v2_path) if fname.endswith('.csv')]\n", "v2_sites = sorted({fname.rsplit('-', 1)[0] for fname in v2_files})\n", "\n", "manual_files = [subdir for subdir in os.listdir(manual_path)\n", " if not os.path.isdir(os.path.join(manual_path, subdir))]\n", "manual_fix_sites = sorted({fname.rsplit('-', 1)[0] for fname in manual_files})" ] }, { "cell_type": "code", "execution_count": null, "id": "e928cfae-5ef1-49bf-a5f4-1662461866e2", "metadata": {}, "outputs": [], "source": [ "def file_making_from_v2(site, df, out_path):\n", " \"\"\"\n", " Outputs complete spectral dataset with manually updated quality flags to 3 CSVs\n", " \n", " Parameters:\n", " -----------\n", " site: string\n", " Contains the v2 file name\n", " df: pandas dataframe\n", " Contains the long-format spectral albedo data from a single day-site\n", " out_path: string\n", " path to output v1 files to\n", " \n", " Returns:\n", " --------\n", " Nothing, but creates three CSV files: (1) all data, (2) all corresponding quality flags and (3) just flag 0 data (passed all QC)\n", " \"\"\"\n", " # making paths and such\n", " asd_path = os.path.join(v2_path, site + '-flags.csv')\n", " \n", " # figure out total position range\n", " if 'Transect' in site:\n", " # add blank positions as needed\n", " if 'R2_main' in site:\n", " pos = pd.DataFrame({'position': range(0, 91, 5)})\n", " elif 'R1_mainB' in site:\n", " pos = pd.DataFrame({'position': range(0, 186, 5)})\n", " elif 'R1_mainA' in site:\n", " pos = pd.DataFrame({'position': range(0, 201, 5)})\n", " elif 'main' in site or 'R1_opticsC' in site or 'R2_tTransect' in site or 'R3_optics' in site:\n", " pos = pd.DataFrame({'position': range(0, 101, 5)})\n", " elif 'raftedIce' in site:\n", " pos = pd.DataFrame({'position': range(15, 101, 5)})\n", " elif 'R2_optics' in site:\n", " pos = pd.DataFrame({'position': range(0, 91, 5)})\n", " elif 'R1_opticsB' in site:\n", " pos = pd.DataFrame({'position': range(5, 66, 5)})\n", "\n", " ## make 3 CSV outputs\n", " \n", " header = pd.read_csv(asd_path,header=None,on_bad_lines='skip')\n", " start_time = header.iloc[3,1]\n", " date = header.iloc[2,1]\n", " date = dt.datetime(int(date[1:5]),int(date[6:8]),int(date[9:11])).date()\n", " lat = header.iloc[6,1]\n", " lon = header.iloc[7,1]\n", " \n", " sky_cond = header.iloc[4,1]\n", " other = str(header.iloc[5,1])\n", " line = header.iloc[0,1]\n", " event_id = header.iloc[1,1]\n", "\n", " header = 'Site, ' + line + '\\n Event ID, ' + str(event_id) + '\\n Date, '+ str(date) + '\\n Start time (UTC), '+ str(start_time) + ' \\n Sky conditions, '+ sky_cond+' \\n Other notes, ' +other+' \\n Ship latitude, '+str(lat)+'\\n Ship longitude, '+str(lon)+'\\n'\n", " header = header.replace(', nan\\n',', NaN\\n')\n", " \n", " df_temp = df.query('variable in [\"albedo\", \"asd_quality_flag\"] ').reset_index()\n", " df_temp = df_temp.pivot(columns='variable', values='value',\n", " index=['time_utc', 'position', 'repetition', 'wavelength', 'surface type', 'surface thickness (cm)', 'notes'])\n", " \n", " # JUST FLAG 0\n", " df_flag0 = df_temp[df_temp['asd_quality_flag'] == 0.0]['albedo'].reset_index()\n", " # make sure all blank columns get included for transects\n", " if 'Transect' in site:\n", " df_flag0 = pos.merge(df_flag0, on='position',how='left')\n", " \n", " df_flag0 = df_flag0.round(3)\n", " df_flag0 = df_flag0.fillna('NaN')\n", " df_flag0 = df_flag0.pivot(columns=['time_utc', 'repetition', 'position', 'surface type', 'surface thickness (cm)', 'notes'],\n", " index='wavelength', values='albedo')\n", " df_flag0 = df_flag0.reindex(range(350, 2501))\n", " df_flag0 = df_flag0.round(3)\n", " \n", " v3_name = os.path.basename(site).replace('-v2', '-v3-flag0.csv') \n", " newfilename = os.path.join(out_path, v3_name)\n", " \n", " with open(newfilename, 'w') as fp:\n", " fp.write(header)\n", " \n", " # write the rest\n", " df_flag0.to_csv(newfilename, header=True,index=True, mode='a', na_rep='NaN')\n", "\n", " ## ALL ALBEDO AND ALL FLAGS\n", " df_albedo = df_temp['albedo'].reset_index()\n", " if 'Transect' in site:\n", " df_albedo = pos.merge(df_albedo, on='position',how='left')\n", " df_albedo = df_albedo.round(3)\n", " df_albedo = df_albedo.fillna('NaN')\n", " df_albedo = df_albedo.pivot(columns=['time_utc', 'repetition', 'position', 'surface type', 'surface thickness (cm)', 'notes'],\n", " index='wavelength', values='albedo')\n", " df_albedo = df_albedo.reindex(range(350, 2501))\n", "\n", " \n", " df_flags = df_temp['asd_quality_flag'].reset_index()\n", " if 'Transect' in site:\n", " df_flags = pos.merge(df_flags, on='position', how='left')\n", " df_flags = df_flags.fillna('NaN')\n", " df_flags = df_flags.pivot(columns=['time_utc', 'repetition', 'position', 'surface type', 'surface thickness (cm)', 'notes'],\n", " index='wavelength', values='asd_quality_flag')\n", " df_flags = df_flags.reindex(range(350, 2501))\n", " df_flags = df_flags.round(3)\n", "\n", " alb_name = os.path.basename(site).replace('-v2', '-v3-albedo.csv')\n", " flag_name = os.path.basename(site).replace('-v2', '-v3-flags.csv')\n", " \n", " newfilename = os.path.join(out_path, alb_name)\n", " \n", " with open(newfilename, 'w') as fp:\n", " fp.write(header)\n", " \n", " # write the rest\n", " df_albedo.to_csv(newfilename, header=True,index=True, mode='a', na_rep='NaN')\n", "\n", " newfilename = os.path.join(out_path, flag_name)\n", " \n", " with open(newfilename, 'w') as fp:\n", " fp.write(header)\n", " \n", " # write the rest\n", " df_flags.to_csv(newfilename, header=True,index=True, mode='a', na_rep='NaN')\n", " \n", " print(\"wrote new versions of \" + site[:-3])" ] }, { "cell_type": "code", "execution_count": null, "id": "93f69863-bfd1-41dd-be80-d302d69cce9e", "metadata": {}, "outputs": [], "source": [ "def v2asd_to_v3(site, out_path):\n", " \"\"\"\n", " Applies manual updates to quality flags from long-formatted spectral data (v0) as needed and deletes any blank rows in created CSVs\n", " \n", " Parameters:\n", " -----------\n", " site: string\n", " Contains the v2 file name\n", " out_path: string\n", " path to output v3 files to\n", " \n", " Returns:\n", " --------\n", " Nothing, but replicates v2 CSV files in place if no manual changes and creates new CSV files with manually updated quality flags for other sites\n", " Returns three CSV files with all data, all corresponding quality flags and just flag 0 data (passed all QC)\n", " \"\"\"\n", " if not site in manual_fix_sites:\n", " for file in glob.glob(os.path.join(v2_path, f\"*{site}*.csv\")):\n", " # Build destination filename (v2 -> v3)\n", " v3_name = file.replace(\"-v2-\", \"-v3-\")\n", " v3_name = v3_name.replace(\"v2_\", \"v3_\")\n", "\n", " with open(file, newline='', encoding='utf-8') as f:\n", " rows = list(csv.reader(f))\n", " rows = [row for row in rows if not all(cell == 'NaN' for cell in row)]\n", " \n", " with open(v3_name, \"w\", newline='', encoding='utf-8') as f:\n", " csv.writer(f).writerows(rows)\n", " print('copied data from', site)\n", " print('')\n", "\n", " else:\n", " # do manual flag updates\n", " # load h5, update v3s, split into files and save\n", " h5_path = os.path.join(v2_path, site[:-3] + '-fullv2.h5')\n", " df = pd.read_hdf(h5_path)\n", "\n", " ## read in manual corrections\n", " corrections_path = os.path.join(manual_path, site + '-flags.csv')\n", " with open(corrections_path, newline='', encoding='utf-8') as f:\n", " reader = csv.reader(f)\n", " for i, row in enumerate(reader, start=1): \n", " ## be warned for some ridiculous reason, rows are 1-indexed and columns are 0 indexed\n", " if i == 10:\n", " repetitions = row[1:]\n", " elif i == 11:\n", " positions = row[1:] \n", " elif i == 16:\n", " if row[0] == \"350.0\" or row[0] == '350':\n", " flags = row[1:]\n", " break\n", " elif i == 17:\n", " flags = row[1:]\n", " break\n", " repetitions = [int(float(x)) if (x != 'nan' and x!= 'NaN') else np.nan for x in repetitions]\n", " if 'Transect' in site:\n", " positions = [int(x) if (x != 'nan' and x!= 'NaN') else np.nan for x in positions]\n", " flags = [int(float(x)) if (x != 'nan' and x!= 'NaN') else np.nan for x in flags]\n", "\n", " ## convert to dataframe\n", " corrections = pd.DataFrame({'position': positions, 'repetition': repetitions, 'value': flags})\n", "\n", " ## temporary dataframe to do corrections with\n", " df = df.reset_index().merge(corrections.rename(columns={'value': 'corr_value'}), on=['position', 'repetition'], how='left')\n", " mask = ((df['variable'] == 'asd_quality_flag') &\n", " (df['value'] != df['corr_value']) & (df['value'] < 6) & (df['value'] != 4))\n", " \n", " df.loc[mask, 'value'] = df.loc[mask, 'corr_value']\n", "\n", " df = df.drop(columns='corr_value')\n", "\n", " file_making_from_v2(site, df, out_path)\n", "\n", " ## same nan row deletion from other case\n", " for file in glob.glob(os.path.join(out_path, f\"*{site[:-3]}*.csv\")):\n", " with open(file, newline='', encoding='utf-8') as f:\n", " rows = list(csv.reader(f))\n", " \n", " rows = [row for row in rows if not all(cell == 'NaN' for cell in row)]\n", "\n", " with open(file, \"w\", newline='', encoding='utf-8') as f:\n", " csv.writer(f).writerows(rows)" ] }, { "cell_type": "code", "execution_count": null, "id": "80443271-a0de-4f6b-aaaf-1a2373ae75b6", "metadata": {}, "outputs": [], "source": [ "# convert data\n", "for i in range(len(v2_sites)):\n", " v2asd_to_v3(v2_sites[i], v3_path)" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3 (ipykernel)", "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.14.3" } }, "nbformat": 4, "nbformat_minor": 5 }