Scheduled Downtime
On Friday 21 April 2023 @ 5pm MT, this website will be down for maintenance and expected to return online the morning of 24 April 2023 at the latest

Found this burried deep in the forum. Did WRF-ARW ever get the gust calculation ?

William.Hatheway

William.Hatheway

Active member
I found some very old code on the WRF User's Forum: https://forum.wrfforum.com/viewtopic.php?f=8&t=948
from Robert Rosumalski.

SUBROUTINE CALGUST(LPBL,ZPBL,GUST)
C$$$ SUBPROGRAM DOCUMENTATION BLOCK
C . . .
C SUBPROGRAM: CALGUST COMPUTE MAX WIND LEVEL
C PRGRMMR: MANIKIN ORG: W/NP2 DATE: 97-03-04
C
C ABSTRACT:
C THIS ROUTINE COMPUTES SURFACE WIND GUST BY MIXING
C DOWN MOMENTUM FROM THE LEVEL AT THE HEIGHT OF THE PBL
C
C
C PROGRAM HISTORY LOG:
C 03-10-15 GEOFF MANIKIN
C 05-03-09 H CHUANG - WRF VERSION
C 05-06-30 R ROZUMALSKI - DYNAMIC MEMORY ALLOCATION AND SMP
C THREAD-SAFE VERSION
C
C USAGE: CALL CALGUST(GUST)
C INPUT ARGUMENT LIST:
C NONE
C
C OUTPUT ARGUMENT LIST:
C GUST - SPEED OF THE MAXIMUM SFC WIND GUST
C
C OUTPUT FILES:
C NONE
C
C SUBPROGRAMS CALLED:
C UTILITIES:
C H2V
C
C LIBRARY:
C COMMON -
C LOOPS
C OPTIONS
C MASKS
C INDX
C
C ATTRIBUTES:
C LANGUAGE: FORTRAN 90
C MACHINE : CRAY C-90
C$$$
C
C
use vrbls3d
use vrbls2d
C
C INCLUDE ETA GRID DIMENSIONS. SET/DERIVE PARAMETERS.
C
!
INCLUDE "params"
C
INCLUDE "CTLBLK.comm"
C
C DECLARE VARIABLES.
C
INTEGER :: LPBL(IM,JM)
REAL :: GUST(IM,JM)
REAL ZPBL(IM,jsta_2l:jend_2u)
C
C
C*****************************************************************************
C START CALMXW HERE.
C
C LOOP OVER THE GRID.
C
DO J=JSTA,JEND
DO I=1,IM
! GUST(I,J) = SPVAL
GUST(I,J) = 0.
ENDDO
ENDDO
C
C ASSUME THAT U AND V HAVE UPDATED HALOS
C
!$omp parallel do
!$omp& private(ie,iw,mxww,u0,v0,wind)
DO 20 J=JSTA_M,JEND_M
DO 20 I=2,IM-1
L=LPBL(I,J)
IF(MODELNAME .EQ. 'NMM')THEN
IE=I+MOD(J+1,2)
IW=I+MOD(J+1,2)-1

USFC=D25*(U10(I,J-1)+U10(IW,J)+
X U10(IE,J)+U10(I,J+1))
VSFC=D25*(V10(I,J-1)+V10(IW,J)+
X V10(IE,J)+V10(I,J+1))
SFCWIND=SQRT(USFC**2 + VSFC**2)
U0 = D25*(U(I,J-1,L)+U(IW,J,L)+
X U(IE,J,L)+U(I,J+1,L))
V0 = D25*(V(I,J-1,L)+V(IW,J,L)+
X V(IE,J,L)+V(I,J+1,L))
WIND=SQRT(U0**2 + V0**2)

ELSE IF(MODELNAME .EQ. 'NCAR')THEN
USFC=U10(I,J)
VSFC=V10(I,J)
SFCWIND=SQRT(USFC**2 + VSFC**2)
U0=U(I,J,L)
V0=V(I,J,L)
WIND=SQRT(U0**2 + V0**2)
END IF
DELWIND=WIND - SFCWIND
ZSFC=FIS(I,J)*GI
DELWIND=DELWIND*(1.0-AMIN1(0.5,ZPBL(I,J)/2000.))
GUST(I,J)=SFCWIND+DELWIND
10 CONTINUE
20 CONTINUE

C END OF ROUTINE.
C
RETURN
END
 
forum.mmm.ucar.edu

Gust_Wind Implementation?

I am upgrading older wrf models to wrf 3.7. One goal was to add windgust to the output parameters. It appears clwrf_gust_wind was added in wrf 3.3.1, but the information available on how to accomplish that seems very limited (from forum and google searches). It appears additional parameters...
forum.mmm.ucar.edu forum.mmm.ucar.edu
 
@dudhia

Code: 
wrfout = Dataset('path/to/wrfout.nc')

frame = wrf.ALL_TIMES

wspd_wdir10 = wrf.getvar(wrfout, 'wspd_wdir10', timeidx=frame)
wspd_wdir = wrf.getvar(wrfout, 'wspd_wdir', timeidx=frame)
pblh = wrf.getvar(wrfout, 'PBLH', timeidx=frame)
height = wrf.getvar(wrfout, 'height_agl', timeidx=frame)
u10 = wrf.getvar(wrfout, 'U10', timeidx=frame)
v10 = wrf.getvar(wrfout, 'V10', timeidx=frame)
wrfout.close()

wspd_pblh = wrf.interplevel(wspd_wdir, height, pblh)                                               # 1.
g_pblh = np.where(pblh <= 1000, pblh, 1000)                                                          # 2.
G = (wspd_wdir10[0] + (wspd_pblh[0] - wspd_wdir10[0])*(1-(g_pblh/2000)))*3.6    # 3. gust in km/h

I have attached the complete script to this email.

_________________________________________________

Download

Matías

1709036270986.png
From the WRF-Python Forum
Credit to Matias Suarez
 
I think different definitions of a windgust exist depending on it's use. I would say a sudden brief increase in windspeed at a certain pressure level. It would be nice to be able to determine the maximum windspeed over all time steps within a certain time period/interval. Have you seen any code for this?
 
Michel said:
I think different definitions of a windgust exist depending on it's use. I would say a sudden brief increase in windspeed at a certain pressure level. It would be nice to be able to determine the maximum windspeed over all time steps within a certain time period/interval. Have you seen any code for this?
Click to expand...
@Michel

Personally no, but it doesn't mean it isn't there.
If we can get some people to test this python code and then show that it works, we can hopefully get it implemented into WRF as a registry file.

Then once it is in there doing what you want should be pretty easy I think
 
William.Hatheway said:
View attachment 14961

I have tested this and the values seem to match up with observations.
Click to expand...
So, this Python code is inspired by two main approaches: Brasseur’s method and the NOAA RUC20 model. Both of these are commonly used for estimating wind gusts near the ground, but they each bring something unique to the table.

Starting with Surface Wind:
First off, we grab the wind speed at 10 meters above the ground (wspd_wdir10). This is our base wind speed at the surface, which is the foundation for predicting gusts. The RUC20 model from NOAA actually does this too, using the 10-meter wind as a starting point since it represents what’s happening right at ground level.

Looking at Wind Speeds Higher Up:
Then, look at wind speeds at different heights within the boundary layer. Brasseur’s method is all about this – it assumes that gusts can be influenced by winds from any level within the PBL (planetary boundary layer) as long as there’s enough turbulence to bring those higher winds down to the surface. So, we grab the winds from each level, which gives us a range of speeds that could potentially contribute to gusts near the ground.

Applying a Height-Based Weight:
Now, this is where the NOAA RUC20 model is utilized. We apply a weight that decreases with height. Up to 1 km, the weight gradually drops from 1.0 to 0.5, and for anything higher, it stays at 0.5. The idea is that winds higher up have less influence on surface gusts, but they’re still somewhat relevant. This weight factor is straight out of the RUC20 model, which scales down the effect of higher-level winds on surface gusts.

Adjusting and Finding the Max Gust:
For each level, we adjust the surface wind by adding the weighted difference between that level’s wind speed and the 10-meter wind. This gives us an adjusted gust speed for each level. Then, we just keep the maximum gust from all the levels, which, again, lines up with Brasseur’s method – he suggests looking at the strongest possible wind from within the PBL as the potential gust.

Why This Works
This code basically blends NOAA RUC20’s height-based scaling with Brasseur’s idea of pulling from all levels within the PBL. By combining these two methods, we’re able to get a solid estimate of what the maximum gust at the surface might look like, especially during turbulent conditions.


# Calculate the wind speed at each level and find the maximum gust speed
max_gust_speed = np.zeros_like(wspd_wdir10) # Assuming wspd_wdir10 is the surface wind speed


# Iterate over each vertical level
for k in range(u.shape[0]):
# Calculate the wind speed at level k

wind_speed_k = np.sqrt(u[k]**2 + v[k]**2)
height_k = height[k]

# Calculate the weight based on height

weight_k = np.where(height_k <= 1000, 1 - height_k / 2000, 0.5)

# Adjust the wind speed by the weighted difference from the surface

adjusted_speed_k = wspd_wdir10 + (wind_speed_k - wspd_wdir10) * weight_k

# Update the maximum gust speed
max_gust_speed = np.maximum(max_gust_speed, adjusted_speed_k)
 
William.Hatheway said:
So, this Python code is inspired by two main approaches: Brasseur’s method and the NOAA RUC20 model. Both of these are commonly used for estimating wind gusts near the ground, but they each bring something unique to the table.

Starting with Surface Wind:
First off, we grab the wind speed at 10 meters above the ground (wspd_wdir10). This is our base wind speed at the surface, which is the foundation for predicting gusts. The RUC20 model from NOAA actually does this too, using the 10-meter wind as a starting point since it represents what’s happening right at ground level.

Looking at Wind Speeds Higher Up:
Then, look at wind speeds at different heights within the boundary layer. Brasseur’s method is all about this – it assumes that gusts can be influenced by winds from any level within the PBL (planetary boundary layer) as long as there’s enough turbulence to bring those higher winds down to the surface. So, we grab the winds from each level, which gives us a range of speeds that could potentially contribute to gusts near the ground.

Applying a Height-Based Weight:
Now, this is where the NOAA RUC20 model is utilized. We apply a weight that decreases with height. Up to 1 km, the weight gradually drops from 1.0 to 0.5, and for anything higher, it stays at 0.5. The idea is that winds higher up have less influence on surface gusts, but they’re still somewhat relevant. This weight factor is straight out of the RUC20 model, which scales down the effect of higher-level winds on surface gusts.

Adjusting and Finding the Max Gust:
For each level, we adjust the surface wind by adding the weighted difference between that level’s wind speed and the 10-meter wind. This gives us an adjusted gust speed for each level. Then, we just keep the maximum gust from all the levels, which, again, lines up with Brasseur’s method – he suggests looking at the strongest possible wind from within the PBL as the potential gust.

Why This Works
This code basically blends NOAA RUC20’s height-based scaling with Brasseur’s idea of pulling from all levels within the PBL. By combining these two methods, we’re able to get a solid estimate of what the maximum gust at the surface might look like, especially during turbulent conditions.


# Calculate the wind speed at each level and find the maximum gust speed
max_gust_speed = np.zeros_like(wspd_wdir10) # Assuming wspd_wdir10 is the surface wind speed


# Iterate over each vertical level
for k in range(u.shape[0]):
# Calculate the wind speed at level k

wind_speed_k = np.sqrt(u[k]**2 + v[k]**2)
height_k = height[k]

# Calculate the weight based on height

