Project

General

Profile

Heat Index Calculation not right

Added by Bryan Kilgore 8 months ago

We are trying to run our script to output grib for heat index based off 2 meter temp and Rh. When we input the script into CDO we are returning a value that is about 3 degrees too high compared to java.

CDO Command:
cdo -setrtomiss,264,265 -expr,' _tF = ((2t - 273.15) * 1.8) + 32; _hi = -42.379 + (2.04901523 * _tF) + (10.14333127 * 2r) - (0.22475541 * _tF * 2r) - (6.83783 * pow(10, -3) * pow(_tF, 2)) - (5.481717 * pow(10, -2) * pow(2r, 2)) + (1.22874 * pow(10, -3) * pow(_tF, 2) * 2r) + (8.5282 * pow(10, -4) * _tF * pow(2r, 2)) - (1.99 * pow(10, -6) * pow(_tF, 2) * pow(2r, 2)); _hi = _tF < 80 ? _tF : _hi; heatx = (_hi - 32) * 0.55555556 + 273.15; ' -invertlat ./temp2m_rh.grib2 ./heat.grib2

When we interigate a point on the grib it is:
2mt = 83.25
rh = 74.11
heat = 94

Run this same formula in java console.
var _tF = 83.25
var r2 = 74.11
console.log(-42.379 + 2.04901523 * _tF + 10.14333127 * r2 - 0.22475541 * _tF * r2 - 6.83783 * (Math.pow(10, -3) * Math.pow(_tF, 2)) - 5.481717 * (Math.pow(10, -2) * Math.pow(r2, 2)) + 1.22874 * (Math.pow(10, -3) * Math.pow(_tF, 2) * r2) + 8.5282 * (Math.pow(10, -4) * _tF * Math.pow(r2, 2)) - 1.99 * (Math.pow(10, -6) * Math.pow(_tF, 2) * Math.pow(r2, 2)))

Output: 90.09733414634763

The java output matches the WPC Heat Index calculator.

https://www.wpc.ncep.noaa.gov/html/heatindex.shtml

