diff --git a/changelog/652.doc.rst b/changelog/652.doc.rst new file mode 100644 index 000000000..5933f6a4c --- /dev/null +++ b/changelog/652.doc.rst @@ -0,0 +1 @@ +Assorted clarifications and improvements to tutorial content. diff --git a/docs/examples/example_cryo_plots.md b/docs/examples/example_cryo_plots.md index 1ec60bd38..65e400c8a 100644 --- a/docs/examples/example_cryo_plots.md +++ b/docs/examples/example_cryo_plots.md @@ -72,7 +72,7 @@ Next we will need to calculate the pixel size in both Longitude and Latitude, wh First get a 2x2 grid of SkyCoord objects: ```{code-cell} ipython3 -space_2 = sp[:2,:2,0].axis_world_coords()[0] # [0] because we only care about space, not time +space_2 = sp[:2,:2,0].axis_world_coords("lon", "lat")[0] # [0] because this always returns a tuple ``` We calculate pixel size by calculating the step in both lat and lon. @@ -84,7 +84,7 @@ dlat = np.abs(space_2[0,1].Ty - space_2[0,0].Ty) aspect = dlon / dlat ``` -Finally, we can plot both our mean-subtracted scan at the peak wavelength and the mean values over all all wavelengths. +Finally, we can plot both our mean-subtracted scan at the peak wavelength and the mean values over all wavelengths. ```{code-cell} ipython3 fig = plt.figure(figsize=(7, 5), layout="constrained") @@ -145,7 +145,7 @@ plt.show() Now we can plot the slit position. We do this by taking the first raster step position of the SP dataset and computing the world coordinates of each pixel along the slit. ```{code-cell} ipython3 -slit_coords = sp[0,:,0].axis_world_coords()[0] # Again, [0] extracts the spatial coordinates and drops time. +slit_coords = sp[0,:,0].axis_world_coords("lon", "lat")[0] # Again, [0] because this always returns a tuple fig = plt.figure(layout="constrained") # Reuse the norm from the first plot @@ -162,7 +162,7 @@ plt.show() Notice how the slit has a larger field of view along the latitude dimension. We can re-plot the image, and crop the extent of the plot back to the extent of the context imager. ```{code-cell} ipython3 -slit_coords = sp[0,:,0].axis_world_coords()[0] # Again, [0] extracts the spatial coordinates and drops time. +slit_coords = sp[0,:,0].axis_world_coords("lon", "lat")[0] # Again, [0] because this always returns a tuple fig = plt.figure(layout="constrained") # Reuse the norm from the first plot diff --git a/docs/examples/reproject_vbi_mosaic.md b/docs/examples/reproject_vbi_mosaic.md index d0ecd6c3c..d837e3e34 100644 --- a/docs/examples/reproject_vbi_mosaic.md +++ b/docs/examples/reproject_vbi_mosaic.md @@ -38,7 +38,7 @@ from dkist.data.sample import VBI_L1_NZJTB ## Obtaining some data -In this example we will use the VBI sample dataset [AJQWW](https://dkist.data.nso.edu/datasetview/AJQWW). +In this example we will use the VBI sample dataset [L1_NZJTB](https://dkist.data.nso.edu/product/L1-NZJTB). If you want to replace this dataset with your own dataset, see {ref}`dkist:howto-guide:sample-data`. Let's load the data with {obj}`dkist.load_dataset`: diff --git a/docs/examples/vbi_extents.md b/docs/examples/vbi_extents.md index 11abaf3fb..e00996afc 100644 --- a/docs/examples/vbi_extents.md +++ b/docs/examples/vbi_extents.md @@ -16,6 +16,7 @@ kernelspec: # Showing the Field of View of VBI on AIA ```{code-cell} ipython3 +:tags: [keep-inputs] import numpy as np import matplotlib.pyplot as plt import astropy.units as u diff --git a/docs/tutorial/1_astropy_and_sunpy.md b/docs/tutorial/1_astropy_and_sunpy.md index 18c080e99..e7f6b272a 100644 --- a/docs/tutorial/1_astropy_and_sunpy.md +++ b/docs/tutorial/1_astropy_and_sunpy.md @@ -101,7 +101,7 @@ Currently, {obj}`astropy.coordinates` supports: * Spatial coordinates via {obj}`astropy.coordinates.SkyCoord` * Spectral coordinates via {obj}`astropy.coordinates.SpectralCoord` -* Stokes profiles via {obj}`astropy.coordinates.StokesCoord` (introduced in astropy 7.0) +* Stokes profiles via {obj}`astropy.coordinates.StokesCoord` (introduced in astropy 5.3) ### Spatial Coordinates @@ -142,6 +142,17 @@ There are few things to notice about the difference between these two `SkyCoord` 1. The `obstime` and `rsun` attributes are still present, but the `observer` attribute isn't. This is because heliographic coordinates are not observer dependent. 1. The `obstime` attribute is still important to transform to other frames, as the heliographic frame needs to know the location of Earth. ++++ + +Sunpy provides more equivalencies that might be useful here, for example if we wish to convert the solar angle values used for this SkyCoord into distance: + +```{code-cell} ipython3 +from sunpy.coordinates.utils import solar_angle_equivalency + +print((hpc1.Tx).to(u.km, equivalencies=solar_angle_equivalency(observer=hpc1.observer))) +print((hpc1.Ty).to(u.km, equivalencies=solar_angle_equivalency(observer=hpc1.observer))) +``` + ### Spectral Coordinates {obj}`astropy.coordinates.SpectralCoord` is a `Quantity`-like object which also holds information about the observer and target coordinates and relative velocities. diff --git a/docs/tutorial/3_dataset_dimensionality.md b/docs/tutorial/3_dataset_dimensionality.md index aa19f4a95..d270e9c5a 100644 --- a/docs/tutorial/3_dataset_dimensionality.md +++ b/docs/tutorial/3_dataset_dimensionality.md @@ -131,7 +131,7 @@ ds.wcs.pixel_axis_names note how these are reversed from one another, we can print them together with: ```{code-cell} ipython3 -for name, length in zip(ds.wcs.pixel_axis_names[::-1], ds.shape): +for name, length in zip(ds.wcs.pixel_axis_names, ds.wcs.pixel_shape): print(f"{name}: {length}") ``` diff --git a/docs/tutorial/exercises/exercises-1.md b/docs/tutorial/exercises/exercises-1.md new file mode 100644 index 000000000..fac4e0c9a --- /dev/null +++ b/docs/tutorial/exercises/exercises-1.md @@ -0,0 +1,328 @@ +--- +jupytext: + text_representation: + extension: .md + format_name: myst + format_version: 0.13 + jupytext_version: 1.17.1 +kernelspec: + display_name: Python 3 (ipykernel) + language: python + name: python3 +--- + ++++ {"editable": true, "slideshow": {"slide_type": "slide"}} + +# Exercise 1 + ++++ {"editable": true, "slideshow": {"slide_type": ""}} + +Load the sample VISP dataset. **Without actually slicing the dataset**, figure out how many files are in the following slices: + +- `ds[0]` +- `ds[:, :200]` +- `ds[:, :200, :200]` +- `ds[:, :, :200, :200]` + +```{code-cell} ipython3 +--- +editable: true +slideshow: + slide_type: skip +--- +from dkist import load_dataset +from dkist.data.sample import VISP_L1_KMUPT + +ds = load_dataset(VISP_L1_KMUPT) +``` + +```{code-cell} ipython3 +--- +editable: true +slideshow: + slide_type: skip +--- +ds.shape +``` + +```{code-cell} ipython3 +--- +editable: true +slideshow: + slide_type: skip +--- +ds.files +``` + +```{code-cell} ipython3 +--- +editable: true +slideshow: + slide_type: skip +--- +ds[0].files +``` + +```{code-cell} ipython3 +--- +editable: true +slideshow: + slide_type: skip +--- +ds[:, :200].files +``` + +```{code-cell} ipython3 +--- +editable: true +slideshow: + slide_type: skip +--- +ds[:, :200, :200].files +``` + +```{code-cell} ipython3 +--- +editable: true +slideshow: + slide_type: skip +--- +ds[:, :, :200, :200].files +``` + ++++ {"editable": true, "slideshow": {"slide_type": "slide"}} + +# Exercise 2 + ++++ + +Find the dataset from this year with the highest number of frames. + +- What is the dataset's Product ID? +- How many frames does the dataset have? +- Is the dataset embargoed? If yes, when does the embargo period end? + ++++ + +**Tips:** + +- You can pass `a.Provider("dkist")` to `Fido.search()` to search only for data provided by DKIST +- You can use `.keys()` or `.colnames` on a `DKISTQueryResponse` object to see all the column names + +```{code-cell} ipython3 +--- +editable: true +slideshow: + slide_type: notes +--- +from sunpy.net import attrs as a, Fido +import dkist.net +``` + +```{code-cell} ipython3 +--- +editable: true +slideshow: + slide_type: skip +--- +# Search for all DKIST data this year +res = Fido.search(a.Provider("dkist"), a.Time("2025-01-01", "2026-01-01")) +res +``` + +```{code-cell} ipython3 +--- +editable: true +slideshow: + slide_type: skip +--- +# Just get the DKIST results table and sort by number of frames, largest first +res2 = res['dkist'] +res2 = res[0] +res2 +``` + +```{code-cell} ipython3 +--- +editable: true +slideshow: + slide_type: skip +--- +res2.sort("Number of Frames", reverse=True) +``` + +```{code-cell} ipython3 +p = res2.copy() +print(p.sort("Number of Frames")) +``` + +```{code-cell} ipython3 +--- +editable: true +slideshow: + slide_type: skip +--- +# Print top result, ie dataset with largest number of frames +res2[:5]["Product ID", "Instrument", "Embargoed", "Number of Frames", "Embargo End Date"] +``` + +```{code-cell} ipython3 +max_frames = res[0]["Number of Frames"].max() +max_frames +``` + +```{code-cell} ipython3 +large = res[0][res[0]["Number of Frames"] == max_frames] +large +``` + +```{code-cell} ipython3 +Fido.fetch(large) +``` + ++++ {"editable": true, "slideshow": {"slide_type": "slide"}} + +# Exercise 3 + ++++ {"editable": true, "slideshow": {"slide_type": ""}} + +Open the VBI sample dataset again. Choose two adjacent tiles and calculate the overlap between them in arcseconds. + +HINTS: +- You will only need the first from of each tile. +- You may find the `SkyCoord.separation()` method helpful + +```{code-cell} ipython3 +--- +editable: true +slideshow: + slide_type: skip +--- +from dkist.data.sample import VBI_L1_NZJTB +from dkist import load_dataset +import astropy.units as u + +ds = load_dataset(VBI_L1_NZJTB) +``` + +```{code-cell} ipython3 +--- +editable: true +slideshow: + slide_type: skip +--- +# We only need the first image of each tile +first_tiles = ds.slice_tiles[0] +``` + +```{code-cell} ipython3 +--- +editable: true +slideshow: + slide_type: skip +--- +# The bottom left tile is 0,0, get the corner pixel on the right side +tile1_bottom_right_pix = first_tiles[0,0].wcs.array_index_to_world(0, 4096) +tile1_bottom_right_pix +``` + +```{code-cell} ipython3 +--- +editable: true +slideshow: + slide_type: skip +--- +# The bottom middle tile is 1,0, get the corner pixel on the left side +tile2_bottom_left_pix = first_tiles[1,0].wcs.array_index_to_world(0, 0) +tile2_bottom_left_pix +``` + +```{code-cell} ipython3 +--- +editable: true +slideshow: + slide_type: skip +--- +(tile1_bottom_right_pix.separation(tile2_bottom_left_pix).to(u.arcsec), + +tile1_bottom_right_pix.Tx - tile2_bottom_left_pix.Tx, + +tile1_bottom_right_pix.Ty - tile2_bottom_left_pix.Ty) +``` + ++++ {"editable": true, "slideshow": {"slide_type": "slide"}} + +# Exercise 4 + ++++ {"editable": true, "slideshow": {"slide_type": ""}} + +Using Fido, find the co-temporal VISP dataset, VBI dataset and AIA image from 19:47 on 2024/04/17. Use `sunpy.visualization.drawing.extent` to plot the FOV of the VBI and VISP datasets onto the AIA image. + +```{code-cell} ipython3 +--- +editable: true +slideshow: + slide_type: skip +--- +# Download one VISP and one VBI metadata ASDF file and a co-temporal AIA image for 19:47 on 2024-04-17 +# Using sunpy.visualization.drawing.extent plot the FOV of the VBi and VISP datasets on the AIA image. +import matplotlib.pyplot as plt + +import sunpy.map +from sunpy.net import Fido, attrs as a +from sunpy.visualization import drawing + +import dkist +import dkist.net +``` + +```{code-cell} ipython3 +--- +editable: true +slideshow: + slide_type: skip +--- +res = Fido.search(a.Time("2024-04-17T19:47:00", "2024-04-17T19:48:00"), + a.Instrument.vbi | a.Instrument.visp | a.Instrument.aia) +# Download the first result from each provider +files = Fido.fetch(res[:,0]) +files.sort() +``` + +```{code-cell} ipython3 +--- +editable: true +slideshow: + slide_type: skip +--- +# Unpack the filenames +vbi_file, visp_file, aia_file = files + +# Load VBI and VISP datasets +vbi, visp = dkist.load_dataset((vbi_file, visp_file)) + +# Load the AIA data +aia = sunpy.map.Map(aia_file) +``` + +```{code-cell} ipython3 +--- +editable: true +slideshow: + slide_type: skip +--- +# Create an axes with the AIA data's coordinate frame +fig = plt.figure(figsize=(8,8)) +ax = fig.add_subplot(projection=aia) + +# Plot the AIA data +aia.plot(axes=ax) + +# Draw the extents of the VBI and VISP data +drawing.extent(ax, vbi[0].wcs, color="C1", label="VBI") +drawing.extent(ax, visp[0,:,0,:].wcs, color="C2", label="VISP") + +plt.legend() + +_ = ax.axis((1000, 2000, 2600, 3400)) +``` diff --git a/docs/tutorial/exercises/exercises-2.md b/docs/tutorial/exercises/exercises-2.md new file mode 100644 index 000000000..959ad9470 --- /dev/null +++ b/docs/tutorial/exercises/exercises-2.md @@ -0,0 +1,97 @@ +--- +jupytext: + text_representation: + extension: .md + format_name: myst + format_version: 0.13 + jupytext_version: 1.17.1 +kernelspec: + display_name: Python 3 (ipykernel) + language: python + name: python3 +--- + ++++ {"editable": true, "slideshow": {"slide_type": "slide"}} + +# Exercise 1 + ++++ {"editable": true, "slideshow": {"slide_type": ""}} + +Take the Cryo-NIRSP CI sample data and download the first 50 frames and save an animation (mp4) of these frames. + +**Tips** + +1. To download more frames of the sample data, you should set the path in the call to `ds.files.download` to download to a different directory. +2. `Dataset.plot` will return a [`mpl-animators.ArrayAnimatorWCS`](https://docs.sunpy.org/projects/mpl-animators/en/stable/api/mpl_animators.ArrayAnimatorWCS.html) object when the dataset is >2D. This object has a `get_animation` method which returns a [`FuncAnimation`](https://matplotlib.org/stable/api/_as_gen/matplotlib.animation.FuncAnimation.html). + +```{code-cell} ipython3 +--- +editable: true +slideshow: + slide_type: skip +--- +import numpy as np +import matplotlib.pyplot as plt +from matplotlib.colors import PowerNorm + +import dkist +from dkist.data.sample import CRYO_L1_TJKGC, CRYO_L1_MSCGD +``` + +```{code-cell} ipython3 +--- +editable: true +slideshow: + slide_type: skip +--- +ci = dkist.load_dataset(CRYO_L1_MSCGD) # Load the context imager (CI) sample dataset +#ci[:50].files.download("~/dkist_data/CRYO_L1_MSCGD") +``` + +## Stage 2 + +```{code-cell} ipython3 +--- +editable: true +slideshow: + slide_type: skip +--- +fig = plt.figure() +ci.data = ci.data.compute() # Let's compute the data here to make things a little easier +vmin, vmax = np.nanpercentile(ci[0].data, [1,99]) +norm = PowerNorm(0.3, vmin=vmin, vmax=vmax) +ax = ci[:50].plot(norm=norm) +``` + +```{code-cell} ipython3 +--- +editable: true +slideshow: + slide_type: skip +--- +ani = ax.get_animation() +``` + +```{code-cell} ipython3 +ani +``` + +```{code-cell} ipython3 +--- +editable: true +slideshow: + slide_type: skip +--- +ani.save("animation.mp4") +``` + +```{code-cell} ipython3 +--- +editable: true +slideshow: + slide_type: skip +--- +from IPython.display import Video + +Video("animation.mp4", embed=True) +``` diff --git a/docs/tutorial/index.md b/docs/tutorial/index.md index a1f94bb33..180611f3d 100644 --- a/docs/tutorial/index.md +++ b/docs/tutorial/index.md @@ -34,4 +34,6 @@ Alternatively we have instructions for running the tutorial as Jupyter Notebooks 5_tileddataset 6_downloading_data 7_visualization +exercises/exercises-1 +exercises/exercises-2 ```