More robust mksurfdata_map logic for determining where to put wetlands

[billsacks]

As noted in #545 , the logic determining where to put wetlands is tied to the mask on the pft raw data file:

https://github.com/ESCOMP/ctsm/blob/a34404419aaa81f5021602963222c50daa10c0f6/tools/mksurfdata_map/src/mksurfdat.F90#L707-L724

I think a more robust way to do this would be: rather than just determining pctland_pft, as is currently done, instead determine a pctland field as the max of the pctland from each pct raw dataset. So, we'd start by setting pctland equal to 0. Then each routine that remaps a pct field (pct pft/cft, pct glacier, pct lake and pct urban) would accept pctland as an inout argument. Each of these routines would compute a my_pctland field, similarly to how pctland_pft is currently computed, and then have code like:

do i = 1, n
   pctland(i) = max(pctland(i), my_pctland(i))
end do

Then we could still turn areas to wetland that are outside the landmask - but we'd just do this for areas outside the union of all PCT fields' landmasks, rather than using only the pctpft landmask.

@dlawrenncar @lawrencepj1 @ekluzek @mvertens does this seem reasonable to you? I don't think we should do this for the cmip6 runs (it could change answers in more than just the Antarctic ice shelves), but I'm considering this as a long-term solution.

[billsacks]

In addition to making the behavior over ice shelves more robust to new pft rawdata files, another area this will make more robust is large water bodies like the Caspian Sea, which we would like to treat as lake rather than wetland.

[billsacks]

Actually, I think we want something a little different from what I first said: Rather than using the landmasks on the various pct special landunit datasets, we should use the actual remapped area - so we change an area to wetland only if it is outside the veg mask and has 0% cover of all special landunits. (I'm not sure whether we want to use the landmask on the pft dataset or use the % cover of the natural veg + crop landunits; that might amount to the same thing.) An example of where this makes a difference: Imagine there is a grid cell that is outside the landmask of most of the raw datasets; it is inside the landmask of the lake raw dataset, but the lake raw dataset says it's 0% lake. In my original algorithm, this would not be converted to wetland, but then its cover type would be indeterminate; with this revision, this point would be converted to wetland.

[dlawrenncar]

I agree that this would be a good long term solution, but that it could potentially change something outside of Antarctica and therefore would be risky to implement now.

…

On Tue, Oct 30, 2018 at 11:45 AM Bill Sacks ***@***.***> wrote: Actually, I think we want something a little different from what I first said: Rather than using the landmasks on the various pct special landunit datasets, we should use the actual remapped area - so we change an area to wetland only if it is outside the veg mask and has 0% cover of all special landunits. (I'm not sure whether we want to use the landmask on the pft dataset or use the % cover of the natural veg + crop landunits; that might amount to the same thing.) An example of where this makes a difference: Imagine there is a grid cell that is outside the landmask of most of the raw datasets; it is inside the landmask of the lake raw dataset, but the lake raw dataset says it's 0% lake. In my original algorithm, this would *not* be converted to wetland, but then its cover type would be indeterminate; with this revision, this point *would* be converted to wetland. — You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub <#553 (comment)>, or mute the thread <https://github.com/notifications/unsubscribe-auth/AUAcVIMwregBwiSWzmy6abbkWOE8w_uxks5uqJArgaJpZM4YCcaa> .

[billsacks]

In discussion with @lawrencepj1 - he's going to try making changes to mksurfdata_map to accomplish the goals of #545 . (This seems significantly easier than trying to change the datasets, and moves us towards the long-term solution.) To first order, this could be as simple as changing this:

https://github.com/ESCOMP/ctsm/blob/527f3262e550c96b28ce78ee307457bf6adb2e46/tools/mksurfdata_map/src/mksurfdat.F90#L711

to something like:

if (pctlnd_pft(n) < 1.e-6_r8 .and. pctgla(n) < 1.e-6_r8) then

However, some thought needs to be given to how this interacts with the error-checks and normalizations that are done later in mksurfdata_map. Two that jump out at me are:

(1) The truncation of percentage fields on the output grid: This currently is done right after the above block of code:

https://github.com/ESCOMP/ctsm/blob/527f3262e550c96b28ce78ee307457bf6adb2e46/tools/mksurfdata_map/src/mksurfdat.F90#L726-L736

If you just do what I suggest above, then I think you can have a situation where a grid cell has very small glacier area, so doesn't trigger the set-to-wetland code, but then the glacier is set to 0, at which point you might have 0% of all special landunits - which I think would result in the grid cell being set to vegetated.

I think the most robust way to fix this would be to move the pctgla truncation (pctgla(n) = float(nint(pctgla(n)))) to above the if (pctlnd_pft(n) < 1.e-6_r8) then block. Then you can change the conditional to:

if (pctlnd_pft(n) < 1.e-6_r8 .and. pctgla(n) <= 0._r8) then

(2) This block in normalizencheck_landuse:

https://github.com/ESCOMP/ctsm/blob/527f3262e550c96b28ce78ee307457bf6adb2e46/tools/mksurfdata_map/src/mksurfdat.F90#L1561-L1568

Consider a grid cell that is outside pftdata_mask and has (say) 50% glacier and no other special landunits. I believe that the code at the start of this subroutine would set the area to 50% glacier, 50% veg, which would cause this error check to fail. This needs to be given some thought, but my initial reaction is that we may want to remove this error check, but make sure that pct_nat_pft and pct_cft are set sensibly even when pftdata_mask is 0 (e.g., being set to 100% bare ground for pct_nat_pft and 100% generic c3 grop for pct_cft); I'm not sure whether that's currently the case.

Edit 2018-11-07: @lawrencepj1 is going to take a more conservative approach for now: see #545 (comment)

[ekluzek]

Fix in release-clm5.0.16