Fix problematic ancillary pixels #44

engelca · 2024-10-25T01:35:41Z

This pull request adds the ability to fix problematic ancillary pixels "on-the-fly" in the replace_landsurface scripts (when handling start dumps).

…lematic pixels.

… put new problematic_pixel code. Moved the mule replace operator to the common_utilities.py script because it is being repeated.

…ll three replace land/surface scripts.

…to the common_utilities.py script

…eadable.

CodeGat · 2024-10-30T05:09:22Z

Sorry @engelca, I'll be away for the next couple weeks and won't be able to review this, so I'll remove my request for review.

atteggiani

Hi @engelca,

I left a some comments and changes related to clarity and performance.

Feel free to ask in the specific comment for any further details about specific changes.

atteggiani · 2024-10-25T14:10:02Z

src/common_utilities.py

+        print('transform')
+        return sources[1]
+
+def replace_in_ff_problematic(f, mf_out, replace, stashcode, canopy_pixels, landsea_pixels):


It would be good to add a docstring to better explain the purpose of the function.

I would also rename to something a bit clearer like replace_problematic_pixels_in_ff

atteggiani · 2024-10-25T14:12:25Z

src/common_utilities.py

+
+    mf_out.fields.append(replace([f, data]))
+
+def problematic_pixels(infile):


I would rename the function to be a bit more specific of its purpose. For example:

get_problematic_pixels

find_problematic_pixels

atteggiani · 2024-10-25T14:18:16Z

src/common_utilities.py

+    current_data = f.get_data()
+    data=current_data.copy()
+
+    if stashcode == 218:


It would be better to capture this number in a constant variable at the beinning of the script (outside the function) that makes it clearer what it stands for.
Something like:

STASH_CANOPY_HEIGHT=218

then change

Suggested change

if stashcode == 218:

if stashcode == STASH_CANOPY_HEIGHT:

atteggiani · 2024-10-25T14:18:52Z

src/common_utilities.py

+            ix=canopy_pixels[j][1]
+            if np.isnan(data[iy,ix]):
+                data[iy,ix]=1.
+    elif stashcode == 33:


Same as above with the hardcoded number for the canopy height stash.

atteggiani · 2024-10-25T14:44:38Z

src/hres_ic.py

@@ -59,13 +59,13 @@ def main():

    # If necessary replace ERA5 land/surface fields with higher-resolution options
    if "era5land" in args.type:
-        replace_landsurface_with_ERA5land_IC.swap_land_era5land(args.mask, args.file, t)
+        replace_landsurface_with_ERA5land_IC.swap_land_era5land(args.mask, args.file, t, fix_problematic_pixels="yes")


Suggested change

replace_landsurface_with_ERA5land_IC.swap_land_era5land(args.mask, args.file, t, fix_problematic_pixels="yes")

replace_landsurface_with_ERA5land_IC.swap_land_era5land(args.mask, args.file, t, fix_problematic_pixels=True)

atteggiani · 2024-10-30T13:52:24Z

src/common_utilities.py

+    lats=orog_data['latitude'].data
+    lons=orog_data['longitude'].data


Suggested change

lats=orog_data['latitude'].data

lons=orog_data['longitude'].data

atteggiani · 2024-10-30T13:55:38Z

src/common_utilities.py

+    # Printing information to the standard output for reporting purposes (so scientists can be aware)
+    npoints,nxy = (canopy_pixels.shape)
+    if npoints>0:
+        for i in range(npoints):
+            print("%.1f, %.1f, Nan canopy"%(lons[canopy_pixels[i,1]],lats[canopy_pixels[i,0]]))
+
+    npoints,nxy = (landsea_pixels.shape)
+    if npoints>0:
+        for i in range(npoints):
+            print("%.1f, %.1f, Misplaced Orography"%(lons[landsea_pixels[i,1]],lats[landsea_pixels[i,0]]))


Suggested change

# Printing information to the standard output for reporting purposes (so scientists can be aware)

npoints,nxy = (canopy_pixels.shape)

if npoints>0:

for i in range(npoints):

print("%.1f, %.1f, Nan canopy"%(lons[canopy_pixels[i,1]],lats[canopy_pixels[i,0]]))

npoints,nxy = (landsea_pixels.shape)

if npoints>0:

for i in range(npoints):

print("%.1f, %.1f, Misplaced Orography"%(lons[landsea_pixels[i,1]],lats[landsea_pixels[i,0]]))

# Printing information to the standard output for reporting purposes (so scientists can be aware)

tmp = missing_canopy.where(missing_canopy.compute(), drop=True)

for lon, lat in zip(tmp.longitude, tmp.latitude):

print(f"{lat:.1f}, {lon:.1f}, NaN Canopy")

tmp = misclassified_land.where(misclassified_land.compute(), drop=True)

for lon, lat in zip(tmp.longitude, tmp.latitude):

print(f"{lat:.1f}, {lon:.1f}, Misplaced Orography")

Also a minor thing, I swapped the printing of 'latitude' and 'longitude' values (placing 'latitude' first), to be more consistent with the usual coordinate specification.

atteggiani · 2024-11-01T06:10:24Z

src/common_utilities.py

