diff --git a/docs/source/handbook/operations/in_flight.md b/docs/source/handbook/operations/in_flight.md index 186c965..6ca5767 100644 --- a/docs/source/handbook/operations/in_flight.md +++ b/docs/source/handbook/operations/in_flight.md @@ -37,7 +37,7 @@ Once the [pre-flight operations](./pre_flight.md) are complete and the tests are 8. If AVAPS is connected to aircraft data (which it should be during flight operations), the meteorological data from the aircraft is taken in automatically. - 9. The dropsonde is now ready. Prepare the sonde for launching. First, disconnect the data cable and put the sonde in the re-conditioning-box. Click on `Continue`. + 9. The dropsonde is now ready. Prepare the sonde for launching. First, disconnect the data cable and put the sonde in the re-conditioning-box. Click on `Continue`. An error message might appear if the data cable is removed before clicking `Continue`. It may be solved by leaving the cable in, clicking `Continue`, then removing the cable. | ![Sonde ready](../../graphics/sonde_ready.png) | | :--------------------------------------------: | diff --git a/docs/source/handbook/operations/post_flight.md b/docs/source/handbook/operations/post_flight.md new file mode 100644 index 0000000..ba4c291 --- /dev/null +++ b/docs/source/handbook/operations/post_flight.md @@ -0,0 +1,6 @@ +# Post-flight Operations + +Data from the dropsondes launched during the flight can be transferred out of the AVAPS computer in two ways: + + 1. By plugging in a USB stick into the USB port on the front side of the HRDS rack. + 2. By connecting to the AVAPS computer remotely via FTP. FileZilla can be used to drag and drop data from AVAPS in the aircraft to a local PC. The IP address and password must be known to access AVAPS remotely. diff --git a/docs/source/howto/create_quicklooks.ipynb b/docs/source/howto/create_quicklooks.ipynb new file mode 100644 index 0000000..98ac46a --- /dev/null +++ b/docs/source/howto/create_quicklooks.ipynb @@ -0,0 +1,267 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "id": "1ff5dc6d-caec-4303-8c22-684792c6cfe5", + "metadata": {}, + "source": [ + "# Some dropsonde quicklook examples" + ] + }, + { + "cell_type": "code", + "execution_count": 1, + "id": "947b5fe5-e348-4d8c-b6fd-de449fc2aa44", + "metadata": {}, + "outputs": [ + { + "ename": "ModuleNotFoundError", + "evalue": "No module named 'halodrops'", + "output_type": "error", + "traceback": [ + "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m", + "\u001b[0;31mModuleNotFoundError\u001b[0m Traceback (most recent call last)", + "Cell \u001b[0;32mIn[1], line 1\u001b[0m\n\u001b[0;32m----> 1\u001b[0m \u001b[38;5;28;01mfrom\u001b[39;00m \u001b[38;5;21;01mhalodrops\u001b[39;00m\u001b[38;5;21;01m.\u001b[39;00m\u001b[38;5;21;01mplotting\u001b[39;00m \u001b[38;5;28;01mimport\u001b[39;00m quicklooks \u001b[38;5;28;01mas\u001b[39;00m ql\n", + "\u001b[0;31mModuleNotFoundError\u001b[0m: No module named 'halodrops'" + ] + } + ], + "source": [ + "from halodrops.plotting import quicklooks as ql" + ] + }, + { + "cell_type": "markdown", + "id": "6b2d5ef7-5111-4978-99dc-36b6a28fed2c", + "metadata": {}, + "source": [ + "### Test plotting for EUREC4A data" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "a3d0cd28-5bce-46d9-abe8-5c45cc32c14e", + "metadata": {}, + "outputs": [], + "source": [ + "import eurec4a\n", + "from datetime import datetime, date\n", + "from collections import defaultdict\n", + "from functools import reduce\n", + "\n", + "# Access EUREC4A catalog\n", + "cat = eurec4a.get_intake_catalog(use_ipfs=False)\n", + "\n", + "# Get flight segments\n", + "meta = eurec4a.get_flight_segments()\n", + "\n", + "segments = [{**s,\n", + " \"platform_id\": platform_id,\n", + " \"flight_id\": flight_id\n", + " }\n", + " for platform_id, flights in meta.items()\n", + " for flight_id, flight in flights.items()\n", + " for s in flight[\"segments\"]\n", + " ]\n", + "\n", + "segments_by_segment_id = {s[\"segment_id\"]: s for s in segments}\n", + "segments_ordered_by_start_time = list(sorted(segments, key=lambda s: s[\"start\"]))\n", + "\n", + "circles_Jan24 = [s[\"segment_id\"]\n", + " for s in segments_ordered_by_start_time\n", + " if \"circle\" in s[\"kinds\"]\n", + " and s[\"start\"].date() == date(2020,1,24)\n", + " and s[\"platform_id\"] == \"HALO\"\n", + " ]\n", + "\n", + "first_circle_Jan24 = circles_Jan24[0]\n", + "\n", + "dropsonde_ids = segments_by_segment_id[first_circle_Jan24][\"dropsondes\"][\"GOOD\"]" + ] + }, + { + "cell_type": "markdown", + "id": "bde5f87f-503d-48f2-9cfc-8a4a8c9010a9", + "metadata": {}, + "source": [ + "Load dropsonde dataset and select first circle of January 24, 2020:" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "5f08a63c-b9e8-4d9e-b48c-ae27983a3920", + "metadata": { + "tags": [] + }, + "outputs": [], + "source": [ + "ds = cat.dropsondes.JOANNE.level3.to_dask()\n", + "\n", + "# Select dropsondes from Jan 24 2020\n", + "mask_sondes_first_circle_Jan24 = reduce(lambda a, b: a | b, [ds.sonde_id==d for d in dropsonde_ids])\n", + "mask_sondes_third_circle_Jan24 = reduce(lambda a, b: a | b, [ds.sonde_id==d for d in dropsonde_ids])\n", + "\n", + "ds_sondes_first_circle_Jan24 = ds.isel(sonde_id=mask_sondes_first_circle_Jan24)\n", + "ds_sondes_third_circle_Jan24 = ds.isel(sonde_id=mask_sondes_third_circle_Jan24)" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "14d71026-054c-410f-8723-b44adcf31f0b", + "metadata": {}, + "outputs": [], + "source": [ + "ds_sondes_first_circle_Jan24" + ] + }, + { + "cell_type": "markdown", + "id": "f01f4463-b18a-4562-ac57-0278a5d27a19", + "metadata": {}, + "source": [ + "Access satellite images:" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "6de0f717-bd26-4eec-8eba-6f43ee98a57a", + "metadata": {}, + "outputs": [], + "source": [ + "satellite_image = ql.get_satellite_data(ds_flight=ds_sondes_first_circle_Jan24,\n", + " satellite_name=\"goes16\", channel = 13, product=\"ABI-L2-CMIPF\")#,\n", + " #extent = (-62,-48,10,30))" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "5118c746-9a84-4323-b3bc-c1b3c2c61793", + "metadata": {}, + "outputs": [], + "source": [ + "satellite_image" + ] + }, + { + "cell_type": "markdown", + "id": "b5834ca1-2c9f-47c7-bf7b-0c8ce62eb897", + "metadata": {}, + "source": [ + "## Plot dropsonde launch locations on a map (over a GOES-16 image)\n", + "\n", + "By default the satellite image will be taken at the average launch time from all dropsondes in the plot. By default the colormap is chosen to show the flight altitude, but it can be any variable with dimensions (`sonde_id`)." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "652d22a6-0ff9-45f8-8b69-8608f18d7ae5", + "metadata": {}, + "outputs": [], + "source": [ + "ql.launch_locations_map(ds_flight=ds_sondes_first_circle_Jan24, save_filepath=\"/home/m/m300931/\", \n", + " satellite_data=satellite_image)\n", + " #extent = (-62,-48,10,30))" + ] + }, + { + "cell_type": "markdown", + "id": "61c0e563-6411-41c1-a21e-d4d23c064537", + "metadata": {}, + "source": [ + "## Plot latitude-time quantities\n", + "\n", + "By default the colormap is chosen to show the flight altitude, but it can be any variable with dimensions (`sonde_id`)." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "5dc93114-76fc-4ce6-88b9-dab6944428f6", + "metadata": {}, + "outputs": [], + "source": [ + "ql.plot_lat_time(ds_flight=ds_sondes_first_circle_Jan24, save_filepath=\"/home/m/m300931/\")" + ] + }, + { + "cell_type": "markdown", + "id": "e7005a0b-8603-4a8e-8306-af5933ce2a5a", + "metadata": {}, + "source": [ + "## Plot vertical profiles\n", + "The variables and number of plots can be adjusted by the user. By default the quicklook shows temperature, potential temperature, relative humidity, and wind speed." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "0d1315af-7d6e-499e-a1ce-2cda46562724", + "metadata": {}, + "outputs": [], + "source": [ + "ql.plot_profiles(ds_flight=ds_sondes_first_circle_Jan24, save_filepath=\"/home/m/m300931/\")" + ] + }, + { + "cell_type": "markdown", + "id": "c54961a5-873d-4091-b798-d72e42d06a45", + "metadata": {}, + "source": [ + "## Plot dropsonde drift" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "805d161e-72bd-4596-948d-75842b6df1e2", + "metadata": {}, + "outputs": [], + "source": [ + "ql.drift_plots(ds_flight=ds_sondes_first_circle_Jan24, save_filepath=\"/home/m/m300931/\")" + ] + }, + { + "cell_type": "markdown", + "id": "31d1e2b6-0bb0-4b03-b94f-cd0619b7abe7", + "metadata": {}, + "source": [ + "## Saving plots in one pdf" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "4b542d28-3a83-4801-b767-316959ede93d", + "metadata": {}, + "outputs": [], + "source": [] + } + ], + "metadata": { + "kernelspec": { + "display_name": "HALO-DROPS", + "language": "python", + "name": "my-kernel" + }, + "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.10.11" + } + }, + "nbformat": 4, + "nbformat_minor": 5 +} diff --git a/src/halodrops/api/plotting.py b/src/halodrops/api/plotting.py new file mode 100644 index 0000000..9c32718 --- /dev/null +++ b/src/halodrops/api/plotting.py @@ -0,0 +1,19 @@ +from halodrops.plotting import quicklooks as ql + + +# Access satellite data +satellite_image = ql.get_satellite_data() + +# Plot launch locations over satellite images +ql.launch_locations_map() + +# Plot longitude/time +ql.plot_lat_time() + +# Plot vertical profiles +ql.plot_profiles() + +# Plot dropsonde drift +ql.drift_plots() + +# Output all quicklooks into a PDF diff --git a/src/halodrops/plotting/quicklooks.py b/src/halodrops/plotting/quicklooks.py new file mode 100644 index 0000000..852ae36 --- /dev/null +++ b/src/halodrops/plotting/quicklooks.py @@ -0,0 +1,369 @@ +from halodrops import sonde +from halodrops.helper import paths +from halodrops.qc import profile + +from importlib import reload +from gogoesgone import processing as pr +from gogoesgone import zarr_access as za + +reload(pr) +reload(za) + +import xarray as xr +import numpy as np +import cartopy as cp +import cartopy.crs as ccrs +from cartopy.mpl.gridliner import LONGITUDE_FORMATTER, LATITUDE_FORMATTER +from cartopy.feature import LAND +import cartopy.feature as cfeature +from matplotlib.offsetbox import AnchoredText +import matplotlib.pyplot as plt +import matplotlib.ticker as mticker +from mpl_toolkits.axes_grid1 import make_axes_locatable +import matplotlib.dates as mdates +from datetime import datetime, date +import s3fs +import pandas as pd + + +def convert_time_to_str(time=None, time_format="%Y%m%d %H:%M:%S"): + + """ + Convert input time into desired string format. + """ + + # Ensure time is in correct format + timestamp = (time - np.datetime64("1970-01-01T00:00:00")) / np.timedelta64(1, "s") + datetime_time = datetime.utcfromtimestamp(timestamp) + str_time = datetime_time.strftime(time_format) + + return str_time + + +def get_mean_launch_time(ds_flight=None, time_format="%Y%m%d %H:%M:%S"): + + """ + Compute mean launch time from all sondes in the dataset. + """ + + mean_time = convert_time_to_str(ds_flight.launch_time.mean().values, time_format) + + return mean_time + + +def get_satellite_data( + satellite_time="mean_launch_time", + time_format="%Y%m%d %H:%M:%S", + ds_flight=None, + satellite_name="goes16", + channel=13, + product="ABI-L2-CMIPF", + extent=None, +): + + """ + Access satellite data for nearest time, map to lon/lat grid, and convert to dataset. + By default use the mean launch time from dropsonde dataset. + """ + + # Get correct time for satellite data. Default is mean of all launches + if satellite_time == "mean_launch_time": + use_time = get_mean_launch_time(ds_flight=ds_flight) + else: + use_time = convert_time_to_str(time=satellite_time) + + # Define default extent based on launch locations, or user-defined extent + if extent is not None: + image_lims = extent + else: + image_lims = ( + ds_flight.lon.min() - 1, + ds_flight.lon.max() + 1, + ds_flight.lat.min() - 1, + ds_flight.lat.max() + 1, + ) + + # Get filepath to satellite data at nearest time. + flist = za.nearest_time_url(use_time, time_format, channel, product, satellite_name) + m = za.get_mapper_from_mzz(flist) + + # Select subset of satellite domain + img = pr.Image(m) + subset = img.subset_region_from_latlon_extents(image_lims, unit="degree") + + return subset + + +def launch_locations_map( + ds_flight=None, + satellite_data=None, + save_filepath="/path/to/save/", + color_coding_var="flight_altitude", + color_coding_cmap="magma", + satellite_time=None, + extent=None, + satellite_cmap="Greys", + satellite_vmin=280, + satellite_vmax=300, +): + + """ + Plot dropsonde launch locations, optionally over satellite images. + """ + + fig = plt.figure(figsize=(10, 8)) + + # Define default extent based on launch locations, or user-defined extent + if extent is not None: + image_lims = extent + else: + image_lims = ( + ds_flight.lon.min() - 1, + ds_flight.lon.max() + 1, + ds_flight.lat.min() - 1, + ds_flight.lat.max() + 1, + ) + + ax = plt.axes(projection=ccrs.PlateCarree()) + ax.coastlines(resolution="50m", linewidth=1.5) + ax.set_extent(image_lims, crs=ccrs.PlateCarree()) + + # Plot satellite image + if satellite_data: + sat_im = satellite_data.CMI.isel(t=0).plot( + ax=ax, + x="lon", + y="lat", + cmap=satellite_cmap, + add_colorbar=True, + cbar_kwargs={ + "pad": 0.05, + "extend": "both", + "aspect": 15, + "shrink": 0.7, + "label": f"{satellite_data.CMI.name} / {satellite_data.CMI.units}", + }, + vmin=satellite_vmin, + vmax=satellite_vmax, + zorder=-1, + transform=ccrs.PlateCarree(), + ) + + # Plot flight path + ax.plot( + ds_flight["lon"].isel(alt=-700), + ds_flight["lat"].isel(alt=-700), + c="red", + linestyle=":", + transform=ccrs.PlateCarree(), + zorder=1, + ) + + # Plot launch locations + im_launches = ax.scatter( + ds_flight["lon"].isel(alt=-700), + ds_flight["lat"].isel(alt=-700), + marker="o", + edgecolor="grey", + s=60, + transform=ccrs.PlateCarree(), + c=ds_flight[color_coding_var], + cmap=color_coding_cmap, + ) + + # Assigning axes ticks + xticks = np.arange(-180, 180, 1) + yticks = np.arange(-90, 90, 1) + + # Setting up the gridlines + gl = ax.gridlines( + crs=ccrs.PlateCarree(), + draw_labels=True, + linewidth=1, + color="gray", + alpha=0.2, + linestyle="--", + ) + gl.xlocator = mticker.FixedLocator(xticks) + gl.ylocator = mticker.FixedLocator(yticks) + gl.xformatter = LONGITUDE_FORMATTER + gl.yformatter = LATITUDE_FORMATTER + gl.xlabel_style = {"size": 10, "color": "k"} + gl.ylabel_style = {"size": 10, "color": "k"} + + # Colorbar + g = fig.colorbar( + im_launches, orientation="horizontal", extend="both", aspect=30, pad=0.05 + ) + g.set_label( + f"{ds_flight[color_coding_var].name} / {ds_flight[color_coding_var].units}", + fontsize=12, + ) + plt.tick_params(labelsize=10) + + # Format time stamp + title_time = convert_time_to_str( + time=satellite_data.t[0].values, time_format="%Y-%m-%d %H:%M" + ) + ax.set_title("") + + plt.title( + f"Sondes {ds_flight.sonde_id.values[0]} to {ds_flight.sonde_id.values[-1]} (Satellite Time = {title_time})", + fontsize=12, + pad=10, + ) + + # Save figure + save_filename = f"{save_filepath}launch-locations-{color_coding_var}-{satellite_data.platform_ID}.png" + plt.savefig(save_filename, dpi=300, bbox_inches="tight") + + +def plot_lat_time( + ds_flight, + color_coding_var="flight_altitude", + color_coding_cmap="magma", + save_filepath="/path/to/save/", +): + + """ + Plot spatio-temporal variation (lat v/s time) of selected variable. + """ + + ax = plt.figure(figsize=(12, 4)) + plt.scatter( + ds_flight["launch_time"].values, + ds_flight["lat"].isel(alt=-700).values, + s=90, + c=ds_flight[color_coding_var], + edgecolor="grey", + cmap=color_coding_cmap, + ) + plt.xlim( + np.min(ds_flight["launch_time"].values) - np.timedelta64(4, "m"), + np.max(ds_flight["launch_time"].values) + np.timedelta64(4, "m"), + ) + plt.gca().spines["right"].set_visible(False) + plt.gca().spines["top"].set_visible(False) + g = plt.colorbar() + g.set_label( + f"{ds_flight[color_coding_var].name} / {ds_flight[color_coding_var].units}", + fontsize=12, + ) + + myFmt = mdates.DateFormatter("%H:%M") + plt.gca().xaxis.set_major_formatter(myFmt) + plt.xlabel("Time / UTC", fontsize=12) + plt.ylabel("Latitude / $\degree$N", fontsize=12) + plt.title( + f"Sondes {ds_flight.sonde_id.values[0]} to {ds_flight.sonde_id.values[-1]}", + fontsize=12, + ) + + save_filename = f"{save_filepath}spatiotemporal-variation-{color_coding_var}.png" + + plt.savefig( + save_filename, + dpi=300, + bbox_inches="tight", + ) + + +def plot_profiles( + ds_flight, + r=["ta", "theta", "rh", "wspd", "wdir"], + r_titles=[ + "T / $\degree$C", + "$\\theta$ / K", + "RH / %", + "Wind speed / ms$^{-1}$", + "Wind direction / $\degree$", + ], + row=1, + col=4, + save_filepath="/path/to/save/", +): + + """ + Plot vertical profiles of specified variables. + """ + + f, ax = plt.subplots(row, col, sharey=True, figsize=(12, 6)) + + for j in range(col): + d = ds_flight[r[j]] + for i in range(1, len(ds_flight["launch_time"]) - 1): + ax[j].plot( + d.isel(sonde_id=i), + ds_flight["alt"] / 1000, + c="grey", + alpha=0.25, + linewidth=0.5, + ) + + ax[j].plot( + np.nanmean(d, axis=0), + ds_flight["alt"] / 1000, + linewidth=3, + c="k", + ) + ax[j].set_xlabel(r_titles[j], fontsize=12) + ax[j].spines["right"].set_visible(False) + ax[j].spines["top"].set_visible(False) + if j == 0: + ax[j].set_ylabel("Altitude / km", fontsize=12) + + plt.suptitle( + f"Sondes {ds_flight.sonde_id.values[0]} to {ds_flight.sonde_id.values[-1]}", + fontsize=12, + ) + + save_filename = f"{save_filepath}vertical-profiles-measured-quantities.png" + + plt.savefig(save_filename, dpi=300, bbox_inches="tight") + + +def drift_plots(ds_flight=None, save_filepath="/path/to/save/"): + + print("Plotting drift in lat and lon...") + + f, ax = plt.subplots(1, 2, sharey=True, figsize=(10, 5)) + + for i in range(len(ds_flight["launch_time"])): + + max_id = np.max(np.where(~np.isnan(ds_flight["lon"].isel(sonde_id=i)))) + + ax[0].plot( + ds_flight["lat"].isel(sonde_id=i) + - ds_flight["lat"].isel(sonde_id=i).isel(alt=max_id), + ds_flight["alt"] / 1000, + linewidth=1.5, + c="grey", + alpha=0.75, + ) + + ax[0].set_xlabel("Drift in Latitude / $\degree$", fontsize=12) + ax[0].set_ylabel("Altitude / km", fontsize=12) + ax[0].spines["right"].set_visible(False) + ax[0].spines["top"].set_visible(False) + + ax[1].plot( + ds_flight["lon"].isel(sonde_id=i) + - ds_flight["lon"].isel(sonde_id=i).isel(alt=max_id), + ds_flight["alt"] / 1000, + linewidth=1.5, + c="grey", + alpha=0.75, + ) + + ax[1].set_xlabel("Drift in Longitude / $\degree$", fontsize=12) + ax[1].spines["right"].set_visible(False) + ax[1].spines["top"].set_visible(False) + + plt.suptitle( + f"Sondes {ds_flight.sonde_id.values[0]} to {ds_flight.sonde_id.values[-1]}", + fontsize=12, + ) + + save_filename = f"{save_filepath}drift-in-lat-lon.png" + + plt.savefig(save_filename, dpi=300, bbox_inches="tight")