We are running following version:
Climate Data Operators version 2.1.1 (https://mpimet.mpg.de/cdo)
System: aarch64-apple-darwin22.3.0
CXX Compiler: /usr/bin/clang++ -std=gnu++17 -pipe -Os -stdlib=libc++ -isysroot/Library/Developer/CommandLineTools/SDKs/MacOSX13.sdk -arch arm64 -pthread
CXX version : Apple clang version 14.0.0 (clang-1400.0.29.202)
C Compiler: /usr/bin/clang -pipe -Os -isysroot/Library/Developer/CommandLineTools/SDKs/MacOSX13.sdk -arch arm64 -pthread -pthread
C version : Apple clang version 14.0.0 (clang-1400.0.29.202)
F77 Compiler: -pipe -Os -m64
Features: 32GB 10threads c++17 pthreads HDF5 NC4/HDF5 OPeNDAP sz udunits2 proj xml2 magics curl fftw3
Libraries: yac/2.6.1 HDF5/1.14.0 proj/8.2.1 xml2/2.10.3 curl/7.87.0 magics/4.5.1
CDI data types: SizeType=size_t
CDI file types: srv ext ieg grb1 grb2 nc1 nc2 nc4 nc4c nc5 nczarr
CDI library version : 2.1.1
ecCodes library version : 2.28.0
NetCDF library version : 4.9.0 of Feb 17 2023 15:11:40 $
HDF5 library version : 1.14.0
exse library version : 1.4.2
FILE library version : 1.9.1

We are using this version because when we upgraded to any version beyond this one other features break.


Replies (10)

RE: Heat Index Calculation not right - Added by Karin Meier-Fleischer 8 months ago

Hi Bryan,

if I run the expression with the given _tF and r2 value I get the correct value.
I used Python for testing:

1. check your Java code result

_tF = 83.25
r2  = 74.11

result = -42.379 + 2.04901523 * _tF + 10.14333127 * r2 - 0.22475541 * \
         _tF * r2 - 6.83783 * (math.pow(10, -3) * math.pow(_tF, 2)) - \
         5.481717 * (math.pow(10, -2) * math.pow(r2, 2)) + 1.22874 * \
         (math.pow(10, -3) * math.pow(_tF, 2) * r2) + 8.5282 * \
         (math.pow(10, -4) * _tF * math.pow(r2, 2)) - 1.99 * (math.pow(10, -6) * \
          math.pow(_tF, 2) * math.pow(r2, 2))

result = 90.09733414634763

2. check cdo result

_tF = str(_tF)
r2  = str(r2)

cmd0 = "'heat= (-42.379 + 2.04901523 * " + _tF + " + 10.14333127 * " + r2 + "- 0.22475541 * " + \
        _tF + " * " + r2 + "- 6.83783 * (pow(10, -3) * pow(" + _tF + ", 2)) - " + \
         "5.481717 * (pow(10, -2) * pow("+r2+", 2)) + 1.22874 * " + \
         "(pow(10, -3) * pow(" + _tF + ", 2) * " + r2 + ") + 8.5282 * " + \
         "(pow(10, -4) * " + _tF + " * pow(" + r2 + ", 2)) - 1.99 * (pow(10, -6) * " + \
         "pow("+ _tF + ", 2) * pow("+ r2 +", 2)))'" 

cdo.infon(input="-expr,"+cmd0+" -chname,topo,heat -topo")  # use -topo as dummy input, it won't be used

​['-1 :       Date     Time   Level Gridsize    Miss :     Minimum        Mean     Maximum : Parameter name',
 '1 : 0001-01-01 00:00:00       0   259200       0 :      90.097      90.097      90.097 : heat']

RE: Heat Index Calculation not right - Added by Bryan Kilgore 8 months ago

If you output the grib does it also give you the same value? We are wondering if there is something wrong with the grib output causing this.

Lat/Lng: -86.317, 39.25

RE: Heat Index Calculation not right - Added by Karin Meier-Fleischer 8 months ago

Send me the cdo command line you are using to compute the above heat, please.

RE: Heat Index Calculation not right - Added by Bryan Kilgore 8 months ago

It is in the original post.

RE: Heat Index Calculation not right - Added by Karin Meier-Fleischer 8 months ago

I asked because it gives me a complete other result.

cdo -setrtomiss,264,265 -expr,' _tF = ((2t - 273.15) * 1.8) + 32; _hi = -42.379 + (2.04901523 * _tF) + (10.14333127 * 2r) - (0.22475541 * _tF * 2r) - (6.83783 * pow(10, -3) * pow(_tF, 2)) - (5.481717 * pow(10, -2) * pow(2r, 2)) + (1.22874 * pow(10, -3) * pow(_tF, 2) * 2r) + (8.5282 * pow(10, -4) * _tF * pow(2r, 2)) - (1.99 * pow(10, -6) * pow(_tF, 2) * pow(2r, 2)); _hi = _tF < 80 ? _tF : _hi; heatx = (_hi - 32) * 0.55555556 + 273.15; ' -invertlat ./temp2m_rh.grib2 ./heat.grib2
cdo -s -outputtab,name,value -remapnn,lon=-86.317_lat=39.25 heat.grib2

Result:

#    name    value 
   heatx 296.6302

RE: Heat Index Calculation not right - Added by Bryan Kilgore 8 months ago

That's why I'm making this post because something isn't right. The commands I posted was straight copy and paste from my terminals. It seems like there is something wrong when it writes out to GRIB.

RE: Heat Index Calculation not right - Added by Karin Meier-Fleischer 8 months ago

The results differ because the given cdo expression and the java expression are not the same. Please, check it, this is up to you.

RE: Heat Index Calculation not right - Added by Bryan Kilgore 8 months ago

Ok let me start back over.

I ran this CDO command:
cdo -setrtomiss,264,265 -expr,' _tF = ((2t - 273.15) * 1.8) + 32; _hi = -42.379 + (2.04901523 * _tF) + (10.14333127 * 2r) - (0.22475541 * _tF * 2r) - (6.83783 * pow(10, -3) * pow(_tF, 2)) - (5.481717 * pow(10, -2) * pow(2r, 2)) + (1.22874 * pow(10, -3) * pow(_tF, 2) * 2r) + (8.5282 * pow(10, -4) * _tF * pow(2r, 2)) - (1.99 * pow(10, -6) * pow(_tF, 2) * pow(2r, 2)); _hi = _tF < 80 ? _tF : _hi; heatx = (_hi - 32) * 0.55555556 + 273.15; ' -invertlat ./temp2m_rh.grib2 ./heat.grib2

I have attached the grib output from the above output.

This has 1 extra part for if the temperature is below 80 degrees. In the following command that isn't necessary since we are specifiying the temperature.
var _tF = 83.25
var r2 = 74.11
console.log(-42.379 + 2.04901523 * _tF + 10.14333127 * r2 - 0.22475541 * _tF * r2 - 6.83783 * (Math.pow(10, -3) * Math.pow(_tF, 2)) - 5.481717 * (Math.pow(10, -2) * Math.pow(r2, 2)) + 1.22874 * (Math.pow(10, -3) * Math.pow(_tF, 2) * r2) + 8.5282 * (Math.pow(10, -4) * _tF * Math.pow(r2, 2)) - 1.99 * (Math.pow(10, -6) * Math.pow(_tF, 2) * Math.pow(r2, 2)))

The Lat/Lng: -86.317, 39.25 that I interrogated this spot in QGIS and received this result.

RE: Heat Index Calculation not right - Added by Karin Meier-Fleischer 8 months ago

What is the difference between heat and heatx? What about 'heatx = (_hi - 32) * 0.55555556 + 273.15'?

RE: Heat Index Calculation not right - Added by Bryan Kilgore 8 months ago

heatx is CDO parameter variable so it codes the output grib as heat index product.
heat in your screenshot is coming from the file name in QGIS.

    (1-10/10)