+        print('transform')
+        return sources[1]
+
+def replace_in_ff_problematic(f, mf_out, replace, stashcode, canopy_pixels, landsea_pixels):


Suggested change

def replace_in_ff_problematic(f, mf_out, replace, stashcode, canopy_pixels, landsea_pixels):

def replace_in_ff_problematic(f, mf_out, replace, missing_canopy, misclassified_land):

stashcode = f.lbuser4

I would avoid passing an argument that can be computed within the function.
stashcode is always going to be f.lbuser4, but since f is an argument of the function, we can compute it inside the function itself.
Note the function call should be modified in the other files after this modification.

For changes in canopy_pixels and landsea_pixels, refer to the related comments below.

atteggiani · 2024-11-01T06:32:20Z

src/common_utilities.py

+            print("%.1f, %.1f, Misplaced Orography"%(lons[landsea_pixels[i,1]],lats[landsea_pixels[i,0]]))
+
+    # returning pixel locations so they can be processed appropriately
+    return canopy_pixels,landsea_pixels


Suggested change

return canopy_pixels,landsea_pixels

return missing_canopy,misclassified_land

This would simplify the logic in the replace_in_ff_problematic function (see related comments above for details).

atteggiani · 2024-11-01T06:40:19Z

src/common_utilities.py

+    current_data = f.get_data()
+    data=current_data.copy()
+
+    if stashcode == 218:
+        for j in range(len(canopy_pixels)):
+            iy=canopy_pixels[j][0]
+            ix=canopy_pixels[j][1]
+            if np.isnan(data[iy,ix]):
+                data[iy,ix]=1.
+    elif stashcode == 33:
+        for j in range(len(landsea_pixels)):
+            data[landsea_pixels[j][0],landsea_pixels[j][1]]=0.


Suggested change

current_data = f.get_data()

data=current_data.copy()

if stashcode == 218:

for j in range(len(canopy_pixels)):

iy=canopy_pixels[j][0]

ix=canopy_pixels[j][1]

if np.isnan(data[iy,ix]):

data[iy,ix]=1.

elif stashcode == 33:

for j in range(len(landsea_pixels)):

data[landsea_pixels[j][0],landsea_pixels[j][1]]=0.

if stashcode == 218:

data = np.where(missing_canopy, 1., f.get_data())

elif stashcode == 33:

data = np.where(misclassified_land, 0., f.get_data())

If we use missing_canopy and misclassified_land (the masks) as input arguments we can simplify the fixing logic without needing loops.

engelca added 5 commits October 25, 2024 09:18

Adding in keyword to "swap" functions to allow for choice to fix prob…

5786851

…lematic pixels.

Creating a new common_utilities.py python script to create a place to…

5a9a278

… put new problematic_pixel code. Moved the mule replace operator to the common_utilities.py script because it is being repeated.

add "common_utilities" to the "problematic_pixels" function call in a…

3fcedd8

…ll three replace land/surface scripts.

add missing libraries and missing replace_in_ff_problematic function …

4841022

…to the common_utilities.py script

Move misplaced python definitions.

6ad51fc

engelca requested review from atteggiani and CodeGat October 25, 2024 01:35

engelca linked an issue Oct 25, 2024 that may be closed by this pull request

Fix problematic ancillary pixels while doing replace land surface in start dump. #43

Open

engelca added 2 commits October 25, 2024 14:07

Remove unnecessary print statement so problematic pixel output more r…

46e4017

…eadable.

Add copyright statement to new common_utilities.py script.

bd98cca

CodeGat removed their request for review October 30, 2024 05:09

atteggiani requested changes Nov 1, 2024

View reviewed changes

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

Fix problematic ancillary pixels #44

Fix problematic ancillary pixels #44

engelca commented Oct 25, 2024

CodeGat commented Oct 30, 2024

atteggiani left a comment

atteggiani Oct 25, 2024

atteggiani Oct 25, 2024

atteggiani Oct 25, 2024

atteggiani Oct 25, 2024

atteggiani Oct 25, 2024

atteggiani Oct 30, 2024

atteggiani Oct 30, 2024

atteggiani Nov 1, 2024

atteggiani Nov 1, 2024

atteggiani Nov 1, 2024


		mf_out.fields.append(replace([f, data]))

		def problematic_pixels(infile):

	replace_landsurface_with_ERA5land_IC.swap_land_era5land(args.mask, args.file, t, fix_problematic_pixels="yes")
	replace_landsurface_with_ERA5land_IC.swap_land_era5land(args.mask, args.file, t, fix_problematic_pixels=True)

		lats=orog_data['latitude'].data
		lons=orog_data['longitude'].data

	def replace_in_ff_problematic(f, mf_out, replace, stashcode, canopy_pixels, landsea_pixels):
	def replace_in_ff_problematic(f, mf_out, replace, missing_canopy, misclassified_land):
	stashcode = f.lbuser4

	return canopy_pixels,landsea_pixels
	return missing_canopy,misclassified_land

Fix problematic ancillary pixels #44

Are you sure you want to change the base?

Fix problematic ancillary pixels #44

Conversation

engelca commented Oct 25, 2024

CodeGat commented Oct 30, 2024

atteggiani left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment