Skip to content

Commit

Permalink
Merge pull request #47 from d70-t/cleanup
Browse files Browse the repository at this point in the history
Some cleanups
  • Loading branch information
tmieslinger authored May 3, 2021
2 parents 532191a + 64dfe9d commit fb6b5ae
Show file tree
Hide file tree
Showing 14 changed files with 206 additions and 208 deletions.
1 change: 1 addition & 0 deletions .gitignore
Original file line number Diff line number Diff line change
@@ -1,3 +1,4 @@
_build
.ipynb_checkpoints
*.swp
.DS_Store
37 changes: 18 additions & 19 deletions how_to_eurec4a/bacardi.md
Original file line number Diff line number Diff line change
Expand Up @@ -17,10 +17,6 @@ The following script exemplifies the access and usage of the Broadband AirCrAft

The dataset is published under [Ehrlich et al. (2021)](https://doi.org/10.25326/160). If you have questions or if you would like to use the data for a publication, please don't hesitate to get in contact with the dataset authors as stated in the dataset attributes `contact` or `author`.

```{code-cell} ipython3
%pylab inline
```

## Get data
* To load the data we first load the EUREC4A meta data catalogue. More information on the catalog can be found [here](https://github.com/eurec4a/eurec4a-intake#eurec4a-intake-catalogue).

Expand All @@ -35,7 +31,9 @@ list(cat.HALO.BACARDI)

* We can further specify the platform, instrument, if applicable dataset level or variable name, and pass it on to dask.

*Note: have a look at the attributes of the xarray dataset `ds` for all relevant information on the dataset, such as author, contact, or citation infromation.*
```{note}
Have a look at the attributes of the xarray dataset `ds` for all relevant information on the dataset, such as author, contact, or citation infromation.
```

```{code-cell} ipython3
ds = cat.HALO.BACARDI.irradiances['HALO-0205'].to_dask()
Expand All @@ -53,24 +51,24 @@ The data from EUREC4A is of 10 Hz measurement frequency and corrected for dynami
We plot the upward and downward irradiances in two panels.

```{code-cell} ipython3
mpl.rcParams['font.size'] = 12
%matplotlib inline
import matplotlib.pyplot as plt
plt.style.use(["./mplstyle/book", "./mplstyle/wide"])
```

fig, (ax1, ax2) = plt.subplots(1,2, figsize=(16,6))
```{code-cell} ipython3
fig, (ax1, ax2) = plt.subplots(1, 2)
ax1.set_prop_cycle(color=['darkblue', 'red'])
for var in ['F_up_solar', 'F_up_terrestrial']:
ds[var].plot(ax=ax1,label= var)
ax1.legend()
ax1.set_ylabel('upward solar and terrestrial irradiance / Wm$^{-2}$')
ax2.set_prop_cycle(color=['grey', 'darkblue', 'skyblue'])
for var in ['F_down_solar_sim', 'F_down_solar', 'F_down_solar_diff']:
ax2.set_prop_cycle(color=['grey', 'skyblue', 'darkblue'])
for var in ['F_down_solar_sim', 'F_down_solar_diff', 'F_down_solar']:
ds[var].plot(ax=ax2, label= var)
ax2.legend()
ax2.set_ylabel('downward solar irradiance / Wm$^{-2}$')
for ax in [ax1, ax2]:
ax.spines['right'].set_visible(False)
ax.spines['top'].set_visible(False)
ax2.set_ylabel('downward solar irradiance / Wm$^{-2}$');
```

The attitude correction of downward solar irradiance does not account for the present cloud situation above HALO. Instead, two data sets, one assuming cloud-free and one assuming overcast (diffuse illumination) conditions, are provided. Depending on the application, the user needs to choose between both data sets. For the downward solar irradiance assuming cloud-free conditions, the data are filtered for turns of HALO, high roll and pitch angles. This filter is not applied for the data assuming overcast/diffuse conditions to provide the full data. However, data during turns of HALO need to be analysed with care. As shown in the example some artifical spikes due to turns are present in the data.
Expand All @@ -81,9 +79,10 @@ The wiggles originate from the about 200km change in location and therewith sol

```{code-cell} ipython3
fig, ax = plt.subplots()
ds.F_down_solar.plot(ax=ax, color = 'darkblue')
ax.set_ylabel('downward solar irradiance \n corrected for cloud-free conditions / Wm$^{-2}$', color = 'darkblue')
ax2=ax.twinx()
ds.sza.plot(ax=ax2, color = 'black')
ax2.set_ylim(110,28)
ds.F_down_solar.plot(ax=ax, color='darkblue')
ax.set_ylabel('downward solar irradiance \n corrected for cloud-free conditions / Wm$^{-2}$',
color='darkblue')
ax2 = ax.twinx()
ds.sza.plot(ax=ax2, color='black')
ax2.set_ylim(110, 28);
```
95 changes: 48 additions & 47 deletions how_to_eurec4a/dropsondes.md
Original file line number Diff line number Diff line change
Expand Up @@ -19,28 +19,21 @@ during EUREC4A - ATOMIC.
More information on the dataset can be found at https://github.com/Geet-George/JOANNE/tree/master/joanne/Level_3#level---3.
If you have questions or if you would like to use the data for a publication, please don't hesitate to get in contact with the dataset authors as stated in the dataset attributes `contact` and `author`.

```{code-cell} ipython3
%pylab inline
```

```{code-cell} ipython3
import logging
from functools import reduce
import eurec4a
```

## Get data
* To load the data we first load the EUREC4A meta data catalogue. More information on the catalog can be found [here](https://github.com/eurec4a/eurec4a-intake#eurec4a-intake-catalogue).

```{code-cell} ipython3
import datetime
import numpy as np
import eurec4a
cat = eurec4a.get_intake_catalog()
list(cat)
```

* We can funrther specify the platform, instrument, if applicable dataset level or variable name, and pass it on to dask.

*Note: have a look at the attributes of the xarray dataset `ds` for all relevant information on the dataset, such as author, contact, or citation infromation.*
```{note}
Have a look at the attributes of the xarray dataset `ds` for all relevant information on the dataset, such as author, contact, or citation infromation.
```

```{code-cell} ipython3
ds = cat.dropsondes.JOANNE.level3.to_dask()
Expand Down Expand Up @@ -93,84 +86,92 @@ dropsonde_ids
We transfer the information from our flight segment selection to the dropsondes data in the xarray dataset.

```{code-cell} ipython3
from functools import reduce
mask_sondes_first_circle_Feb05 = reduce(lambda a, b: a | b, [ds.sonde_id==d
for d in dropsonde_ids])
ds_sondes_first_circle_Feb05 = ds.isel(sounding=mask_sondes_first_circle_Feb05)
```

## Plots
You can get a list of available variables in the dataset from `ds.variables.keys()`
*Note: fetching the data and displaying it might take a few seconds*

```{note}
fetching the data and displaying it might take a few seconds
```

```{code-cell} ipython3
import matplotlib as mpl
%matplotlib inline
import matplotlib.pyplot as plt
mpl.rcParams['font.size'] = 14
plt.style.use("./mplstyle/book")
```

* Figure 1: Temperature and relative humidity as stored in the xarray dataset
### Temperature and relative humidity as stored in the xarray dataset

```{code-cell} ipython3
fig, (ax0, ax1) = plt.subplots(1, 2, figsize=(12,7))
ds_sondes_first_circle_Feb05.ta.T.plot(ax=ax0, cmap="magma")
ds_sondes_first_circle_Feb05.rh.T.plot(ax=ax1, cmap="BrBG")
fig, (ax0, ax1) = plt.subplots(1, 2)
ds_sondes_first_circle_Feb05.ta.transpose("alt", "sounding").plot(ax=ax0, cmap="magma")
ds_sondes_first_circle_Feb05.rh.transpose("alt", "sounding").plot(ax=ax1, cmap="BrBG")
None
```

* Figure 2: Temperature and relative humidity profiles.
### Temperature and relative humidity profiles.
The temperature profiles are colored according to their launch time, while relative humidity profiles are colored according to their integrated water vapour in the measured column, i.e. precipitable water.

```{code-cell} ipython3
def dt64_to_dt(dt64):
return datetime.datetime.utcfromtimestamp((dt64 - np.datetime64('1970-01-01T00:00:00'))
/ np.timedelta64(1, 's'))
epoch = np.datetime64('1970-01-01T00:00:00')
second = np.timedelta64(1, 's')
return datetime.datetime.utcfromtimestamp(int((dt64 - epoch) / second))
```

```{code-cell} ipython3
fig, (ax0, ax1) = plt.subplots(1, 2, figsize=(12,7))
fig, (ax0, ax1) = plt.subplots(1, 2)
y = ds_sondes_first_circle_Feb05.alt
x0 = ds_sondes_first_circle_Feb05.ta
ax0.set_prop_cycle(color=plt.get_cmap("viridis")(np.linspace(0, 1, len(dropsonde_ids))))
ax0.plot(x0.T, y.data[:, np.newaxis])
x0 = ds_sondes_first_circle_Feb05.ta.transpose("alt", "sounding")
ax0.set_prop_cycle(color=plt.cm.viridis(np.linspace(0, 1, len(dropsonde_ids))))
ax0.plot(x0, y.data[:, np.newaxis])
ax0.set_xlabel(f"{x0.long_name} / {x0.units}")
ax0.set_ylabel(f"{y.name} / m")
ax0.legend([dt64_to_dt(d).strftime("%H:%M:%S")
for d in ds_sondes_first_circle_Feb05.launch_time],
title=x0.launch_time.name)
x1 = ds_sondes_first_circle_Feb05.rh
x1 = ds_sondes_first_circle_Feb05.rh.transpose("alt", "sounding")
c = ds_sondes_first_circle_Feb05.PW
ax1.set_prop_cycle(color=mpl.cm.BrBG([int(i) for i in (c - c.min())#cividis_r
ax1.set_prop_cycle(color=plt.cm.BrBG([int(i) for i in (c - c.min())#cividis_r
/ (c - c.min()).max() * 255]))
p = ax1.plot(x1.T, y.data[:, np.newaxis])
p = ax1.plot(x1, y.data[:, np.newaxis])
ax1.set_xlabel(f"{x1.long_name} / {x1.units}")
ax1.set_ylabel(f"{y.name} / m")
ax1.legend(ds_sondes_first_circle_Feb05.PW.values,
title=ds_sondes_first_circle_Feb05.PW.standard_name)
for ax in [ax0, ax1]:
ax.spines['right'].set_visible(False)
ax.spines['top'].set_visible(False)
fig.suptitle('Dropsondes from 1st circle an February 5', fontsize=18)
None
```

* Figure 3: wind speed variations throughout February 5
### wind speed variations throughout February 5

```{code-cell} ipython3
mask_sondes_Feb05 = [dt64_to_dt(t).date() == datetime.date(2020,2,5)
for t in ds.launch_time.values]
mask_sondes_Feb05 = ds.launch_time.astype("<M8[D]") == np.datetime64("2020-02-05")
ds_sondes_Feb05 = ds.isel(sounding=mask_sondes_Feb05)
```

```{code-cell} ipython3
fig, ax = plt.subplots(figsize=(12,4))
p = ax.pcolormesh(ds_sondes_Feb05.wspd.T, vmin=0, vmax=20)#var.T
plt.colorbar(p, ax=ax, label=f"{ds_sondes_Feb05.wspd.long_name} / {ds_sondes_Feb05.wspd.units}")
xticks = np.arange(0, ds_sondes_Feb05.sounding.size, int(ds_sondes_Feb05.sounding.size/10))
ax.set_xticks(xticks)
ax.set_xticklabels([dt64_to_dt(t).strftime("%H:%M") for t in ds.launch_time.values[xticks]])
ax.set_xlabel(f"{ds_sondes_Feb05.launch_time.long_name}")
yticks = np.arange(0, ds_sondes_Feb05.alt.size, int(ds_sondes_Feb05.alt.size/5))
ax.set_yticks(yticks)
ax.set_yticklabels(ds_sondes_Feb05.alt.values[yticks])
ax.set_ylabel(f"{ds_sondes_Feb05.alt.name} / m")
with plt.style.context("mplstyle/wide"):
fig, ax = plt.subplots()
p = ax.pcolormesh(ds_sondes_Feb05.wspd.transpose("alt", "sounding"), vmin=0, vmax=20)
plt.colorbar(p, ax=ax, label=f"{ds_sondes_Feb05.wspd.long_name} / {ds_sondes_Feb05.wspd.units}")
xticks = np.arange(0, ds_sondes_Feb05.sounding.size, int(ds_sondes_Feb05.sounding.size/10))
ax.set_xticks(xticks)
ax.set_xticklabels([dt64_to_dt(t).strftime("%H:%M") for t in ds.launch_time.values[xticks]])
ax.set_xlabel(f"{ds_sondes_Feb05.launch_time.long_name}")
yticks = np.arange(0, ds_sondes_Feb05.alt.size, int(ds_sondes_Feb05.alt.size/5))
ax.set_yticks(yticks)
ax.set_yticklabels(ds_sondes_Feb05.alt.values[yticks])
ax.set_ylabel(f"{ds_sondes_Feb05.alt.name} / m")
None
```
65 changes: 30 additions & 35 deletions how_to_eurec4a/flight_tracks_leaflet.md
Original file line number Diff line number Diff line change
Expand Up @@ -20,45 +20,30 @@ This notebook shows how to plot the flight tracks of the HALO aircraft. First, w
First, we'll import pylab and the EUREC4A data catalog.

```{code-cell} ipython3
%pylab inline
import eurec4a
cat = eurec4a.get_intake_catalog()
```

### Restructuring the dataset
When inspecting the contents of the currently available BAHAMAS quicklook dataset, we see that it does not handle coordinates according to the ideas of [CF conventions](http://cfconventions.org/Data/cf-conventions/cf-conventions-1.8/cf-conventions.html), which makes using the data a lot more difficult.
### Inspecting the dataset

```{code-cell} ipython3
ds = cat.HALO.BAHAMAS.QL['HALO-0122'].to_dask()
ds = cat.HALO.BAHAMAS.PositionAttitude['HALO-0122'].to_dask()
ds
```

In order to fix that a bit, we could come up with a little function which transforms the dataset we get into a dataset we would like to have. Of course the proper way of handling this would be to actually fix the dataset, but this way of doing it may still be illustrative.
Another advantage of this preparatory step is that amount of data in the dataset can be reduced, such that when we later load the datasets early, less data will be fetched.

```{code-cell} ipython3
def fix_halo_ql(ds):
import xarray as xr
ds = ds.rename({"tid":"time"})
return xr.Dataset({
"time": ds.TIME,
"lat": ds.IRS_LAT,
"lon": ds.IRS_LON,
"alt": ds.IRS_ALT,
})
```

This at least has a `time` coordinate.

```{code-cell} ipython3
dsfixed = fix_halo_ql(ds)
dsfixed
```

Just to have a better visual impression, we can create a quick overview plot of the data:

```{code-cell} ipython3
plt.plot(dsfixed.lon, dsfixed.lat)
%matplotlib inline
import numpy as np
import matplotlib.pyplot as plt
plt.style.use("./mplstyle/book")
plt.plot(ds.lon, ds.lat)
center_lat = float(ds.lat.mean())
plt.gca().set_aspect(1./np.cos(np.deg2rad(center_lat)))
plt.xlabel("longitude / °")
plt.ylabel("latitude / °")
plt.show()
```

Expand All @@ -82,20 +67,29 @@ def simplify_dataset(ds, tolerance):
We can now use that algorithm to generate a simplified version of the dataset with a tolerance of $10^{-5}$ degrees, which otherwise looks the same to the previous version.

```{code-cell} ipython3
dssimplified = simplify_dataset(dsfixed, 1e-5)
dssimplified = simplify_dataset(ds, 1e-5)
dssimplified
```

We can now compare those two tracks side by side while keeping a look at the number of data points required.

```{code-cell} ipython3
fig, (ax1, ax2) = subplots(1, 2, figsize=(15,4))
ax1.plot(dsfixed.lon, dsfixed.lat)
ax1.set_title(f"{len(dsfixed.time)} data points")
fig, (ax1, ax2) = plt.subplots(1, 2)
ax1.plot(ds.lon, ds.lat)
ax1.set_title(f"{len(ds.time)} data points")
ax2.plot(dssimplified.lon, dssimplified.lat)
ax2.set_title(f"{len(dssimplified.time)} data points")
for ax in (ax1, ax2):
ax.set_aspect(1./np.cos(np.deg2rad(center_lat)))
ax.set_xlabel("longitude / °")
ax.set_ylabel("latitude / °")
ax.label_outer()
plt.show()
ratio = len(dssimplified.time) / len(dsfixed.time)
ratio = len(dssimplified.time) / len(ds.time)
print(f"compression ratio: {ratio*100:.2f} %")
```

Expand Down Expand Up @@ -140,10 +134,10 @@ from multiprocessing.pool import ThreadPool
pool = ThreadPool(20)
def get_dataset(flight_id):
ds = cat.HALO.BAHAMAS.QL[flight_id].to_dask()
return flight_id, fix_halo_ql(ds).load()
ds = cat.HALO.BAHAMAS.PositionAttitude[flight_id].to_dask()
return flight_id, ds.load()
full_tracks = dict(pool.map(get_dataset, cat.HALO.BAHAMAS.QL))
full_tracks = dict(pool.map(get_dataset, cat.HALO.BAHAMAS.PositionAttitude))
```

We still have to simplify the dataset, which is done here:
Expand All @@ -156,6 +150,7 @@ tracks = {flight_id: simplify_dataset(ds, 1e-5)
Let's also quickly grab some colors from a matplotlib colorbar, such that we can show all tracks in individual colors:

```{code-cell} ipython3
import matplotlib
colors = [matplotlib.colors.to_hex(c)
for c in plt.cm.inferno(np.linspace(0, 1, len(tracks)))]
```
Expand Down
19 changes: 7 additions & 12 deletions how_to_eurec4a/hamp_comparison.md
Original file line number Diff line number Diff line change
Expand Up @@ -12,20 +12,12 @@ kernelspec:
---

# HAMP comparison

This notebook shows a comparison of cloud features as detected by different parts of the Halo Microwave Package (HAMP) and includes data from WALES and specMACS for additional context.
Another application of the HAMP data is included in the chapter {doc}`unified`.

```{code-cell} ipython3
import xarray
import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.axes_grid1.inset_locator import inset_axes, mark_inset
import glob
import eurec4a
%pylab inline
cat = eurec4a.get_intake_catalog()
```

Expand Down Expand Up @@ -179,10 +171,13 @@ radar_y = radar_y - wgs_correction

## Plot timeseries
```{code-cell} ipython3
%matplotlib inline
import matplotlib.pyplot as plt
plt.style.use("./mplstyle/book")
fig, (ax3, ax2, ax1, ax,) = plt.subplots(
nrows=4, figsize=(10, 8), sharex=True,
nrows=4, sharex=True,
gridspec_kw=dict(height_ratios=[1, 1, .75, 0.5]),
constrained_layout=True
)
ax.plot(wales.time, wales.cloud_mask + 0.2, '.', label='WALES', markersize=3)
Expand Down
Loading

0 comments on commit fb6b5ae

Please sign in to comment.