Skip to content

Conversation

@keltonhalbert
Copy link

Work in progress for #1273, using the modifications provided by @EdwardColon-NOAA via email.

I currently do not have access to WCOSS2, so I am trying to install dependencies locally in order to confirm that it compiles without error. In the mean-time, I wanted to go ahead and provide the starting point, and ask a few clarifying questions (@BenjaminBlake-NOAA) :

  1. In MISCLN.f, I just want to double-check and make sure that the T1D and TBND arrays are potential temperature, and not air temperature. It seems like they should be, based on the code, but I could not find a definitive answer and wanted to make absolutely sure.

  2. Right now, the mixed-layer and surface-based parcel calculations are overwriting the lowest model grid level values. Is that something that we want to do? It may not ultimately be of significant consequence, but there is another way this can be approached in which the parcel starting P/T/Q is used to compute the LCL, and then the parcel value below the LCL has a constant virtual potential temperature. It would effectively preserve the values in the 1D arrays from the model. This would likely require significant modifications, though, since it looks like the LCL is an integer index rather than a computed value. My gut says "probably not", but I wanted to pose the question.

  3. Do these modifications need to be made to CALCAPE2? I'm not well versed in the UPP code and am not exactly sure what this does differently, and for what purpose. I do see some stuff about effective inflow layer calculations 0-3 km CAPE, which could be sensitive to this change as well.

Lastly, I noticed that in the MLCAPE code the pressure values are averaged.

UPP/sorc/ncep_post.fd/MISCLN.f

Lines 3092 to 3106 in 8f6caa9

DO J=JSTA,JEND
DO I=ISTA,IEND
EGRID1(I,J) = -H99999
EGRID2(I,J) = -H99999
LB2(I,J) = (LVLBND(I,J,1) + LVLBND(I,J,2) + &
LVLBND(I,J,3))/3
P1D(I,J) = (PBND(I,J,1) + PBND(I,J,2) + PBND(I,J,3))/3
T1D(I,J) = (TBND(I,J,1) + TBND(I,J,2) + TBND(I,J,3))/3
Q1D(I,J) = (QBND(I,J,1) + QBND(I,J,2) + QBND(I,J,3))/3
ENDDO
ENDDO
!
DPBND = 0.
CALL CALCAPE(ITYPE,DPBND,P1D,T1D,Q1D,LB2,EGRID1, &
EGRID2,EGRID3,EGRID4,EGRID5)

I could be misunderstanding what the code is doing (potentially related to my question about the T1D array), but canonically at SPC, the mixed-layer parcel is lifted from the shelter pressure value. The 100mb mean layer potential temperature is used as a proxy for boundary-layer mixing processes, but is still lifted from the surface. I wrote the code to maintain the existing functionality, so if we want to keep it that way we can... but also figured I'd mention the difference just in case. I'm also not sure what LB2 is and whether this average needs modification, too.

@keltonhalbert keltonhalbert marked this pull request as draft July 22, 2025 17:14
@BenjaminBlake-NOAA
Copy link
Collaborator

@keltonhalbert Thanks for working on this! I'm tagging @jaymes-kenyon here as he may be more familiar with the UPP CAPE routines than I am, but here are my two cents:

  1. It looks like the T1D and TBND arrays are actually for air temperature and not potential temperature. T1D and TBND are used as inputs to the CALPOT routine, which calculates potential temperature. CALPOT is called within MISCLN.f.

  2. I am also thinking no, but it would be good to get confirmation on this.

  3. I believe the same modifications should also be made to CALCAPE2 within UPP_PHYSICS.f.

  4. I'm not sure what approach regarding the MLCAPE averaging is best, but perhaps we should stick with the current approach for now.

@keltonhalbert
Copy link
Author

1. It looks like the T1D and TBND arrays are actually for air temperature and not potential temperature.  T1D and TBND are used as inputs to the CALPOT routine, which calculates potential temperature.  CALPOT is called within MISCLN.f.

Ahhhh -- well that is definitely a difference worth noting. The averaging process being in T/P space instead of THETA space could result in some differences, particularly in models with greater grid spacing in the lowest levels. We can keep it in T/P space, but it would probably be more appropriate to average in potential temperature space.

I can make the modifications to CALCAPE2 then!