weight_k = np.where(height_k <= 1000, 1 - height_k / 2000, 0.5)

# Adjust the wind speed by the weighted difference from the surface

adjusted_speed_k = wspd_wdir10 + (wind_speed_k - wspd_wdir10) * weight_k

# Update the maximum gust speed
max_gust_speed = np.maximum(max_gust_speed, adjusted_speed_k)
Click to expand...
I'm sharing this python code so that other meteorologists can test it and see how it compares to surface observations.
 
I don't think this wind gust calculation proposed by Rosumalski is included in WRF.
If the option "output_diagnostics" is turned on, WRF will output maximum winds (U10MAx and V10MAX) at 10-m height. However, this is apparently different to the "wind gust" concept he proposed.
 
Ming Chen said:
I don't think this wind gust calculation proposed by Rosumalski is included in WRF.
If the option "output_diagnostics" is turned on, WRF will output maximum winds (U10MAx and V10MAX) at 10-m height. However, this is apparently different to the "wind gust" concept he proposed.
Click to expand...
Its something I think could be added to the WRF framework if it turns out to be valid. I can provide the source material if needed.

@Ming Chen
 
William.Hatheway said:
Its something I think could be added to the WRF framework if it turns out to be valid. I can provide the source material if needed.

@Ming Chen
Click to expand...
@Ming Chen

I made an error here's the source material and the new code.

Code: 
SUBROUTINE CALGUST(LPBL,ZPBL,GUST)
C$$$ SUBPROGRAM DOCUMENTATION BLOCK
C . . .
C SUBPROGRAM: CALGUST COMPUTE MAX WIND LEVEL
C PRGRMMR: MANIKIN ORG: W/NP2 DATE: 97-03-04
C
C ABSTRACT:
C THIS ROUTINE COMPUTES SURFACE WIND GUST BY MIXING
C DOWN MOMENTUM FROM THE LEVEL AT THE HEIGHT OF THE PBL
C
C
C PROGRAM HISTORY LOG:
C 03-10-15 GEOFF MANIKIN
C 05-03-09 H CHUANG - WRF VERSION
C 05-06-30 R ROZUMALSKI - DYNAMIC MEMORY ALLOCATION AND SMP
C THREAD-SAFE VERSION
C
C USAGE: CALL CALGUST(GUST)
C INPUT ARGUMENT LIST:
C NONE
C
C OUTPUT ARGUMENT LIST:
C GUST - SPEED OF THE MAXIMUM SFC WIND GUST
C
C OUTPUT FILES:
C NONE
C
C SUBPROGRAMS CALLED:
C UTILITIES:
C H2V
C
C LIBRARY:
C COMMON -
C LOOPS
C OPTIONS
C MASKS
C INDX
C
C ATTRIBUTES:
C LANGUAGE: FORTRAN 90
C MACHINE : CRAY C-90
C$$$
C
C
use vrbls3d
use vrbls2d
C
C INCLUDE ETA GRID DIMENSIONS. SET/DERIVE PARAMETERS.
C
!
INCLUDE "params"
C
INCLUDE "CTLBLK.comm"
C
C DECLARE VARIABLES.
C
INTEGER :: LPBL(IM,JM)
REAL :: GUST(IM,JM)
REAL ZPBL(IM,jsta_2l:jend_2u)
C
C
C*****************************************************************************
C START CALMXW HERE.
C
C LOOP OVER THE GRID.
C
DO J=JSTA,JEND
DO I=1,IM
! GUST(I,J) = SPVAL
GUST(I,J) = 0.
ENDDO
ENDDO
C
C ASSUME THAT U AND V HAVE UPDATED HALOS
C
!$omp parallel do
!$omp& private(ie,iw,mxww,u0,v0,wind)
DO 20 J=JSTA_M,JEND_M
DO 20 I=2,IM-1
L=LPBL(I,J)
IF(MODELNAME .EQ. 'NMM')THEN
IE=I+MOD(J+1,2)
IW=I+MOD(J+1,2)-1

