You signed in with another tab or window. Reload to refresh your session.You signed out in another tab or window. Reload to refresh your session.You switched accounts on another tab or window. Reload to refresh your session.Dismiss alert
Copy file name to clipboardExpand all lines: 06-raster-intro.md
+73-17Lines changed: 73 additions & 17 deletions
Display the source diff
Display the rich diff
Original file line number
Diff line number
Diff line change
@@ -173,7 +173,9 @@ downloading high-resolution images when only quick previews are required.
173
173
174
174
Overviews are often computed using powers of 2 as down-sampling (or zoom) factors. So, typically, the first level
175
175
overview (index 0) corresponds to a zoom factor of 2, the second level overview (index 1) corresponds to a zoom factor
176
-
of 4, and so on. Here, we open the third level overview (index 2, zoom factor 8) and check that the resolution is about 80 m:
176
+
of 4, and so on. According to [GDAL documentation](https://gdal.org/en/stable/drivers/raster/cog.html#drivers/raster/cog-co-RESAMPLING), the default resampling method for non-paletted images is `CUBIC`, for paletted images it is `NEAREST`.
177
+
178
+
Here, we open the third level overview (index 2, zoom factor 8) and check that the resolution is about 80 m:
177
179
178
180
```python
179
181
import rioxarray
@@ -217,6 +219,44 @@ More options can be consulted [here](https://docs.xarray.dev/en/v2024.02.0/gener
217
219
218
220
:::
219
221
222
+
:::challenge
223
+
## Exercise: visualize the NIR band
224
+
225
+
There are other bands in the downloaded Sentinel-2 scene in `data/sentinel2/` directory. In this exercise, let's look in to the NIR (near-infrared) band `data/sentinel2/nir.tif`. Please perform the following steps:
226
+
227
+
1. Load the NIR band using `rioxarray.open_rasterio()`, what is the dimension of the NIR band? What is the resolution?
228
+
2. Based on the dimension and resolution, load the NIR band again with the appropriate `overview_level` for visualization.
229
+
3. Use the `plot` method to visualize the NIR band. Do you see interesting patterns in the image? What do you think they are?
230
+
231
+
::::solution
232
+
233
+
We can use the same techniques that we have used for the red band to explore the NIR band, then inspect the dimension and resolution.
So it turned out that the NIR band has the same dimension and resolution as the red band. Therefore, we can load the NIR band with `overview_level=2` for visualization.
{alt="NIR band 80m resolution"}
254
+
255
+
As we can see in the image, there is a pattern of low values in NIR band on Rhodes Island, which is likely to be the burned area. We will further explore this in [episode 9](/episodes/09-raster-calculations.md).
256
+
257
+
::::
258
+
:::
259
+
220
260
## View Raster Coordinate Reference System (CRS) in Python
221
261
Another information that we're interested in is the CRS, and it can be accessed with `.rio.crs`. We introduced the concept of a CRS in [an earlier
222
262
episode](03-crs.md). Now we will see how features of the CRS appear in our data file and what
@@ -278,15 +318,7 @@ crs.area_of_use
278
318
AreaOfUse(west=24.0, south=0.0, east=30.0, north=84.0, name='Between 24°E and 30°E, northern hemisphere between equator and 84°N, onshore and offshore. Belarus. Bulgaria. Central African Republic. Democratic Republic of the Congo (Zaire). Egypt. Estonia. Finland. Greece. Latvia. Lesotho. Libya. Lithuania. Moldova. Norway. Poland. Romania. Russian Federation. Sudan. Svalbard. Türkiye (Turkey). Uganda. Ukraine.')
279
319
```
280
320
281
-
:::challenge
282
-
## Exercise: find the axes units of the CRS
283
-
What units are our data in? See if you can find a method to examine this information using `help(crs)` or `dir(crs)`
284
321
285
-
::::solution
286
-
`crs.axis_info` tells us that the CRS for our raster has two axis and both are in meters.
287
-
We could also get this information from the attribute `rhodes_red_80.rio.crs.linear_units`.
288
-
::::
289
-
:::
290
322
291
323
### Understanding pyproj CRS Summary
292
324
Let's break down the pieces of the `pyproj` CRS summary. The string contains all of the individual CRS elements that Python or another GIS might need, separated into distinct sections, and datum.
@@ -517,24 +549,48 @@ rhodes_visual.plot.imshow()
517
549
518
550
{alt="true-color image overview"}
519
551
520
-
Note that the `DataArray.plot.imshow()` function makes assumptions about the shape of the input DataArray, that since it has three channels, the correct colormap for these channels is RGB. It does not work directly on image arrays with more than 3 channels. One can replace one of the RGB channels with another band, to make a false-color image.
552
+
One can also specify the size of the plot, the aspect ratio, and use the `robust` option to set the color limits:
{alt="raster plot with correct aspect ratio"}
559
+
560
+
521
561
522
562
:::challenge
523
-
## Exercise: set the plotting aspect ratio
524
-
As seen in the figure above, the true-color image is stretched. Let's visualize it with the right aspect ratio. You can use the [documentation](https://xarray.pydata.org/en/stable/generated/xarray.DataArray.plot.imshow.html) of `DataArray.plot.imshow()`.
563
+
## Exercise: what is the correct order of the color channels?
525
564
526
-
::::solution
527
-
Since we know the height/width ratio is 1:1 (check the `rio.height` and `rio.width` attributes), we can set the aspect ratio to be 1. For example, we can choose the size to be 5 inches, and set `aspect=1`. Note that according to the [documentation](https://xarray.pydata.org/en/stable/generated/xarray.DataArray.plot.imshow.html) of `DataArray.plot.imshow()`, when specifying the `aspect` argument, `size` also needs to be provided.
565
+
We just visualized the true-color image `rhodes_visual`, which has three color channels: red, green, and blue. Apparently, the `plot.imshow()` function makes assumptions on the order of the channels in `rhodes_visual` to plot a true-color image. Can you figure out this order?
566
+
567
+
There are multiple ways to do this.
568
+
569
+
Hints:
570
+
571
+
1. You can use the [documentation](https://xarray.pydata.org/en/stable/generated/xarray.DataArray.plot.imshow.html) of `DataArray.plot.imshow()`, or
572
+
573
+
2. You can experiment this by yourself by mannually replacing values in one channel and see how the image changes.
574
+
575
+
:::solution
576
+
577
+
The order is red, green, blue.
578
+
579
+
To replace values of one channel, for example you can do:
528
580
529
581
```python
530
-
rhodes_visual.plot.imshow(size=5, aspect=1)
582
+
# The following code generates a red-ish image
583
+
# Therefore we know that the first channel is red.
584
+
rhodes_visual_modified = rhodes_visual.copy()
585
+
rhodes_visual_modified.values[0,:,:] =255# set the values to maximum
531
586
```
532
587
533
-
{alt="raster plot with correct aspect ratio"}
534
-
535
588
::::
536
589
:::
537
590
591
+
In conclusion, the `DataArray.plot.imshow()` function makes assumptions about the shape of the input DataArray, and the order of these channels. The correct colormap for these channels is RGB. It does not work directly on image arrays with more than 3 channels. One can also replace one of the RGB channels with another, to make a false-color image.
592
+
593
+
538
594
:::keypoints
539
595
-`rioxarray` and `xarray` are for working with multidimensional arrays like pandas is for working with tabular data.
540
596
-`rioxarray` stores CRS information as a CRS object that can be converted to an EPSG code or PROJ4 string.
0 commit comments