@EdwardColon-NOAA
Copy link
Contributor

EdwardColon-NOAA commented Jul 22, 2025 via email

@WenMeng-NOAA WenMeng-NOAA added the enhancement New feature or request label Jul 23, 2025
@keltonhalbert
Copy link
Author

(1) I am in agreement with Ben regarding the T1D and TBND array units which would complicate things a bit in the SCAPE and MUCAPE calculation.
(2) The direct substitution of T2M and Q2M is erroneous so, no, it needs to be changed.
(3) CALCAPE2 would need to reflect the changes in CALCAPE, and (4) the MLCAPE averaging is simplistic and may require some weighting of components but is the more important of the tasks.

I have removed the calculations of potential temperature derived from the shelter fields, and have instead just done a direct read from the shelter arrays into the appropriate variable. Something that was missing was also the setting of the shelter pressure... the short variable names in FORTRAN make it hard for me to know for sure, but I think I set it to the correct variable PKL based on looking at the surrounding code. If someone could help confirm, that would be appreciated. I also copied there adjustments into CALCAPE2.

A few additional questions and observations as I continue to familiarize myself with the code...

  1. The shelter fields are currently only set when ITYPE = 1 is set. I saw a comment that says for this type, P1D, T1D, and Q1D are "dummy arrays"... not entirely sure what this means here? The documentation says that this is for finding the maximum THETA-E layer, and it reads from P1D, T1D, etc. Either way, anything with this ITYPE should now include the contributions from the 2-meter fields. I think/hope...
  2. When ITYPE=2, manual "intervention" is required. This is currently set up for mixed layer calculations, but there are others listed that may need to have this change included for consistency. See below, with line numbers corresponding to MISCLN.f.
  • ITYPE = 1: (these are covered by changes to CALCAPE[2])

    • Most-unstable CAPE/CINH | L 3242
    • Significant Tornado Parameter (fixed-layer) | L4172
  • ITYPE = 2: (These have to be covered by changes external to CALCAPE[2])

    • Best boundary layer CAPE/CINH | L2077
    • Mixed-layer CAPE/CINH | L3107
    • Effective Inflow layer | L3477
    • 0-3 km CAPE/CINH | L3602
    • mixed-layer LFC (part of STP calc -- explicitely sets virtual T??) | L4297

Something else I noticed is that in the call to CALCAPE for the virtual mixed-layer LFC (L4297), it explicitly makes the call to the virtual temperature correction for the mixed-layer, where the other routines do not. However, within CALCAPE, I see calls to TVIRTUAL as well. My assumption was that CALCAPE is already doing the virtual temperature correction -- is this not the case? I'm a bit confused and could definitely use some clarification on what's happening there.

I went ahead and updated the mixed-layer LFC calculation to include the 2-meter fields with the virtual temperature correction. However, it would be good to have some clarification on why the other calculations do not have or need these corrections. Mixed-layer, most-unstable, and surface-based parcel temperature definitions all need to use the virtual temperature consistently to be correct/accurate.

Are you compiling UPP on Jet?

Unfortunately I do not have access to Jet, and at one point I was given a WCOSS user account. However, that was earlier in the spring, and I never could get Tectia installed on my workstation due to staffing issues. Now, I am waiting for SPC to be given a license for the new software, SecureCRT, in order to access. So, in the mean time, I am using Spack to set up a local dev environment for UPP that I can at least test compilation on -- though I may have to employ assistance to test on actual output until I get on one of these systems.

@EdwardColon-NOAA
Copy link
Contributor

EdwardColon-NOAA commented Jul 24, 2025 via email

@keltonhalbert
Copy link
Author

Hi @EdwardColon-NOAA,

Thank you for some of the clarifications, and your willingness to be an intermediate compiler/tester!

Regarding your question about the 2-meter virtual temperature fields, this is the static link to the lines in MISCLN.f:

UPP/sorc/ncep_post.fd/MISCLN.f

Lines 4285 to 4292 in 0ac3cfa

