akankshabalha
commented
10 months ago Hi Mattia,
Could you please help me in this.
Closed akankshabalha closed 9 months ago
akankshabalha
commented
10 months ago Hi Mattia,
Could you please help me in this.
malmans2
commented
10 months ago Hi, other users requested support before you did, so I've been working on that.
The spatial resolution of both the datasets (noaa: 4km and proba: 1km) is different. Hence, I tried to regrid noaa dataset at similar resolution as that of proba. It provided me the regridded noaa map but with a warning (mentioned in the nb). Could you please help me in resolving the warning.
You already had this problem. The solution is thoroughly explained here: https://github.com/bopen/c3s-eqc-toolbox-template/issues/121#issuecomment-1829789739
In proba map, I want to estimate the area covered by pixels having albedo values > 0.1625 as well as < 0.1625. Similarly, in regridded noaa map, I want to estimate the area covered by pixels having albedo values > 0.215 as well as < 0.215. The purpose is to determine how much lower albedo values area have decreased/increased over the time-period.
I don't see the code to do that in the notebook. What's the issue? If you don't know how to mask, here is a couple of examples:
da_masked_1 = da.where(da > vmin)
da_masked_2 = da.where(da < vmax)
da_masked_3 = da.where((da > vmin) & (da < vmax))Hi Mattia,
Thank you for your response.
I am able to resolve the regridding issue and also calculate the area. But there seems to be a problem in the calculation of area probably due to grid size (I am not sure of the reason though).
As per my jupyter notebook, following are the results of area calculation. But, it cannot be correct as the sum of forest + non-forest area in 1990 is not equal to that of 2018, however, the same size of study area is chosen for year 1990 as well as 2018.
Could you please have a look at my jupyter notebooks and help me in resolving this miscalculation.
Forest area in 1990 (calculated using noaa_regrid) - 30 square units. Non-forest area in 1990 (calculated using noaa_regrid) - 24 square units.
Forest area in 2018 (calculated using proba) - 229 square units. Non-forest area in 2018 (calculated using proba) - 497 square units.
akankshabalha
commented
10 months ago malmans2
commented
10 months ago Hi, I don't fully follow. If you compare the number of pixels of variables with different resolution, the total number of pixels is going to be different. Also, in your code you were not interpolating the data because you were using an output grid identical to the input grid.
Try this snippet, I think it does what you need (maps + stats). You always refer to "area", but I don't see the code to actually compute area from lat lon. I added it, so both the number of pixels and area are displayed.
def plot_and_print_area_statistics(da, threshold):
resolution = dict(da.sizes)
print(f"{threshold=}")
print(f"{resolution=}")
grid_cell_area = diagnostics.grid_cell_area(da)
masks = {
"Non-Forest": da > threshold,
"Forest": da < threshold,
}
for label, mask in masks.items():
n_pixels = mask.sum()
total_area = grid_cell_area.where(mask).sum(keep_attrs=True)
if total_area.attrs["units"] == "m2":
with xr.set_options(keep_attrs=True):
total_area *= 1.0e-6
total_area.attrs["units"] = "km2"
print(f"\n{label}:")
print(f"\tNumber of pixels: {n_pixels.values}")
print(f"\tTotal area: {total_area.values:.2f}{total_area.attrs['units']}")
colors = [
"blue",
"skyblue",
"cyan",
"lightgreen",
"yellow",
"orange",
"red",
"maroon",
]
custom_cmap = mcolors.LinearSegmentedColormap.from_list("CustomColorMap", colors)
return plot.projected_map(da, cmap=custom_cmap, show_stats=False)
da_albedo_proba = ds_albedo_maps_proba["AL_BH_BB"]
_ = plot_and_print_area_statistics(da_albedo_proba, threshold=0.1625)
plt.title("proba")
plt.show()
da_albedo_noaa = ds_albedo_maps_noaa["AL_BH_BB"]
_ = plot_and_print_area_statistics(da_albedo_noaa, threshold=0.215)
plt.title("noaa")
plt.show()
da_albedo_noaa_regrid = diagnostics.regrid(
da_albedo_noaa, grid_out=da_albedo_proba, method="nearest_s2d"
)
_ = plot_and_print_area_statistics(da_albedo_noaa_regrid, threshold=0.215)
plt.title("noaa regridded")
plt.show()Is this issue solved?
akankshabalha
commented
9 months ago Yes, it is solved now. Thank you.
Notebook description
In the attached notebook, I have downloaded and transform the albedo satellite dataset for noaa_11 and proba. I want to elaborate the nb code to achieve following objectives:
1) The spatial resolution of both the datasets (noaa: 4km and proba: 1km) is different. Hence, I tried to regrid noaa dataset at similar resolution as that of proba. It provided me the regridded noaa map but with a warning (mentioned in the nb). Could you please help me in resolving the warning.
2) In proba map, I want to estimate the area covered by pixels having albedo values > 0.1625 as well as < 0.1625. Similarly, in regridded noaa map, I want to estimate the area covered by pixels having albedo values > 0.215 as well as < 0.215. The purpose is to determine how much lower albedo values area have decreased/increased over the time-period. Albedo_Forests.zip
Notebook link or upload
No response
Anything else we need to know?
No response
Environment