USFC=D25*(U10(I,J-1)+U10(IW,J)+
X U10(IE,J)+U10(I,J+1))
VSFC=D25*(V10(I,J-1)+V10(IW,J)+
X V10(IE,J)+V10(I,J+1))
SFCWIND=SQRT(USFC**2 + VSFC**2)
U0 = D25*(U(I,J-1,L)+U(IW,J,L)+
X U(IE,J,L)+U(I,J+1,L))
V0 = D25*(V(I,J-1,L)+V(IW,J,L)+
X V(IE,J,L)+V(I,J+1,L))
WIND=SQRT(U0**2 + V0**2)

ELSE IF(MODELNAME .EQ. 'NCAR')THEN
USFC=U10(I,J)
VSFC=V10(I,J)
SFCWIND=SQRT(USFC**2 + VSFC**2)
U0=U(I,J,L)
V0=V(I,J,L)
WIND=SQRT(U0**2 + V0**2)
END IF
DELWIND=WIND - SFCWIND
ZSFC=FIS(I,J)*GI
DELWIND=DELWIND*(1.0-AMIN1(0.5,ZPBL(I,J)/2000.))
GUST(I,J)=SFCWIND+DELWIND
10 CONTINUE
20 CONTINUE

C END OF ROUTINE.
C
RETURN
END

1731062878784.png

Key References​

  1. Benjamin, S. G., et al. (2004). An hourly assimilation–forecast cycle: The RUC. Monthly Weather Review, 132(2), 495–518.
    • This paper documents the structure and science behind the RUC model, including the methods for diagnosing surface wind gusts by mixing down PBL winds.
  2. Benjamin, S. G., et al. (2016). A North American hourly assimilation and model forecast cycle: The Rapid Refresh. Monthly Weather Review, 144(4), 1669–1694.
    • This study describes the RAP model, which evolved from RUC, and includes refinements in diagnosing surface wind gusts using PBL heights and layer influences.
  3. Panofsky, H. A., & Dutton, J. A. (1984). Atmospheric Turbulence: Models and Methods for Engineering Applications. John Wiley & Sons.
    • Provides background on turbulence and momentum transfer in the atmosphere, key concepts used in PBL wind gust parameterization.
  4. Stull, R. B. (1988). An Introduction to Boundary Layer Meteorology. Springer.
    • This book explains boundary layer dynamics, including how wind speeds vary with height and can be used to estimate surface gusts based on the PBL.

Python: 
    # Calculate the 10-meter wind speed
    wspd_wdir10 = np.sqrt((u10**2) + (v10**2))
    wspd_wdir = np.sqrt((u**2) + (v**2))
    pblh = wrf.getvar(ncfile, 'PBLH', meta=True)
    height = wrf.getvar(ncfile, 'height_agl', meta=True, units="m")

    # Initialize the maximum gust speed array (same shape as surface wind speed)
    surface_gust_speed = np.zeros_like(wspd_wdir10)  # Assuming wspd_wdir10 is the surface wind speed

    # Iterate over each vertical level
    for k in range(u.shape[0]):
        # Calculate the wind speed at level k
        wind_speed_k = np.sqrt(u[k]**2 + v[k]**2)
        height_k = np.minimum(height[k], 1000)  # Cap height at 1000 meters
    
        # Calculate the weight based on capped height
        weight_k = (1 - (height_k / 2000))  # For heights > 1000m, weight_k will be 0.5
    
        # Adjust the wind speed by the weighted difference from the surface
        surface_gust_speed = wspd_wdir10 + ((wind_speed_k - wspd_wdir10) * weight_k)
      
    # Apply Gaussian smoothing and convert to mph
    gusts = gaussian_filter(surface_gust_speed * 2.23694, sigma=1)  # Convert m/s to mph


made using wrf-python
 
William.Hatheway said:
@Ming Chen

I made an error here's the source material and the new code.