P1D(I,J) = (PBND(I,J,1) + PBND(I,J,2) + PBND(I,J,3) + &
PSHLTR(I, J))/4
T1D(I,J) = (TVIRTUAL(TBND(I,J,1),QBND(I,J,1)) + &
TVIRTUAL(TBND(I,J,2),QBND(I,J,2)) + &
TVIRTUAL(TBND(I,J,3),QBND(I,J,3)) + &
TVIRTUAL(TSHLTR(I,J), QSHLTR(I,J))/4
Q1D(I,J) = (QBND(I,J,1) + QBND(I,J,2) + QBND(I,J,3) + &
QSHLTR(I,J))/4

This is the section that is potentially "double correcting" for the virtual temperature when computing the mixed-layer LFC. It appears to then call CALCAPE2 with this configuration.

@keltonhalbert
Copy link
Author

If we are in agreement on this, I have removed the calls to TVIRTUAL outside of CALCAPE[2]. It appears this was the only instance of that function being called in MISCLN.f, so I removed it from the import at the top of the file.

The following calls to CALCAPE[2] have been updated to include the 2-meter values:

  • Mixed-layer CAPE/CINH | L3107
  • 0-3 km CAPE/CINH | L3602
  • mixed-layer LFC (part of STP calc -- no longer sets virtual T) | L4297

The two remaining calculations that need addressing:

  • Best boundary layer CAPE/CINH | L2077
  • Effective Inflow layer | L3477

Best Boundary Layer CAPE

https://github.com/keltonhalbert/UPP/blob/581d81acdbeac1e08b1fa1d33e85663476823944/sorc/ncep_post.fd/MISCLN.f#L2059-L2078

Effective Inflow Layer

https://github.com/keltonhalbert/UPP/blob/581d81acdbeac1e08b1fa1d33e85663476823944/sorc/ncep_post.fd/MISCLN.f#L3459-L3478

Both of these are loops over a vertical coordinate, with the best boundary layer being over NBND and the effective inflow layer being over LM. Neither of these pose easy ability to add the shelter fields to the loop indices, but there are probably a few options. I'm having a hard time deciphering the "best boundary layer CAPE", but I am assuming it is taking the maximum value... so there should be a way to compare the shelter field CAPE with the model levels and take the maximum. The effective inflow layer would be a bit more tricky, since it is also using the equilibrium level value based on the effective inflow criteria.

@keltonhalbert
Copy link
Author

After many trials and tribulations, I finally managed to get the UPP and its dependencies compiled and installed. Turned out, the version of CRTM spack defaults to does not work, but setting it to the version in the spack.yaml file for CI worked.

I can confirm that my changes at least compile. I can't speak to the quality of the output, but it at least builds.

@EdwardColon-NOAA
Copy link
Contributor

EdwardColon-NOAA commented Jul 28, 2025 via email

@EdwardColon-NOAA
Copy link
Contributor

EdwardColon-NOAA commented Jul 29, 2025 via email

@keltonhalbert
Copy link
Author

@EdwardColon-NOAA,

Since the changes are relatively lightweight, I can apply them to a fork of release/3drtma if need be. The last thing that needs to be figured out is how to include the shelter fields in the vertical loops for the effective inflow layer and best boundary layer CAPE. I have been working forecast shifts this week so I haven't come up with a solution yet, but if and when that is the case, can work to apply it to the other branch as well.

@EdwardColon-NOAA
Copy link
Contributor

EdwardColon-NOAA commented Jul 30, 2025 via email

@BenjaminBlake-NOAA
Copy link
Collaborator

@keltonhalbert Work is in progress on updating the 3DRTMA UPP. After this work is complete we will revisit your CAPE/CIN PR. Thanks for your patience!

@BenjaminBlake-NOAA BenjaminBlake-NOAA changed the title WIP: CAPE/CINH using 2-m fields CAPE/CINH using 2-m fields: develop branch Jan 27, 2026
@BenjaminBlake-NOAA BenjaminBlake-NOAA moved this from On Hold to Draft in PRs to Process Jan 28, 2026
P1D(I,J) = (PBND(I,J,1) + PBND(I,J,2) + PBND(I,J,3) + &
PSHLTR(I,J))/4
T1D(I,J) = (TBND(I,J,1) + TBND(I,J,2) + TBND(I,J,3) + &
TSHLTR(I,J))/4
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Thanks @GangZhao-NOAA for catching a glitch: TSHLTR inside UPP is potential temperature. I suggest changing it into 2-m T as follows:

TSHLTR(I,J) --> TSHLTR(I,J)*(PSHLTR(I,J)/P1000)**CAPA

Tagging @GangZhao-NOAA @BenjaminBlake-NOAA @keltonhalbert

Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@GangZhao-NOAA Thanks for catching this - @WenMeng-NOAA I'll get this updated tomorrow. At some point we'll also want to get this updated in the 3DRTMA release branch, perhaps combined with a future PR that Gang will submit.

@WenMeng-NOAA
Copy link
Collaborator

@BenjaminBlake-NOAA There are a syntax errors as below, which might cause CI failure:
'''
Error: Symbol ‘p1000’ at (1) has no IMPLICIT type; did you mean ‘p00k’?
/home/runner/work/UPP/UPP/UPP/sorc/ncep_post.fd/UPP_PHYSICS.f:701:59:
Error: Symbol ‘p1000’ at (1) has no IMPLICIT type; did you mean ‘p00k’?
make[2]: *** [sorc/ncep_post.fd/CMakeFiles/upp.dir/build.make:1625: sorc/ncep_post.fd/CMakeFiles/upp.dir/UPP_PHYSICS.f.o] Error 1
'''

@BenjaminBlake-NOAA
Copy link
Collaborator

@WenMeng-NOAA Thanks, looks like I need to declare P1000 - I'll get this fixed.

@WenMeng-NOAA
Copy link
Collaborator

The baseline of gfs test has the following changes:

892:2486066395:CAPE:surface:rpn_corr=0.973987:rpn_rms=245.415
893:2488576822:CIN:surface:rpn_corr=0.744611:rpn_rms=37.2207
1044:3071206478:CAPE:90-0 mb above ground:rpn_corr=0.987647:rpn_rms=94.6981
1045:3073362134:CIN:90-0 mb above ground:rpn_corr=0.741131:rpn_rms=33.1882
1046:3074565713:CAPE:255-0 mb above ground:rpn_corr=0.972574:rpn_rms=245.248
1047:3076869436:CIN:255-0 mb above ground:rpn_corr=0.687106:rpn_rms=35.5124
1048:3077818705:PLPL:255-0 mb above ground:rpn_corr=0.969716:rpn_rms=4173

@WenMeng-NOAA WenMeng-NOAA linked an issue Feb 4, 2026 that may be closed by this pull request
@BenjaminBlake-NOAA BenjaminBlake-NOAA marked this pull request as ready for review February 4, 2026 16:57
@BenjaminBlake-NOAA BenjaminBlake-NOAA added Baseline Change The baselines of the UPP regression tests are changed. Ready for Review This PR is ready for code review. labels Feb 4, 2026
@BenjaminBlake-NOAA
Copy link
Collaborator

@keltonhalbert After discussing with @WenMeng-NOAA , we are going to move forward with merging the work in this PR. For most affected applications (except 3DRTMA) we will leave capecin_2m=.false. for now until further testing/evaluation has been completed. I know there is a desire to add these changes to GFS and RRFS. Thank you for your significant contribution to UPP!

@WenMeng-NOAA With the current settings, there are only baseline changes to the 3DRTMA test (capecin_2m=.true.). There are baseline changes for RRFS but they can be ignored because I know you already updated the baselines for my RRFS PR, which we can process first. GEFSv12 and RRFS_IFI_missing were the only RTs that did not have baseline changes when capecin_2m=.true.

@BenjaminBlake-NOAA
Copy link
Collaborator

The CAPE/CIN modifications in this PR result in changes to the following UPP RTs when the capecin_2m namelist option is enabled (all tests except gefsv12 and rrfs_ifi_missing):
3drtma, aqm, dafs, gcafs, gefsv13, gfs, hafs, hrrr, hrrr_ifi, mpas, mpas_hfip, nmmb, rap, rrfs, rrfs_ifi, sfs

The variables that are changed with each test are outlined in this Google Doc (accessible to anyone with a NOAA account)

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

Baseline Change The baselines of the UPP regression tests are changed. enhancement New feature or request Ready for Review This PR is ready for code review.

Projects

Status: Draft

Development

Successfully merging this pull request may close these issues.

Updating UPP CAPE/CIN calculations to include 2-m fields

4 participants