Air temperature incomparable across MOST and RSL simulations

General
2

Hi @Harro, welcome to the forum, and thanks for such a precise diagnosis. You have found the right culprit, and it is a limitation of how SUEWS defines its “2 m” diagnostics rather than anything in your setup.

The mechanism is written up in my reply on the tree-density thread, so I will not repeat it here. The short form: under RSL the diagnostics are read at 2 m and 10 m above ground, inside the canopy; under MOST they are read at 2 + z_d + z_{0m} and 10 + z_d + z_{0m}, which is the standard MOST convention and puts “2 m” near canopy top. The two schemes therefore report different quantities, and a city map that mixes them is not comparable by construction. The relevant lines are in suews_phys_rslprof.f95 at 2026.6.5, L618 to L628, identical in 2026.4.3.

Two details specific to your case:

  • The automatic switch uses only 0.1 < \mathrm{PAI} < 0.68 and z_H > 2 m; the FAI criterion mentioned in the docs is not active in the code. That is why your steps sit exactly at those thresholds.
  • The switch itself is not used in the energy balance, but with roughness_sublayer_level: basic the diagnosed T_2 feeds back into anthropogenic heat and surface conductance, which is the small flux step you see. Set roughness_sublayer_level: none if you want the fluxes strictly scheme-independent.

Options for a comparable city-wide field. None is perfect, so I would rather lay them out than pick for you:

  1. Force RSL everywhere (roughness_sublayer: rst). Consistent 2 m above ground, but the Harman and Finnigan profile is then applied outside the PAI range it was derived for, and for very low canopies the code clamps z_H to 2 m, so the low-PAI grids are extrapolations.
  2. Force MOST everywhere (roughness_sublayer: most). One scheme, but the sampling height still varies grid to grid with z_d + z_{0m}.
  3. Write out the profile and interpolate yourself. Add RSL to model.control.output.groups (txt output; in the Python API the RSL variables are always in the results DataFrame). The group carries the 30-level z, U, T and q profiles plus zd_RSL, z0_RSL, zH_RSL, UStar_heat and TStar_RSL, so you can read every grid at one common height and report which height that is. One constraint: under MOST the profile starts just above z_d + z_{0m}, so a common height above ground only exists if it clears the largest z_d + z_{0m} in your domain, or if you combine this with option 1. A common height above z_d + z_{0m} works if a canopy-top quantity serves your purpose.

For a city map, option 3 is the most defensible, and whichever you choose I would carry the surface energy fluxes and surface temperature alongside, since they do not have this reference-height ambiguity.

On our side, the docs should state plainly that the MOST diagnostics are referenced to z_d + z_{0m}, and the output should make the actual sampling height visible. How best to do that is still open, and this thread is a good place to work it out.