Code: 
SUBROUTINE CALGUST(LPBL,ZPBL,GUST)
C$$$ SUBPROGRAM DOCUMENTATION BLOCK
C . . .
C SUBPROGRAM: CALGUST COMPUTE MAX WIND LEVEL
C PRGRMMR: MANIKIN ORG: W/NP2 DATE: 97-03-04
C
C ABSTRACT:
C THIS ROUTINE COMPUTES SURFACE WIND GUST BY MIXING
C DOWN MOMENTUM FROM THE LEVEL AT THE HEIGHT OF THE PBL
C
C
C PROGRAM HISTORY LOG:
C 03-10-15 GEOFF MANIKIN
C 05-03-09 H CHUANG - WRF VERSION
C 05-06-30 R ROZUMALSKI - DYNAMIC MEMORY ALLOCATION AND SMP
C THREAD-SAFE VERSION
C
C USAGE: CALL CALGUST(GUST)
C INPUT ARGUMENT LIST:
C NONE
C
C OUTPUT ARGUMENT LIST:
C GUST - SPEED OF THE MAXIMUM SFC WIND GUST
C
C OUTPUT FILES:
C NONE
C
C SUBPROGRAMS CALLED:
C UTILITIES:
C H2V
C
C LIBRARY:
C COMMON -
C LOOPS
C OPTIONS
C MASKS
C INDX
C
C ATTRIBUTES:
C LANGUAGE: FORTRAN 90
C MACHINE : CRAY C-90
C$$$
C
C
use vrbls3d
use vrbls2d
C
C INCLUDE ETA GRID DIMENSIONS. SET/DERIVE PARAMETERS.
C
!
INCLUDE "params"
C
INCLUDE "CTLBLK.comm"
C
C DECLARE VARIABLES.
C
INTEGER :: LPBL(IM,JM)
REAL :: GUST(IM,JM)
REAL ZPBL(IM,jsta_2l:jend_2u)
C
C
C*****************************************************************************
C START CALMXW HERE.
C
C LOOP OVER THE GRID.
C
DO J=JSTA,JEND
DO I=1,IM
! GUST(I,J) = SPVAL
GUST(I,J) = 0.
ENDDO
ENDDO
C
C ASSUME THAT U AND V HAVE UPDATED HALOS
C
!$omp parallel do
!$omp& private(ie,iw,mxww,u0,v0,wind)
DO 20 J=JSTA_M,JEND_M
DO 20 I=2,IM-1
L=LPBL(I,J)
IF(MODELNAME .EQ. 'NMM')THEN
IE=I+MOD(J+1,2)
IW=I+MOD(J+1,2)-1

USFC=D25*(U10(I,J-1)+U10(IW,J)+
X U10(IE,J)+U10(I,J+1))
VSFC=D25*(V10(I,J-1)+V10(IW,J)+
X V10(IE,J)+V10(I,J+1))
SFCWIND=SQRT(USFC**2 + VSFC**2)
U0 = D25*(U(I,J-1,L)+U(IW,J,L)+
X U(IE,J,L)+U(I,J+1,L))
V0 = D25*(V(I,J-1,L)+V(IW,J,L)+
X V(IE,J,L)+V(I,J+1,L))
WIND=SQRT(U0**2 + V0**2)

ELSE IF(MODELNAME .EQ. 'NCAR')THEN
USFC=U10(I,J)
VSFC=V10(I,J)
SFCWIND=SQRT(USFC**2 + VSFC**2)
U0=U(I,J,L)
V0=V(I,J,L)
WIND=SQRT(U0**2 + V0**2)
END IF
DELWIND=WIND - SFCWIND
ZSFC=FIS(I,J)*GI
DELWIND=DELWIND*(1.0-AMIN1(0.5,ZPBL(I,J)/2000.))
GUST(I,J)=SFCWIND+DELWIND
10 CONTINUE
20 CONTINUE

C END OF ROUTINE.
C
RETURN
END

View attachment 16239

