/* * Copyright (c) 2019 by Alfonso Crisci, Marco Morabito e Alessandro Messeri
* * @authors Alfonso Crisci, Marco Morabito e Alessandro Messeri
* * Corresponet Email: a.crisci@ibe.cnr.it
* * This program is free software; you can redistribute it and/or
* * modify it under the terms of the GNU General Public License
* * as published by the Free Software Foundation; either version 2
* * of the License, or (at your option) any later version,
* * provided that any use properly credits the author.
* * This program is distributed in the hope that it will be useful,
* * but WITHOUT ANY WARRANTY; without even the implied warranty of
* * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* * GNU General Public License for more details at http://www.gnu.org * * */
/*
* Define global constants for functions.
*/
var Patm = 101325.0; // pascal
var CpAir= 1004.0;
var CpWat= 4186.0;
var CpVap= 1805.0;
var Hfg= 2501000.0;
var RAir= 287.055;
var TKelConv= 273.15;
var PI = Math.PI,
pi = Math.PI,
twopi = 2 * pi,
rad = PI / 180,
radians = PI / 180,
degrees = 180 / PI,
J1970 = 2440588,
J2000 = 2451545;
/*
* Format functions.
*/
function OneDec(c)
{
return Math.round(10 * c) / 10;
}
function TwoDec(c)
{
return Math.round(100 * c) / 100;
}
function ThreeDec(c)
{
return Math.round(1000 * c) / 1000;
}
function FourDec(c)
{
return Math.round(10000 * c) / 10000;
}
function scientificNotation(c,e)
{
return c.toPrecision(e);
}
/** * @(#)pnorm.js and qnorm.js
* * Copyright (c) 2000 by Sundar Dorai-Raj
* * @author Sundar Dorai-Raj
* * Email: sdoraira@vt.edu
* * This program is free software; you can redistribute it and/or
* * modify it under the terms of the GNU General Public License
* * as published by the Free Software Foundation; either version 2
* * of the License, or (at your option) any later version,
* * provided that any use properly credits the author.
* * This program is distributed in the hope that it will be useful,
* * but WITHOUT ANY WARRANTY; without even the implied warranty of
* * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* * GNU General Public License for more details at http://www.gnu.org * * */
// "Applied Statistics Algorithms" (1985)
// P. Griffiths and I. D. Hill, editor
function pnorm(z,tail) {
// Algorithm AS66 Applied Statistics (1973) vol 22 no.3
// Computes P(Z<z)
if( tail === undefined ) { tail = false;};
z=parseFloat(z);
var ltone= 7.0;
var utzero= 18.66;
con= 1.28;
a1 = 0.398942280444;
a2 = 0.399903438504;
a3 = 5.75885480458;
a4 = 29.8213557808;
a5 = 2.62433121679;
a6 = 48.6959930692;
a7 = 5.92885724438;
b1 = 0.398942280385;
b2 = 3.8052e-8;
b3 = 1.00000615302;
b4 = 3.98064794e-4;
b5 = 1.986153813664;
b6 = 0.151679116635;
b7 = 5.29330324926;
b8 = 4.8385912808;
b9 = 15.1508972451;
b10= 0.742380924027;
b11= 30.789933034;
b12= 3.99019417011;
if(z<0) {
tail=!tail;
z=-z;
}
if(z<=ltone || tail && z<=utzero) {
y=0.5*z*z;
if(z>con) {
alnorm=b1*Math.exp(-y)/(z-b2+b3/(z+b4+b5/(z-b6+b7/(z+b8-b9/(z+b10+b11/(z+b12))))));
}
else {
alnorm=0.5-z*(a1-a2*y/(y+a3-a4/(y+a5+a6/(y+a7))));
}
}
else {
alnorm=0;
}
if(!tail) alnorm=1-alnorm;
return(alnorm);
}
function qnorm(p) {
// ALGORITHM AS 111, APPL.STATIST., VOL.26, 118-121, 1977.
// Computes z=invNorm(p)
p=parseFloat(p);
split=0.42;
a0= 2.50662823884;
a1=-18.61500062529;
a2= 41.39119773534;
a3=-25.44106049637;
b1= -8.47351093090;
b2= 23.08336743743;
b3=-21.06224101826;
b4= 3.13082909833;
c0= -2.78718931138;
c1= -2.29796479134;
c2= 4.85014127135;
c3= 2.32121276858;
d1= 3.54388924762;
d2= 1.63706781897;
q=p-0.5;
if(Math.abs(q)<=split) {
r=q*q;
ppnd=q*(((a3*r+a2)*r+a1)*r+a0)/((((b4*r+b3)*r+b2)*r+b1)*r+1);
}
else {
r=p;
if(q>0) r=1-p;
if(r>0) {
r=Math.sqrt(-Math.log(r));
ppnd=(((c3*r+c2)*r+c1)*r+c0)/((d2*r+d1)*r+1);
if(q<0) ppnd=-ppnd;
}
else {
ppnd=0;
}
}
return(ppnd);
}
/**
* Air Pressure a elevation
*
* @param {number} press,topo
* @return {number}
*
*/
function pheight(press,topo)
{
var pressalt= press * Math.pow((1-(2.25577*(0.00001)*topo)),5.25588);
return(pressalt)
}
/**
* Error function
*
* @param {number} x
* @return {number}
*
*/
function erf(x) {x=parseFloat(x);
return(2 * pnorm(x * Math.sqrt(2)) - 1);}
/**
* Wind reduction in boundary layer.
* @param {number} x,ref,fin
* @return {number}
*
*/
function reducewind(x,ref,fin) {
if( ref === undefined ) { tresh = 10;};
if( fin === undefined ) { fin = 2;};
return(x*1/(Math.log(ref/0.01)/Math.log(fin/0.01)));
}
///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
// Math & Trigonometric functions.
function radToDeg(angleRad)
{
return (180.0 * angleRad / Math.PI);
}
function degToRad(angleDeg)
{
return (Math.PI * angleDeg / 180.0);
}
function compass_16(direction)
{
var dir = ["N", "NNE", "NE", "ENE", "E", "ESE", "SE",
"SSE", "S", "SSW", "SW", "WSW", "W", "WNW",
"NW", "NNW"];
var dirint=(((direction+11.25)/22.5));
return(dir[parseInt(dirint % 16)]);
}
function compass_8(direction)
{ var dir = ["N", "NE", "E", "SE",
"S", "SW", "W", "NW"];
var dirint=(((direction+22.5)/45));
return(dir[parseInt(dirint % 8)]);
}
///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
// Math & Trigonometric functions.
String.prototype.toDate = function(format)
{
var normalized = this.replace(/[^a-zA-Z0-9]/g, '-');
var normalizedFormat= format.toLowerCase().replace(/[^a-zA-Z0-9]/g, '-');
var formatItems = normalizedFormat.split('-');
var dateItems = normalized.split('-');
var monthIndex = formatItems.indexOf("mm");
var dayIndex = formatItems.indexOf("dd");
var yearIndex = formatItems.indexOf("yyyy");
var hourIndex = formatItems.indexOf("hh");
var minutesIndex = formatItems.indexOf("ii");
var secondsIndex = formatItems.indexOf("ss");
var today = new Date();
var year = yearIndex>-1 ? dateItems[yearIndex] : today.getFullYear();
var month = monthIndex>-1 ? dateItems[monthIndex]-1 : today.getMonth()-1;
var day = dayIndex>-1 ? dateItems[dayIndex] : today.getDate();
var hour = hourIndex>-1 ? dateItems[hourIndex] : today.getHours();
var minute = minutesIndex>-1 ? dateItems[minutesIndex] : today.getMinutes();
var second = secondsIndex>-1 ? dateItems[secondsIndex] : today.getSeconds();
return new Date(year,month,day,hour,minute,second);
};
function sun_data(strtime,lat,lon,parameter)
{
var datetime=new Date(strtime);
var udtTimedHours=datetime.getHours() - 0;
var udtTimedMinutes =datetime.getMinutes() - 0;
var udtTimedSeconds = datetime.getSeconds() - 0;
var udtLocationdLongitude = lon - 0;
var udtLocationdLatitude = lat - 0;
var dEarthMeanRadius = 6371.01;
var dAstronomicalUnit = 149597890;
var dDecimalHours = udtTimedHours + (udtTimedMinutes + udtTimedSeconds / 60) / 60;
var dJulianDate = datetime.valueOf() / (1000.0 *60.0*60.0 * 24.0) - 0.5 + J1970;
var dElapsedJulianDays = dJulianDate - 2451545;
var dOmega = 2.1429 - 0.0010394594 * dElapsedJulianDays;
var dMeanLongitude = 4.895063 + 0.017202791698 * dElapsedJulianDays;
var dMeanAnomaly = 6.24006 + 0.0172019699 * dElapsedJulianDays;
var dEclipticLongitude = dMeanLongitude + 0.03341607 * Math.sin(dMeanAnomaly) + 3.4894E-4 * Math.sin(2 * dMeanAnomaly) - 1.134E-4 - 2.03E-5 * Math.sin(dOmega);
var dEclipticObliquity = 0.4090928 - 6.214E-9 * dElapsedJulianDays + 3.96E-5 * Math.cos(dOmega);
var dSin_EclipticLongitude = Math.sin(dEclipticLongitude);
var dY = Math.cos(dEclipticObliquity) * dSin_EclipticLongitude;
var dX = Math.cos(dEclipticLongitude);
var dRightAscension = Math.atan2(dY, dX);
0 > dRightAscension && (dRightAscension += twopi);
var dDeclination = Math.asin(Math.sin(dEclipticObliquity) * dSin_EclipticLongitude);
var dGreenwichMeanSiderealTime = 6.6974243242 + 0.0657098283 * dElapsedJulianDays + dDecimalHours;
var dLocalMeanSiderealTime = (15 * dGreenwichMeanSiderealTime +udtLocationdLongitude) * rad;
var dHourAngle = dLocalMeanSiderealTime - dRightAscension;
var dLatitudeInRadians = udtLocationdLatitude * rad;
var dCos_Latitude = Math.cos(dLatitudeInRadians);
var dSin_Latitude = Math.sin(dLatitudeInRadians);
var dCos_HourAngle = Math.cos(dHourAngle);
var udtSunCoordinatesdZenithAngle = Math.acos(dCos_Latitude * dCos_HourAngle * Math.cos(dDeclination) + Math.sin(dDeclination) * dSin_Latitude);
dY = -Math.sin(dHourAngle);
dX = Math.tan(dDeclination) * dCos_Latitude - dSin_Latitude * dCos_HourAngle;
var udtSunCoordinatesdAzimuth = Math.atan2(dY,dX);
0 > udtSunCoordinatesdAzimuth && (udtSunCoordinatesdAzimuth += twopi);
udtSunCoordinatesdAzimuth /= rad;
var dParallax = dEarthMeanRadius / dAstronomicalUnit * Math.sin(udtSunCoordinatesdZenithAngle);
udtSunCoordinatesdZenithAngle = (udtSunCoordinatesdZenithAngle + dParallax) / rad;
var azimuth = udtSunCoordinatesdAzimuth;
var zenith = udtSunCoordinatesdZenithAngle;
var elevation = 90 - udtSunCoordinatesdZenithAngle;
if (elevation > 85.0) {var refractionCorrection = 0.0;}
else {
var te = Math.tan (degToRad(elevation));
if (elevation > 5.0)
{var refractionCorrection = 58.1 / te - 0.07 / (te*te*te) + 0.000086 / (te*te*te*te*te);}
else if (elevation > -0.575)
{var refractionCorrection = 1735.0 + elevation * (-518.2 + elevation * (103.4 + elevation * (-12.79 + elevation * 0.711) ) );}
else {var refractionCorrection = -20.774 / te;}
refractionCorrection = refractionCorrection / 3600.0;
}
var solarZen = zenith - refractionCorrection;
if ( parameter == "azimuth") {return(azimuth)}
else if ( parameter == "zenith") {return(zenith)}
else if ( parameter == "solarZenith") {return(solarZen)}
else if ( parameter == "altitude") {return(elevation)}
else if ( parameter == "declination") {return(dDeclination*(180/PI))}
else if ( parameter == "JD") {return(dJulianDate)}
else { return("Parameter not indicated!")};
}
/**
* Given t air temperature (Celsius degrees), rh relative humidity (%) , wind speed in m per second, global solar radiation in Watt on square meter
* tg globometric temperature gives, solar zenith in radians, pair Air Pressure in millibar (hPa),
* alb_sfc mean albedo surface, ratio diffuse and directed solar radiation and irad if radiation is computed in
* in heat globe realease the function computes
* Wet-bulb globe temperature (WBGT) index following Liljegren scheme .
* @param {number}t,rh,wind,solar,zenith,pair,alb_sfc,fdir,irad.
* @return {number}
* @customfunction
*/
function wbgt_liljegren(t,rh,wind,solar,zenith,tg,pair,alb_sfc,fdir,irad,diam_globe)
{
var wbgt;
if( pair === undefined ) {pair = 1013.25;};
if( alb_sfc === undefined ) {alb_sfc = 0.4;};
if( fdir === undefined ) { fdir = 0.8;};
if( irad === undefined ) {irad = 1;};
if( diam_globe === undefined ) {diam_globe=0.0508;};
if( tg === undefined ) { tg =Tglob_sphere(t,rh,wind,solar,zenith,pair,alb_sfc,fdir,diam_globe)};
var tw = natural_wetbulb(t,rh,wind,solar,zenith,pair,alb_sfc,fdir,irad);
wbgt = 0.7*tw+0.2*tg+0.1*t;
return wbgt;
}
/**
* Given t air temperature (Celsius), rh relative humidity (%) and Tg globometric temperature gives Wet-bulb globe temperature (WBGT) index indoor.
* @param {number} t,rh,tg,pa
* @return {number}
* @customfunction
*/
function wbgt_stull(t,rh,tg)
{
var wbgt;
var tw = wetbulb_stull(t,rh)
wbgt = 0.7*tw+0.2*tg+0.1*t;
return wbgt;
}
/**
* Given tg globometric of sphere .
* @param {number} t,rh,speed,solar,zenith,pair,alb_sfc,fdir,diam
* @return {number}
* @customfunction
*/
function Tglob_sphere(t,rh,speed,solar,zenith,pair,alb_sfc,fdir,diam)
{
if( solar < 5 ) { return(t)};
if(zenith === undefined ) {zenith = 0.0001;};
if(pair === undefined ) { pair =1013.25;};
if(alb_sfc === undefined ) {alb_sfc = 0.4;};
if(fdir === undefined ) { fdir = 0.8;};
if(zenith <= 0 ){ zenith = 0.0000000001;};
if(zenith > 1.57 ){ zenith = 1.57;};
var converge,cza,dT,Tref,h,Tglobe,tmp_globe;
zenith=zenith/(180/Math.PI);
converge = 0.05;
cza = Math.cos(zenith);
var emis_air=emis_atm(t,rh);
var speedmin = 0.1;
var alb_globe = 0.05;
var emis_globe = 0.95;
var emis_sfc = 0.999;
var stefanb = 0.0000000567;
var Tair = t + 273.15;
var Tsfc = t + 273.15;
var Tglobe_prev = t + 273.15;
var iter=1;
do {
Tref = 0.5 * (Tglobe_prev + Tair);
h = h_sphere_in_air(Tref, speed, pair,speedmin, diam);
tmp_globe= 0.5 * (emis_air * Math.pow( Tair,4) + Math.pow(Tsfc,4)) -
h / (emis_globe * stefanb) * (Tglobe_prev - Tair) +
solar / ((2 * emis_globe * stefanb) * (1 - alb_globe) * (fdir * (1 / (2 * cza) - 1) + 1 + alb_sfc));
Tglobe =Math.pow(tmp_globe ,0.25);
dT = Tglobe - Tglobe_prev;
if(Math.abs(dT) < converge) { Tglobe = Tglobe - 273.15;break;} else {Tglobe_prev =Tglobe_prev+0.01;}
iter=iter+1;
} while ( iter < 1000);
if ( iter===1000) {Tglobe=-999};
return(Tglobe);
}
/**
* Given air temperature (Celsius), relative humidity (%) gives natural wetbulb in celsius degrees.
* Formulation from Stull 2011, Journal of Applied Meteorology and Climatology.
* @param {number} t,rh
* @return {number}
* @customfunction
*/
function wetbulb_stull(t,rh)
{
c = new Array(6);
c[0] =0.151977;
c[1] =8.313659;
c[2] =1.676331;
c[3] =0.00391838;
c[4] =0.023101;
c[5] =4.686035;
var wetbulb = t * Math.atan(c[0] * Math.sqrt(rh + c[1]))
+ Math.atan(t + rh) - Math.atan(rh - c[2])
+ c[3] * (Math.pow(rh,3/2)) * Math.atan(c[4] * rh) - c[5];
return(wetbulb)
}
/**
* Given air temperature (Celsius), relative humidity (%) and pressure ( pa) gives natural wetbulb in Ceslius degrees.
* Brice and HALL vapor pressure https://www.weather.gov/epz/wxcalc_rh
* @param {number} t, rh, pair
* @return {number}
* @customfunction
*/
function wetbulb_brice_hall(t,rh,pair)
{
if (pair == undefined) {pair = 1013.25};
var Ewguess,Eguess,wetbulb,cursign;
var Twguess = 0;
var incr = 1;
var previoussign = 1;
var Edifference = 1;
var E2 = pvap(t,rh);
outerloop:
while (Math.abs(Edifference) > 0.005)
{
Ewguess = 6.112 * Math.exp((17.67 * Twguess) / (Twguess + 243.5));
Eguess = Ewguess - (pair) * (t - Twguess) * 0.00066 * (1 + (0.00115 * Twguess));
Edifference = E2 - Eguess;
if (Edifference == 0)
{
break outerloop;
} else {
if (Edifference < 0)
{
cursign = -1;
if (cursign != previoussign)
{
previoussign = cursign;
incr = incr/10;
} else {
incr = incr;
}
} else {
cursign = 1;
if (cursign != previoussign)
{
previoussign = cursign;
incr = incr/10;
} else {
incr = incr;
}
}
}
Twguess = Twguess + incr * previoussign;
}
wetbulb = Twguess;
return wetbulb;
}
/**
* Given t air temperature (Celsius degrees), rh relative humidity (%) , wind speed in m per second,
* global solar radiation in watt on square meters, solar zenith in radians, pair Air Pressure in millibar (hPa),
* alb_sfc mean albedo surface, ratio diffuse and directed solar radiation and irad if radiation is computed in
* in heat globe realease the function computes
* natural wet-bulb temperature index following Liljegren scheme .
* @param {number} t,rh,wind,solar,zenith,pair,alb_sfc,fdir,irad
* @return {number}
* @customfunction
*/
function natural_wetbulb(t,rh,wind,solar,zenith,pair,alb_sfc,fdir,irad)
{
if( solar < 10 ) { solar == 0};
if( zenith === undefined ) {zenith = 0.00001;};
if( pair === undefined ) {pair = 1013.25;};
if( alb_sfc === undefined ) {alb_sfc = 0.4;};
if( fdir === undefined ) { fdir = 0.8;};
if( irad === undefined ) {irad = 1;};
zenith=zenith/(180/Math.PI);
var converge,cza,dT,Tref,h,Twb,Fatm,density,eair;
converge = 0.05;
if(zenith <= 0 ){ zenith = 0.0000000001;};
if(zenith > 1.57 ){ zenith = 1.57;};
var speedmin = 0.1;
var emis_sfc = 0.999;
var stefanb=0.000000056696;
/* heat capaticy of dry air at constant pressure */
var cp=1003.5 ;
var m_air=28.97;
var m_h2o=18.015;
var r_gas=8314.34;
var r_air=r_gas / m_air;
var ratio=cp * m_air/ m_h2o;
var Pr=cp / (cp + (1.25 * r_air));
// Wick constants
var emis_wick=0.95 // emissivity
var alb_wick=0.4 // albedo
var diam_wick=0.007 // diameter (in m)
var len_wick=0.0254 // length (in m)
// Globe constants
var Tair = t + 273.15;
var Tsfc = t + 273.15;
var Twb_prev = dewpoint(t,rh) + 273.15;
var emis_air= emis_atm(t,rh);
eair = (rh/100) * esat(Tair);
density = (pair * 100) / (Tair * r_air);
var iter=1;
do {
Tref = 0.5 * (Twb_prev + Tair);
// Radiative heating term
Fatm = stefanb * emis_wick * (0.5 * (emis_air * Math.pow(Tair, 4)
+ emis_sfc * Math.pow(Tsfc, 4) - Math.pow(Twb_prev, 4))
+ (1 - alb_wick) * solar * ((1 - fdir) * (1 + 0.25 * diam_wick / len_wick)
+ (Math.tan(zenith) / 3.1416) + 0.25 * diam_wick / len_wick) * fdir + alb_sfc);
// Schmidt number
var Sc = viscosity(Tair) / (density * diffusivity(Tref, pair));
// Calculate the convective heat transfer coefficient for a long cylinder in cross flow
h = h_cylinder_in_air(Twb_prev, wind, pair,speedmin, diam_wick);
// Calculate the saturation vapor pressure (hPa) over liquid water
var ewick = esat(Twb_prev);
// Calculate the heat of evaporation, J/(kg K)
var evap = h_evap(Twb_prev);
Twb = Tair - (evap/ratio) * ((ewick - eair) / (pair - ewick)) * Math.pow((Pr / Sc),0.56) + Fatm / h * irad;
dT = Twb - Twb_prev;
if (Math.abs(dT) < converge) { Twb = Twb - 273.15; break;} else {Twb_prev=Twb_prev+0.01;}
iter=iter+1;
} while ( iter < 1000);
if ( iter===1000) {Twb=-999}
return(Twb);
}
/**
* lost_productivity is a function to estimate population heat exposure and impacts on working people
* @param {number} wbgt,tresh
* @return {number}
* Estimating population heat exposure and impacts on working people in conjunction with climate change
* Tord Kjellstrom & Chris Freyberg & Bruno Lemke & Matthias Otto & David Briggs
* Int J Biometeorol (2018) 62:291-306 DOI 10.1007/s00484-017-1407-0
*/
function lost_productivity(wbgt,tresh) {
if( tresh === undefined ) { tresh = 33.5;};
return(0.5*(1+ erf((wbgt-tresh)/(3.75*Math.sqrt(2))))*100)
}
function BMW(h,w)
{return(w/(Math.pow((h/100),2)));}
/**
* Given body surface area by using antropometric features.
* @param {number} h,w
* @return {number}
* @customfunction
*/
function BSA (h,w)
{ return( Math.pow(h/100,0.725)*(0.20247*(Math.pow(w,0.425)))); }
/**
* Given met rate level.
* @param {number} BSA, isolevel
* @return {number}
* @customfunction
*/
function met_rate(BSA, isolevel) {
if ( isolevel === undefined) {isolevel = 1 };
return(BSA*(isolevel*50));
}
/**
* Given met and clo level gives acclimated tresholds of wbgt.
* @param {number} met,clolevel
* @return {number}
* @customfunction
*/
function rel_acclimatized(met,clolevel) {
if ( clolevel === undefined) {isolevel = 1 };
return(56.7-(11.5*Math.LN10(met))-clolevel);
}
/**
* Given met and clo level gives unacclimated tresholds of wbgt.
* @param {number} met,clolevel
* @return {number}
* @customfunction
*/
function ral_unacclimatized(met,clolevel) {
if ( clolevel === undefined) {isolevel = 1 };
return(59.9-(14.1*Math.LN10(met))-clolevel);
}
/* https://www.rapidtables.com/convert/color/hex-to-rgb.html
RISCHIO = 80% NESSUN RISCHIO (GREEN) rgb(0,255,0)
80% < RISCHIO = 100% ATTENZIONE (YELLOW) rgb(255,255,0)
100% < RISCHIO = 120% ALLARME (ORANGE) rgb(255,165,0)
RISCHIO > 120% EMERGENZA (RED) rgb(255,0,0)
*/
/**
* Given t air temperature (Celsius), rh relative humidity (%) and tg globometric temperature gives heat risk level by using WBGT index.
* @param {number} t,rh,tg,tresh
* @return {number}
* @customfunction
*/
function heat_risk_text_level(t,rh,tg,tresh) {
if ( tresh === undefined) {tresh = 28.5 };
var risk =wbgt(t,rh,tg)/tresh;
var level_list=["NESSUN RISCHIO","ATTENZIONE","ALLARME","EMERGENZA"];
var class;
if ( risk <= 0.8) { class = 1;}
else if (risk <= 1) { class = 2;}
else if (risk <= 1.2) { class = 3;}
else { class = 4};
return(level_list[class]);
}
/**
* Given t air temperature (Celsius), rh relative humidity (%) and tg globometric temperature gives heat risk color level by using WBGT index.
* @param {number} t,rh,tg,tresh
* @return {number}
* @customfunction
*/
function heat_risk_color_level(wbgt,tresh) {
if ( tresh === undefined) {tresh = 28.5 };
var risk =wbgt/tresh;
var colorcode_list=["green","yellow","orange","red"];
var class;
if ( risk <= 0.8) { class = 1;}
else if (risk <= 1) { class = 2;}
else if (risk <= 1.2) { class = 3;}
else { class = 4};
return(colorcode_list[class]);
}
/**
* Gives heat color risk level in hex format by using WBGT index.
* @param {number} t,rh,tg,tresh
* @return {number}
* @customfunction
*/
function heat_risk_hexrgb_level(wbgt,tresh) {
if ( tresh === undefined) {tresh = 28.5 };
var risk= wbgt/tresh;
var colorcode_hex=["#00ff00","#ffff00","#ffa500","#ff0000"];
var class;
if ( risk <= 0.8) { class = 1;}
else if (risk <= 1) { class = 2;}
else if (risk <= 1.2) { class = 3;}
else { class = 4};
return(colorcode_hex[class]);
}
/**
* Given t air temperature (Celsius), rh relative humidity (%) and tg globometric temperature gives heat risk value by using WBGT index.
* @param {number} t,rh,tg,tresh
* @return {number}
* @customfunction
*/
function heat_risk_index_level(wbgt,tresh) {
if ( tresh === undefined) {tresh = 28.5 };
var risk =(wbgt/tresh)*100;
return(risk);
}
/**
* Given a air temperature t (Celsius degree), globe temperature (Celsius degree), wind speed in m per second and diam of glbe in meters compute the mean radiant temperature (Celsius degree) ;
* Reference; BSL, page 23.;
*
* @param {number} t, tg, wind, diam_glob
* @return {number}
* @customfunction
*/
function mrt_globe(t, tg, wind, diam_globe,emis_globe)
{ if ( diam_globe === undefined) {diam_globe = 0.05;} ; // in meters.
if ( emis_globe === undefined) {emis_globe = 0.97;} ;
var stefanb = 0.0000000567;
return Math.pow((Math.pow(tg + 273.15, 4) + ((1.1 * Math.pow(10,8) * Math.pow(wind,0.6)) /(emis_globe*Math.pow(diam_globe,0.4))) * (tg - t)), 0.25)- 273.15;
}
/**
* Given air temperature t (Celsius degrees) and relative humidity (%) gives deficit of saturarion in hPa.
*
* @param {number} t,rh
* @return {number}
* @customfunction
*/
function deficitsat(t,rh)
{
var pws = PVS(t);
var pae=rh/100*pws;
return (pws-pae);
}
/**
* Given air temperature t (Celsius degrees) and relative humidity (%) gives saturation vapor pressure in hPa.
*
* @param {number} t,rh
* @return {number}
* @customfunction
*/
function es(t)
{
var es_air, tk,i;
g = new Array(8);
g[0] =-2.8365744E3;
g[1] =-6.028076559E3;
g[2] =1.954263612E1;
g[3] =-2.737830188E-2;
g[4] =1.6261698E-5;
g[5] =7.0229056E-10;
g[6] =-1.8680009E-13;
g[7] = 2.7150305;
tk = t+273.15;
es_air = g[7]*Math.log(tk);
for ( var i = 0; i < 7; i ++ )
{es_air = es_air+g[i]*Math.pow(tk,i-2)};
var es = Math.exp(es_air)*0.01;
return es;
}
/**
* Given air temperature (Celsius), relative humidity (%) and pressure ( pa) gives saturation vapor pressure in hPa.
*
* @param {number} t,rh
* @return {number}
* @customfunction
*/
function pvap(t,rh)
{
var E;
E = es(t) * (rh/100);
return E;
}
/**
* Given air temperature t (Celsius) function gives Saturation Vapor Pressure in hectoPascal (hPa)
*
* @param {number} t
* @return {number}
* @customfunction
*/
function PVS(t)
{
t=t+273.16;
var minT = 173; // -100 Deg. C.
var maxT = 678;
var noresp = -9999;
var pvsresp;
if (t < minT || t> maxT)
{return noresp;}
else if (t<273.15)
{pvsresp=pvsIce(t)/100;}
else
{pvsresp=pvsWater(t)/100;}
return(TwoDec(pvsresp));
}
/*
* Saturation Vapor Pressure formula for range -100..0 Deg. C.
* Hardy B, 1998,"ITS-90 Formulations for Vapor Pressure, Frostpoint Temperature,Dewpoint Temperature, and Enhancement Factors in the Range 100 to +100 C".
* Proceedings of the Third International Symposium on Humidity & Moisture",Teddington, London, England, April 1998
*/
function pvsIce(T)
{
var k0 = -5.8666426e3;
var k1 = 2.232870244e1;
var k2 = 1.39387003e-2;
var k3 = -3.4262402e-5;
var k4 = 2.7040955e-8;
var k5 = 6.7063522e-1;
var lnP = k0/T + k1 + (k2 + (k3 + (k4*T))*T)*T + k5*Math.log(T);
return Math.exp(lnP);
}
/**
* Given air temperature T (Celsius) function gives Saturation Vapor Pressure. Dimension of outcomes in kPascal (kPa)
* Saturation Vapor Pressure formula for range 273..678 Deg. K.
* Equation (30) in Section 8.1 "The Saturation-Pressure Equation (Basic Equation),Erlangen, Germany, September 1997.
* @param {number} T
* @return {number}
* @customfunction
*/
function pvsWater(T)
{
var n1 = 0.11670521452767e4;
var n6 = 0.14915108613530e2;
var n2 = -0.72421316703206e6;
var n7 = -0.48232657361591e4;
var n3 = -0.17073846940092e2;
var n8 = 0.40511340542057e6;
var n4 = 0.12020824702470e5;
var n9 = -0.23855557567849;
var n5 = -0.32325550322333e7;
var n10 = 0.65017534844798e3;
var th = T+n9/(T-n10);
var A = (th+n1)*th+n2;
var B = (n3*th+n4)*th+n5;
var C = (n6*th+n7)*th+n8;
var p = 2*C/(-B+Math.sqrt(B*B-4*A*C));
p *= p;
p *= p;
return p*1e6;
}
/**
* Given a air temperature t (Celsius) and air relative humidity (RH) and formula gives the Dew Point in Celsius degrees.
*
* @param {number} t,RH
* @return {number}
* @customfunction
*/
function dewpoint(t,rh,formula) {
if (formula == undefined) {formula == "NOAA"}
var RHD = rh / 100;
if (formula == "Paroscientific")
{return 237.3 * (Math.log(RHD) / 17.27 + t / (237.3 + t)) / (1 - Math.log(RHD) / 17.27 - t / (237.3 + t));}
else if (formula == "Sonntag")
{return 243.12 * (Math.log(RHD) / 17.62 + t / (243.12 + t)) / (1 - Math.log(RHD) / 17.62 - t / (243.12 + t));}
else
{return 243.5 * (Math.log(RHD) / 17.67 + t / (243.5 + t)) / (1 - Math.log(RHD) / 17.67 - t / (243.5 + t));}
}
function h_evap(tk){
var evap = (313.15 - tk) / 30 * (-71100) + 2407300;
return(evap);
}
function emis_atm(t,rh)
{
// Reference; Oke (2nd edition), page 373.;
var e = (rh*0.01) * esat(t+273.15);
return 0.575 * Math.pow(e,0.143);
}
/**
* Calculate the saturation vapor pressure (mb) over liquid water given the temperature (K).;
* Reference Buck's (1981) approximation (eqn 3) of Wexler's (1976) formulae.;
*
* @param {number} t
* @return {number}
* @customfunction
*/
function esat (tk)
{
var esat= 6.1121 * Math.exp(17.502 * (tk - 273.15) / (tk - 32.18));
return 1.004 * esat;
}
/**
* Given a air temperature t (Celsius) compute the viscosity of air, kg/(m s) given temperature, K;
* Reference; BSL, page 23.;
*
* @param {number} t
* @return {number}
* @customfunction
*/
function viscosity(t)
{
var omega = ((t+273.15) / 97 - 2.9) / 0.4 * (-0.034) + 1.048;
return 0.0000026693 * Math.pow((28.97 * t),0.5) / (Math.pow(3.617,2) * omega);
}
/**
* Given a air temperature tk (Kelvin) and air pressure in millibar(hPA) compute the heat diffusivity in air;
* Reference; BSL, page 23.;
*
* @param {number} t
* @return {number}
* @customfunction
*/
function diffusivity(tk, pair) {
if( pair === undefined ) {pair=1013.25;};
var pcrit13 = Math.pow(36.4 * 218, 1 / 3);
var tcrit512 = Math.pow(132 * 647.3,(5 / 12));
var Tcrit12 = Math.pow((132 * 647.3), 0.5);
var Mmix = Math.pow((1 / 28.97 + 1 / 18.015),0.5);
var diffusivity = 0.000364 * Math.pow((tk / Tcrit12),2.334) * pcrit13 * tcrit512 * Mmix / (pair / 1013.25) * 0.0001;
return(diffusivity);
}
// Purpose: to calculate the convective h the heat tranfer coefficient for flow around a sphere.; t in Kelvin degrees;
// Reference : Bird, Stewart, && Lightfoot (BSL), page 409.;
function h_sphere_in_air(tk,speed,pair,speedmin,diam_globe)
{
if( diam_globe === undefined ) {diam_globe=0.0508;};
if( speedmin === undefined ) {speedmin=0.1;};
if( pair === undefined ) {pair=1013.25;};
var Rair = 8314.34 / 28.97;
var Pr = 1003.5 / (1003.5 + 1.25 * Rair);
var thermal_con = (1003.5 + 1.25 * 8314.34 / 28.97) * viscosity(tk-273.15);
var density = (pair * 100) / (Rair * tk); // kg/m3;
if(speed < speedmin ){speed = speedmin};
var Re = (speed * density * diam_globe)/ viscosity(tk-273.15);
var Nu = 2 + 0.6 * Math.pow(Re,0.5) * Math.pow(Pr, 0.3333);
return (Nu * thermal_con) / diam_globe; // W/(m2 K);
}
// Purpose: to calculate the convective heat tranfer coefficient for heat flow around a cylinder; t in K
// Reference : Bird, Stewart, && Lightfoot (BSL), page 409.;
function h_cylinder_in_air(tk,speed,pair,speedmin,diam_wick){
if( diam_wick === undefined ) { diam_wick=0.007;};
if( speedmin === undefined ) { speedmin=0.1;};
if( pair === undefined ) { pair=1013.25;};
var m_air = 28.97;
var r_gas = 8314.34;
var r_air = r_gas / m_air;
var cp = 1003.5; // heat capaticy at constant pressure of dry air
var Pr = cp / (cp + (1.25 * r_air));
// Calculate the thermal conductivity of air, W/(m K)
var thermal_con = (cp + 1.25 * r_gas / m_air) * viscosity(tk-273.15);
// Density of the air
var density = pair * 100 / (r_air * tk);
if (speed < speedmin) {speed = speedmin};
// Reynolds number
var Re = speed * density * diam_wick / viscosity(tk);
// Nusselt number
var Nu = 0.281 * Math.pow(Re,0.6) * Math.pow(Pr, 0.44);
// Convective heat transfer coefficient in W/(m2 K) for a long cylinder in cross flow
var h_cylinder_in_air = Nu * thermal_con / diam_wick ;
return(h_cylinder_in_air);
}
/**
* Given t, air temperature (Celsius degree)
* rh, relative humidity (%) Used only this way to input humidity level
* wind, relative air velocity (m/s)
* trad, mean radiant temperature (Celsius degree)
* M, metabolic rate (met)
* W, external work, normally around 0 (met)
* clo, clothing (clo) compute the Predicted Mean Vote for indoor moderate thermal environments.
* Reference; ISO7730
*
* @param {number} t,rh,wind,trad,M,W,clo
* @return {number}
* @customfunction
*/
function PMV_ISO7730(t,rh,wind,trad,M,W,clo)
{
var pa, icl, mw, fcl, hcf, taa, tra, tcla, p1, p2, p3, p4,
p5, xn, xf, eps, hcn, hc, tcl, hl1, hl2, hl3, hl4, hl5, hl6,
ts, pmv, ppd, n;
pa = rh * 10 * Math.exp(16.6536 - 4030.183 / (t + 235));
icl = 0.155 * clo; //thermal insulation of the clothing in M2K/W
mw = M - W; //internal heat production in the human body
if (icl <= 0.078) {fcl = 1 + (1.29 * icl)} else {fcl = 1.05 + (0.645 * icl)};
// heat transf. coeff. by forced convection
hcf = 12.1 * Math.sqrt(wind);
taa = t + 273;
tra = trad + 273;
tcla = taa + (35.5 - t) / (3.5 * icl + 0.1);
p1 = icl * fcl;
p2 = p1 * 3.96;
p3 = p1 * 100;
p4 = p1 * taa;
p5 = 308.7 - 0.028 * mw + p2 * Math.pow(tra / 100, 4);
xn = tcla / 100;
xf = tcla / 50;
eps = 0.00015;
n = 0;
while (Math.abs(xn - xf) > eps) {
xf = (xf + xn) / 2;
hcn = 2.38 * Math.pow(Math.abs(100.0 * xf - taa), 0.25);
if (hcf > hcn) hc = hcf;
else hc = hcn;
xn = (p5 + p4 * hc - p2 * Math.pow(xf, 4)) / (100 + p3 * hc);
++n;
if (n > 150) {
return('Max iterations exceeded');
}
}
tcl = 100 * xn - 273;
// heat loss diff. through skin
hl1 = 3.05 * 0.001 * (5733 - (6.99 * mw) - pa);
// heat loss by sweating
if (mw > 58.15) {hl2 = 0.42 * (mw - 58.15)} else {hl2 = 0};
// latent respiration heat loss
hl3 = 1.7 * 0.00001 * M * (5867 - pa);
// dry respiration heat loss
hl4 = 0.0014 * M * (34 - t);
// heat loss by radiation R
hl5 = 3.96 * fcl * (Math.pow(xn, 4) - Math.pow(tra / 100, 4));
// heat loss by convection C
hl6 = fcl * hc * (tcl - t);
ts = 0.303 * Math.exp(-0.036 * M) + 0.028;
pmv = ts * (mw - hl1 - hl2 - hl3 - hl4 - hl5 - hl6);
return(TwoDec(pmv));
}
/**
* Given air temperature (Celsius), relative humidity (%), wind velocity (m/sec), direct beam short-wavelength radiation (W/mq)
* and sunelev sun elevation angle (degrees) gives Steadman outdoor index.
*
* @param {number} t, rh, wind, rshort, sunelev
* @return {number}
* @customfunction
*/
function steadman_outdoor(t,rh,wind,rad,sunelev)
{
var steadman_outdoor_sun=9999;
if (rh > 100.1 || rh < 0.0) {return 9999};
if (t > 100.0 || t < -100.0) {return 9999};
if (rad < 50) { return(steadman_outdoor_shade(t,rh,wind))}
else {
var ee = (rh/1000.0)*(6.105*Math.exp((t*17.27)/(237.7+t)));
var q_glob = 0.56*(0.386-(0.0032*sunelev))*rad + 0.224*(0.1*rad)+ 0.028*rad- 150.0*(0.38-0.16*(Math.pow(ee,0.5)));
if (q_glob > 0.0) {steadman_outdoor_sun = t+3.48*(ee)-0.7*wind +0.7*q_glob/(wind+10.0)-4.25};
};
return TwoDec(steadman_outdoor_sun);
}
/**
* Given air temperature (Celsius), rh relative humidity (%), wind velocity (m/sec) gives Steadman outdoor shade index.
*
* @param {number} t,rh,wind
* @return {number}
* @customfunction
*/
function steadman_outdoor_shade(t,rh,wind)
{
var steadman_outdoor_shade,e;
steadman_outdoor_shade = 9999;
if (rh > 100.1 || rh < 0.0)
{return 9999}
else if (wind > 130.0 || wind < 0.0)
{return 9999}
else if (t > 100.0 || t < -100.0)
{return 9999}
else
{
e = (rh/100.0)*(6.105*Math.exp((t*17.27)/(237.7+t)));
steadman_outdoor_shade = t+(0.33*e)-(0.7*wind)-4.0;
};
return TwoDec(steadman_outdoor_shade);
}
/**
* Given air temperature (Celsius), relative humidity (%) gives Steadman indoor index.
*
* @param {number} t,rh
* @return {number}
* @customfunction
*/
function steadman_indoor(t,rh)
{
var steadman_indoor,e;
steadman_indoor = -9999;
if (rh > 100.1 || rh < 0.0) {return -9999}
else if (t > 100.0 || t < -100.0) {return -9999}
else
{
e = (rh/100.0)*(6.105*Math.exp((t*17.27)/(237.7+t)));
steadman_indoor = -2.56+(0.89*t)+(0.382*e);
};
return TwoDec(steadman_indoor);
}
/**
* Given air temperature (Celsius Degrees), air relative humidity (%), wind velocity ( mpersec) and mean radiant temperature (trad) in Celsius Degrees
* gives 10 UTCI class.
* @param {number} t,rh,wind,trad
* @return {number}
* @customfunction
*/
function utci_class10(t,rh,wind,trad)
{
var utci_v,utci_c;
utci_v=UTCI(t,rh,wind,trad);
if ( utci_v > 46.0)
{utci_c=10.0;}
else if (utci_v > 38.0 && utci_v <= 46.0)
{utci_c=9.0;}
else if (utci_v > 32.0 && utci_v <= 38.0)
{utci_c=8.0;}
else if (utci_v > 26.0 && utci_v <= 32.0)
{utci_c=7.0;}
else if (utci_v > 9.0 && utci_v <= 26.0)
{utci_c=6.0;}
else if (utci_v > 0.0 && utci_v <= 9.0)
{utci_c=5.0;}
else if (utci_v > -13.0 && utci_v <= 0)
{utci_c=4.0;}
else if (utci_v > -27.0 && utci_v <= -13.0)
{utci_c=3.0;}
else if (utci_v > -40.0 && utci_v <= -27.0)
{utci_c=2.0;}
else if (utci_v <= -40.0)
{utci_c=1.0;}
else if (utci_v == 9999)
{utci_c=9999};
return utci_c;
}
/**
* Given air temperature (Celsius Degrees), air relative humidity (%), wind velocity ( mpersec) and mean radiant temperature (trad) in Celsius Degrees
* gives 7 UTCI class.
* @param {number} t,rh,wind,trad
* @return {number}
* @customfunction
*/
function utci_class7(t,rh, wind,trad)
{
var utci_v,utci_c;
utci_v=UTCI(t,rh,wind,trad);
if ( utci_v > 46.0)
{utci_c=7.0;}
else if (utci_v > 38.0 && utci_v <= 46.0)
{utci_c=7.0;}
else if (utci_v > 32.0 && utci_v <= 38.0)
{utci_c=7.0;}
else if (utci_v > 26.0 && utci_v <= 32.0)
{utci_c=6.0;}
else if (utci_v > 16.0 && utci_v <= 26.0)
{utci_c=6.0;}
else if (utci_v > 0.0 && utci_v <= 16.0)
{utci_c=5.0;}
else if (utci_v > -13.0 && utci_v <= 0)
{utci_c=4.0;}
else if (utci_v > -27.0 && utci_v <= -13.0)
{utci_c=3.0;}
else if (utci_v > -40.0 && utci_v <= -27.0)
{utci_c=2.0;}
else if (utci_v <= -40.0)
{utci_c=1.0;}
else if (utci_v == 9999)
{utci_c=9999};
return utci_c;
}
/**
* Given air temperature (t Celsius degrees), relative humidity (%), wind velocity (m/sec)
* and mean radiant temperature (tmrt in Celsius degree)
* gives the Universal Thermal Climate Index (Celsius degree).
*
* @param {number} t, rh, wind, tmrt
* @return {number}
* @customfunction
*/
function UTCI(ta,rh,wind,tmrt)
{
var ta,pa,va, e, dtm,i;
e = es(ta)/10; // use vapour pressure in kPa
pa = (e*rh/100.0);
va = wind;
if ( va < 0.51) {va=0.5;};
if ( va > 17) {va=17;};
dtm = tmrt - ta;
utci = new Array(210);
utci[0]=ta;
utci[1]=6.07562052E-01;
utci[2]=-2.27712343E-02*ta;
utci[3]=8.06470249E-04*ta*ta;
utci[4]=-1.54271372E-04*ta*ta*ta;
utci[5]=-3.24651735E-06*ta*ta*ta*ta;
utci[6]=7.32602852E-08*ta*ta*ta*ta*ta;
utci[7]=1.35959073E-09*ta*ta*ta*ta*ta*ta;
utci[8]=-2.25836520E-00*va;
utci[9]=8.80326035E-02*ta*va;
utci[10]=2.16844454E-03*ta*ta*va;
utci[11]=-1.53347087E-05*ta*ta*ta*va;
utci[12]=-5.72983704E-07*ta*ta*ta*ta*va;
utci[13]=-2.55090145E-09*ta*ta*ta*ta*ta*va;
utci[14]=-7.51269505E-01*va*va;
utci[15]=-4.08350271E-03*ta*va*va;
utci[16]=-5.21670675E-05*ta*ta*va*va;
utci[17]=1.94544667E-06*ta*ta*ta*va*va;
utci[18]=1.14099531E-08*ta*ta*ta*ta*va*va;
utci[19]=1.58137256E-01*va*va*va;
utci[20]=-6.57263143E-05*ta*va*va*va;
utci[21]=2.22697524E-07*ta*ta*va*va*va;
utci[22]=-4.16117031E-08*ta*ta*ta*va*va*va;
utci[23]=-1.27762753E-02*va*va*va*va;
utci[24]=9.66891875E-06*ta*va*va*va*va;
utci[25]=2.52785852E-09*ta*ta*va*va*va*va;
utci[26]=4.56306672E-04*va*va*va*va*va;
utci[27]=-1.74202546E-07*ta*va*va*va*va*va;
utci[28]=-5.91491269E-06*va*va*va*va*va*va;
utci[29]=3.98374029E-01*dtm;
utci[30]=1.83945314E-04*ta*dtm;
utci[31]=-1.73754510E-04*ta*ta*dtm;
utci[32]=-7.60781159E-07*ta*ta*ta*dtm;
utci[33]=3.77830287E-08*ta*ta*ta*ta*dtm;
utci[34]=5.43079673E-10*ta*ta*ta*ta*ta*dtm;
utci[35]=-2.00518269E-02*va*dtm;
utci[36]=8.92859837E-04*ta*va*dtm;
utci[37]=3.45433048E-06*ta*ta*va*dtm;
utci[38]=-3.77925774E-07*ta*ta*ta*va*dtm;
utci[39]=-1.69699377E-09*ta*ta*ta*ta*va*dtm;
utci[40]=1.69992415E-04*va*va*dtm;
utci[41]=-4.99204314E-05*ta*va*va*dtm;
utci[42]=2.47417178E-07*ta*ta*va*va*dtm;
utci[43]=1.07596466E-08*ta*ta*ta*va*va*dtm;
utci[44]=8.49242932E-05*va*va*va*dtm;
utci[45]=1.35191328E-06*ta*va*va*va*dtm;
utci[46]=-6.21531254E-09*ta*ta*va*va*va*dtm;
utci[47]=-4.99410301E-06*va*va*va*va*dtm;
utci[48]=-1.89489258E-08*ta*va*va*va*va*dtm;
utci[49]=8.15300114E-08*va*va*va*va*va*dtm;
utci[50]=7.55043090E-04*dtm*dtm;
utci[51]=-5.65095215E-05*ta*dtm*dtm;
utci[52]=-4.52166564E-07*ta*ta*dtm*dtm;
utci[53]=2.46688878E-08*ta*ta*ta*dtm*dtm;
utci[54]=2.42674348E-10*ta*ta*ta*ta*dtm*dtm;
utci[55]=1.54547250E-04*va*dtm*dtm;
utci[56]=5.24110970E-06*ta*va*dtm*dtm;
utci[57]=-8.75874982E-08*ta*ta*va*dtm*dtm;
utci[58]=-1.50743064E-09*ta*ta*ta*va*dtm*dtm;
utci[59]=-1.56236307E-05*va*va*dtm*dtm;
utci[60]=-1.33895614E-07*ta*va*va*dtm*dtm;
utci[61]=2.49709824E-09*ta*ta*va*va*dtm*dtm;
utci[62]=6.51711721E-07*va*va*va*dtm*dtm;
utci[63]=1.94960053E-09*ta*va*va*va*dtm*dtm;
utci[64]=-1.00361113E-08*va*va*va*va*dtm*dtm;
utci[65]=-1.21206673E-05*dtm*dtm*dtm;
utci[66]=-2.18203660E-07*ta*dtm*dtm*dtm;
utci[67]=7.51269482E-09*ta*ta*dtm*dtm*dtm;
utci[68]=9.79063848E-11*ta*ta*ta*dtm*dtm*dtm;
utci[69]=1.25006734E-06*va*dtm*dtm*dtm;
utci[70]=-1.81584736E-09*ta*va*dtm*dtm*dtm;
utci[71]=-3.52197671E-10*ta*ta*va*dtm*dtm*dtm;
utci[72]=-3.36514630E-08*va*va*dtm*dtm*dtm;
utci[73]=1.35908359E-10*ta*va*va*dtm*dtm*dtm;
utci[74]=4.17032620E-10*va*va*va*dtm*dtm*dtm;
utci[75]=-1.30369025E-09*dtm*dtm*dtm*dtm;
utci[76]=4.13908461E-10*ta*dtm*dtm*dtm*dtm;
utci[77]=9.22652254E-12*ta*ta*dtm*dtm*dtm*dtm;
utci[78]=-5.08220384E-09*va*dtm*dtm*dtm*dtm;
utci[79]=-2.24730961E-11*ta*va*dtm*dtm*dtm*dtm;
utci[80]=1.17139133E-10*va*va*dtm*dtm*dtm*dtm;
utci[81]=6.62154879E-10*dtm*dtm*dtm*dtm*dtm;
utci[82]=4.03863260E-13*ta*dtm*dtm*dtm*dtm*dtm;
utci[83]=1.95087203E-12*va*dtm*dtm*dtm*dtm*dtm;
utci[84]=-4.73602469E-12*dtm*dtm*dtm*dtm*dtm*dtm;
utci[85]=5.12733497E-00*pa;
utci[86]=-3.12788561E-01*ta*pa;
utci[87]=-1.96701861E-02*ta*ta*pa;
utci[88]=9.99690870E-04*ta*ta*ta*pa;
utci[89]=9.51738512E-06*ta*ta*ta*ta*pa;
utci[90]=-4.66426341E-07*ta*ta*ta*ta*ta*pa;
utci[91]=5.48050612E-01*va*pa;
utci[92]=-3.30552823E-03*ta*va*pa;
utci[93]=-1.64119440E-03*ta*ta*va*pa;
utci[94]=-5.16670694E-06*ta*ta*ta*va*pa;
utci[95]=9.52692432E-07*ta*ta*ta*ta*va*pa;
utci[96]=-4.29223622E-02*va*va*pa;
utci[97]=5.00845667E-03*ta*va*va*pa;
utci[98]=1.00601257E-06*ta*ta*va*va*pa;
utci[99]=-1.81748644E-06*ta*ta*ta*va*va*pa;
utci[100]=-1.25813502E-03*va*va*va*pa;
utci[101]=-1.79330391E-04*ta*va*va*va*pa;
utci[102]=2.34994441E-06*ta*ta*va*va*va*pa;
utci[103]=1.29735808E-04*va*va*va*va*pa;
utci[104]=1.29064870E-06*ta*va*va*va*va*pa;
utci[105]=-2.28558686E-06*va*va*va*va*va*pa;
utci[106]=-3.69476348E-02*dtm*pa;
utci[107]=1.62325322E-03*ta*dtm*pa;
utci[108]=-3.14279680E-05*ta*ta*dtm*pa;
utci[109]=2.59835559E-06*ta*ta*ta*dtm*pa;
utci[110]=-4.77136523E-08*ta*ta*ta*ta*dtm*pa;
utci[111]=8.64203390E-03*va*dtm*pa;
utci[112]=-6.87405181E-04*ta*va*dtm*pa;
utci[113]=-9.13863872E-06*ta*ta*va*dtm*pa;
utci[114]=5.15916806E-07*ta*ta*ta*va*dtm*pa;
utci[115]=-3.59217476E-05*va*va*dtm*pa;
utci[116]=3.28696511E-05*ta*va*va*dtm*pa;
utci[117]=-7.10542454E-07*ta*ta*va*va*dtm*pa;
utci[118]=-1.24382300E-05*va*va*va*dtm*pa;
utci[119]=-7.38584400E-09*ta*va*va*va*dtm*pa;
utci[120]=2.20609296E-07*va*va*va*va*dtm*pa;
utci[121]=-7.32469180E-04*dtm*dtm*pa;
utci[122]=-1.87381964E-05*ta*dtm*dtm*pa;
utci[123]=4.80925239E-06*ta*ta*dtm*dtm*pa;
utci[124]=-8.75492040E-08*ta*ta*ta*dtm*dtm*pa;
utci[125]=2.77862930E-05*va*dtm*dtm*pa;
utci[126]=-5.06004592E-06*ta*va*dtm*dtm*pa;
utci[127]=1.14325367E-07*ta*ta*va*dtm*dtm*pa;
utci[128]=2.53016723E-06*va*va*dtm*dtm*pa;
utci[129]=-1.72857035E-08*ta*va*va*dtm*dtm*pa;
utci[130]=-3.95079398E-08*va*va*va*dtm*dtm*pa;
utci[131]=-3.59413173E-07*dtm*dtm*dtm*pa;
utci[132]=7.04388046E-07*ta*dtm*dtm*dtm*pa;
utci[133]=-1.89309167E-08*ta*ta*dtm*dtm*dtm*pa;
utci[134]=-4.79768731E-07*va*dtm*dtm*dtm*pa;
utci[135]=7.96079978E-09*ta*va*dtm*dtm*dtm*pa;
utci[136]=1.62897058E-09*va*va*dtm*dtm*dtm*pa;
utci[137]=3.94367674E-08*dtm*dtm*dtm*dtm*pa;
utci[138]=-1.18566247E-09*ta*dtm*dtm*dtm*dtm*pa;
utci[139]=3.34678041E-10*va*dtm*dtm*dtm*dtm*pa;
utci[140]=-1.15606447E-10*dtm*dtm*dtm*dtm*dtm*pa;
utci[141]=-2.80626406E-00*pa*pa;
utci[142]=5.48712484E-01*ta*pa*pa;
utci[143]=-3.99428410E-03*ta*ta*pa*pa;
utci[144]=-9.54009191E-04*ta*ta*ta*pa*pa;
utci[145]=1.93090978E-05*ta*ta*ta*ta*pa*pa;
utci[146]=-3.08806365E-01*va*pa*pa;
utci[147]=1.16952364E-02*ta*va*pa*pa;
utci[148]=4.95271903E-04*ta*ta*va*pa*pa;
utci[149]=-1.90710882E-05*ta*ta*ta*va*pa*pa;
utci[150]=2.10787756E-03*va*va*pa*pa;
utci[151]=-6.98445738E-04*ta*va*va*pa*pa;
utci[152]=2.30109073E-05*ta*ta*va*va*pa*pa;
utci[153]=4.17856590E-04*va*va*va*pa*pa;
utci[154]=-1.27043871E-05*ta*va*va*va*pa*pa;
utci[155]=-3.04620472E-06*va*va*va*va*pa*pa;
utci[156]=5.14507424E-02*dtm*pa*pa;
utci[157]=-4.32510997E-03*ta*dtm*pa*pa;
utci[158]=8.99281156E-05*ta*ta*dtm*pa*pa;
utci[159]=-7.14663943E-07*ta*ta*ta*dtm*pa*pa;
utci[160]=-2.66016305E-04*va*dtm*pa*pa;
utci[161]=2.63789586E-04*ta*va*dtm*pa*pa;
utci[162]=-7.01199003E-06*ta*ta*va*dtm*pa*pa;
utci[163]=-1.06823306E-04*va*va*dtm*pa*pa;
utci[164]=3.61341136E-06*ta*va*va*dtm*pa*pa;
utci[165]=2.29748967E-07*va*va*va*dtm*pa*pa;
utci[166]=3.04788893E-04*dtm*dtm*pa*pa;
utci[167]=-6.42070836E-05*ta*dtm*dtm*pa*pa;
utci[168]=1.16257971E-06*ta*ta*dtm*dtm*pa*pa;
utci[169]=7.68023384E-06*va*dtm*dtm*pa*pa;
utci[170]=-5.47446896E-07*ta*va*dtm*dtm*pa*pa;
utci[171]=-3.59937910E-08*va*va*dtm*dtm*pa*pa;
utci[172]=-4.36497725E-06*dtm*dtm*dtm*pa*pa;
utci[173]=1.68737969E-07*ta*dtm*dtm*dtm*pa*pa;
utci[174]=2.67489271E-08*va*dtm*dtm*dtm*pa*pa;
utci[175]=3.23926897E-09*dtm*dtm*dtm*dtm*pa*pa;
utci[176]=-3.53874123E-02*pa*pa*pa;
utci[177]=-2.21201190E-01*ta*pa*pa*pa;
utci[178]=1.55126038E-02*ta*ta*pa*pa*pa;
utci[179]=-2.63917279E-04*ta*ta*ta*pa*pa*pa;
utci[180]=4.53433455E-02*va*pa*pa*pa;
utci[181]=-4.32943862E-03*ta*va*pa*pa*pa;
utci[182]=1.45389826E-04*ta*ta*va*pa*pa*pa;
utci[183]=2.17508610E-04*va*va*pa*pa*pa;
utci[184]=-6.66724702E-05*ta*va*va*pa*pa*pa;
utci[185]=3.33217140E-05*va*va*va*pa*pa*pa;
utci[186]=-2.26921615E-03*dtm*pa*pa*pa;
utci[187]=3.80261982E-04*ta*dtm*pa*pa*pa;
utci[188]=-5.45314314E-09*ta*ta*dtm*pa*pa*pa;
utci[189]=-7.96355448E-04*va*dtm*pa*pa*pa;
utci[190]=2.53458034E-05*ta*va*dtm*pa*pa*pa;
utci[191]=-6.31223658E-06*va*va*dtm*pa*pa*pa;
utci[192]=3.02122035E-04*dtm*dtm*pa*pa*pa;
utci[193]=-4.77403547E-06*ta*dtm*dtm*pa*pa*pa;
utci[194]=1.73825715E-06*va*dtm*dtm*pa*pa*pa;
utci[195]=-4.09087898E-07*dtm*dtm*dtm*pa*pa*pa;
utci[196]=6.14155345E-01*pa*pa*pa*pa;
utci[197]=-6.16755931E-02*ta*pa*pa*pa*pa;
utci[198]=1.33374846E-03*ta*ta*pa*pa*pa*pa;
utci[199]=3.55375387E-03*va*pa*pa*pa*pa;
utci[200]=-5.13027851E-04*ta*va*pa*pa*pa*pa;
utci[201]=1.02449757E-04*va*va*pa*pa*pa*pa;
utci[202]=-1.48526421E-03*dtm*pa*pa*pa*pa;
utci[203]=-4.11469183E-05*ta*dtm*pa*pa*pa*pa;
utci[204]=-6.80434415E-06*va*dtm*pa*pa*pa*pa;
utci[205]=-9.77675906E-06*dtm*dtm*pa*pa*pa*pa;
utci[206]=8.82773108E-02*pa*pa*pa*pa*pa;
utci[207]=-3.01859306E-03*ta*pa*pa*pa*pa*pa;
utci[208]=1.04452989E-03*va*pa*pa*pa*pa*pa;
utci[209]=2.47090539E-04*dtm*pa*pa*pa*pa*pa;
utci[210]=1.48348065E-03*pa*pa*pa*pa*pa*pa;
var total = 0;
for (i = 0, n = utci.length; i < n; ++i)
{
total = total+utci[i];
};
if ( ta < -50.0 || ta > 60.0 ) {total=9999};
if ( (tmrt < ta-30.0) || (tmrt > ta+70.0 )) {total=9999};
if ( rh <= 0.0 || rh >= 100.0 ) {total=9999};
return TwoDec(total);
}
/**
* Given air temperature (Celsius), relative humidity (%) give a heat index in Celsius degree.
* References: http://www.wpc.ncep.noaa.gov/html/heatindex.shtml [2] https://en.wikipedia.org/wiki/Heat_index [3] http://www.srh.noaa.gov/images/ffc/pdf/ta_htindx.PDF
*
* @param {number} t,rh
* @return {number}
* @customfunction
*/
function heatindex(t,rh)
{
var tf, tf2, ur2, hif;
tf = C2F(t);
tf2 = Math.pow(tf, 2.0);
ur2 = Math.pow(rh, 2.0);
hif = -42.379 + 2.04901523 * tf + 10.1433127 * rh - 0.22475541 *tf * rh
- 6.83783 * 0.001* tf2 - 5.481717 * 0.01* ur2 +1.22874 * 0.001* tf2* rh
+ 8.5282 * 0.0001* tf * ur2 -1.99 * 0.000001* tf2* ur2;
if (t < 44 & t > 27 & rh < 13 )
{
return (TwoDec(F2C(hif-((13-rh)/4)*Math.sqrt((17-Mat.abs(tf-95.))/17))));
}
else if (t < 31 & t > 27 & rh > 85 ) {
return (TwoDec(F2C(hif-((rh-85)/10) * ((87-tf)/5))));
}
if (t > 27)
{
return (TwoDec(F2C(hif)));
}
else
{return(TwoDec(F2C(0.5 * (tf + 61.0 + ((tf-68.0) *1.2) + (rh*0.094) ))))};
}
///////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////