Key References​

  1. Benjamin, S. G., et al. (2004). An hourly assimilation–forecast cycle: The RUC. Monthly Weather Review, 132(2), 495–518.
    • This paper documents the structure and science behind the RUC model, including the methods for diagnosing surface wind gusts by mixing down PBL winds.
  2. Benjamin, S. G., et al. (2016). A North American hourly assimilation and model forecast cycle: The Rapid Refresh. Monthly Weather Review, 144(4), 1669–1694.
    • This study describes the RAP model, which evolved from RUC, and includes refinements in diagnosing surface wind gusts using PBL heights and layer influences.
  3. Panofsky, H. A., & Dutton, J. A. (1984). Atmospheric Turbulence: Models and Methods for Engineering Applications. John Wiley & Sons.
    • Provides background on turbulence and momentum transfer in the atmosphere, key concepts used in PBL wind gust parameterization.
  4. Stull, R. B. (1988). An Introduction to Boundary Layer Meteorology. Springer.
    • This book explains boundary layer dynamics, including how wind speeds vary with height and can be used to estimate surface gusts based on the PBL.

Python: 
    # Calculate the 10-meter wind speed
    wspd_wdir10 = np.sqrt((u10**2) + (v10**2))
    wspd_wdir = np.sqrt((u**2) + (v**2))
    pblh = wrf.getvar(ncfile, 'PBLH', meta=True)
    height = wrf.getvar(ncfile, 'height_agl', meta=True, units="m")

    # Initialize the maximum gust speed array (same shape as surface wind speed)
    surface_gust_speed = np.zeros_like(wspd_wdir10)  # Assuming wspd_wdir10 is the surface wind speed

    # Iterate over each vertical level
    for k in range(u.shape[0]):
        # Calculate the wind speed at level k
        wind_speed_k = np.sqrt(u[k]**2 + v[k]**2)
        height_k = np.minimum(height[k], 1000)  # Cap height at 1000 meters
   
        # Calculate the weight based on capped height
        weight_k = (1 - (height_k / 2000))  # For heights > 1000m, weight_k will be 0.5
   
        # Adjust the wind speed by the weighted difference from the surface
        surface_gust_speed = wspd_wdir10 + ((wind_speed_k - wspd_wdir10) * weight_k)
     
    # Apply Gaussian smoothing and convert to mph
    gusts = gaussian_filter(surface_gust_speed * 2.23694, sigma=1)  # Convert m/s to mph


made using wrf-python
Click to expand...
I think there is another error in this I will review again
 
RalphIPTony said:
Great topic for discussion, thanks for bringing this up!
Click to expand...
Python: 
# Retrieve the planetary boundary layer height (PBLH) from the netCDF file
    # PBLH represents the top height of the atmospheric boundary layer
    pblh = wrf.getvar(ncfile, 'PBLH', meta=True)

    # Retrieve the height above ground level at different atmospheric levels
    # This provides information on the altitude of each vertical level
    height = wrf.getvar(ncfile, 'height_agl', meta=True, units="m")

    # Define the capped height as the minimum of the PBLH and 1000 meters
    capped_pblh = np.minimum(pblh, 1000)

    # Interpolate the wind speed to the top of the planetary boundary layer or 1000m, whichever is lower
    # This corresponds to vPBL in the wind gust parameterization formula
    wspd_pblh = wrf.interplevel(wspd_wdir, height, capped_pblh)
  
    # Calculate the weight based on the capped height
    # This weight decreases linearly with height, reducing the influence of wind at higher levels
    weight_k = (1 - (capped_pblh / 2000))  # For heights ≥1000m, weight_k becomes 0.5
  
    # Adjust the surface gust speed based on the parameterization formula
    # This combines the surface wind speed (10-meter level) with the weighted difference between
    # the wind speed at the top of the PBL and the surface wind speed
    surface_gust_speed = wspd_wdir10 + ((wspd_pblh - wspd_wdir10) * weight_k)

I think this works better
 
You must log in or register to reply here.
Top