Commit 7cd89161 authored by Uwe Schulzweida's avatar Uwe Schulzweida
Browse files

ECA update

parent 37971e07
......@@ -163,7 +163,7 @@ doc/tex/mod/Vardup -text
doc/tex/mod/Vargen -text
doc/tex/mod/Vertint -text
doc/tex/mod/Vertstat -text
doc/tex/mod/Wcl -text
doc/tex/mod/Wct -text
doc/tex/mod/Wind -text
doc/tex/mod/Writegrid -text
doc/tex/mod/Writerandom -text
......@@ -280,7 +280,6 @@ src/Splitrec.c -text
src/Splittime.c -text
src/Splityear.c -text
src/Subtrend.c -text
src/Tchill.c -text
src/Templates.c -text
src/Test.c -text
src/Timpctl.c -text
......@@ -293,6 +292,7 @@ src/Vargen.c -text
src/Varrms.c -text
src/Vertint.c -text
src/Vertstat.c -text
src/Wct.c -text
src/Wind.c -text
src/Writegrid.c -text
src/Writerandom.c -text
......
2007-01-02 Ralf Quast <ralf.quast@brockmann-consult.de>
* rename tchill to wct
* split eca_strwind into eca_strwin, eca_strbre, eca_strgal and eca_hurr
2006-12-14 Uwe Schulzweida <Uwe.Schulzweida@zmaw.de>
* using CDI library version 1.0.5
......
......@@ -75,6 +75,7 @@ Operator catalog:
Seltime selyear Select years
Seltime selseas Select seasons
Seltime seldate Select dates
Seltime selsmon Select single month
Selbox sellonlatbox Select a longitude/latitude box
Selbox selindexbox Select an index box
-------------------------------------------------------------
......@@ -376,45 +377,45 @@ Operator catalog:
Vardup varmul Multiply variables
Rotuv rotuvb Backward rotation
Mastrfu mastrfu Mass stream function
Hi hi Humidity index
Tchill tchill Windchill temperature
Hi hi Humidity index (C)
Wct wct Windchill temperature (C)
-------------------------------------------------------------
ECA indices
-------------------------------------------------------------
EcaCdd eca_cdd Consecutive dry days
EcaCfd eca_cfd Consecutive frost days
EcaCsu eca_csu Consecutive summer days
EcaCwd eca_cwd Consecutive wet days
EcaCwdi eca_cwdi Cold wave duration index
EcaCwfi eca_cwfi Cold-spell days
EcaCdd eca_cdd Consecutive dry days index per time period
EcaCfd eca_cfd Consecutive frost days index per time period
EcaCsu eca_csu Consecutive summer days index per time period
EcaCwd eca_cwd Consecutive wet days index per time period
EcaCwdi eca_cwdi Cold wave duration index wrt mean of reference period
EcaCwfi eca_cwfi Cold-spell days index wrt 10th percentile of reference period
EcaEtr eca_etr Intra-period extreme temperature range
EcaFd eca_fd Frost days
EcaFdns eca_fdns Frost days where no snow
EcaGsl eca_gsl Growing season length
EcaHd eca_hd Heating degree days
EcaHwdi eca_hwdi Heat wave duration index
EcaHwfi eca_hwfi Warm-spell days
EcaId eca_id Ice days
EcaR10mm eca_r10mm Heavy precipitation days
EcaR20mm eca_r20mm Very heavy precipitation days
EcaFd eca_fd Frost days index per time period
EcaFdns eca_fdns Frost days where no snow index per time period
EcaGsl eca_gsl Growing season length index
EcaHd eca_hd Heating degree days per time period
EcaHwdi eca_hwdi Heat wave duration index wrt mean of reference period
EcaHwfi eca_hwfi Warm spell days index wrt 90th percentile of reference period
EcaId eca_id Ice days index per time period
EcaR10mm eca_r10mm Heavy precipitation days index per time period
EcaR20mm eca_r20mm Very heavy precipitation days index per time period
EcaR75p eca_r75p Moderate wet days wrt 75th percentile of reference period
EcaR75ptot eca_r75ptot Precipitation fraction due to R75p days
EcaR90p eca_r90p Very wet days wrt 90th percentile of reference period
EcaR90ptot eca_r90ptot Precipitation fraction due to R90p days
EcaR75ptot eca_r75ptot Precipitation percent due to R75p days
EcaR90p eca_r90p Wet days wrt 90th percentile of reference period
EcaR90ptot eca_r90ptot Precipitation percent due to R90p days
EcaR95p eca_r95p Very wet days wrt 95th percentile of reference period
EcaR95ptot eca_r95ptot Precipitation fraction due to R95p days
EcaR95ptot eca_r95ptot Precipitation percent due to R95p days
EcaR99p eca_r99p Extremely wet days wrt 99th percentile of reference period
EcaR99ptot eca_r99ptot Precipitation fraction due to R99p days
EcaRr1 eca_rr1 Wet days
EcaRx1day eca_rx1day Highest one-day precipitation amount
EcaRx5day eca_rx5day Highest five-day precipitation amount
EcaSdii eca_sdii Simple daily intensity index
EcaStrwind eca_strwind Strong wind days
EcaSu eca_su Summer days
EcaTg10p eca_tg10p Cold days wrt 10th percentile of reference period
EcaTg90p eca_tg90p Warm days wrt 90th percentile of reference period
EcaTn10p eca_tn10p Cold nights wrt 10th percentile of reference period
EcaTn90p eca_tn90p Warm nights wrt 90th percentile of reference period
EcaTr eca_tr Tropical nights
EcaTx10p eca_tx10p Cold days wrt 10th percentile of reference period
EcaTx90p eca_tx90p Warm days wrt 90th percentile of reference period
EcaR99ptot eca_r99ptot Precipitation percent due to R99p days
EcaRr1 eca_rr1 Wet days index per time period
EcaRx1day eca_rx1day Highest one day precipitation amount per time period
EcaRx5day eca_rx5day Highest five-day precipitation amount per time period
EcaSdii eca_sdii Simple daily intensity index per time period
EcaStrwin eca_strwin Strong wind days index per time period
EcaSu eca_su Summer days index per time period
EcaTg10p eca_tg10p Cold days percent wrt 10th percentile of reference period
EcaTg90p eca_tg90p Warm days percent wrt 90th percentile of reference period
EcaTn10p eca_tn10p Cold nights percent wrt 10th percentile of reference period
EcaTn90p eca_tn90p Warm nights percent wrt 90th percentile of reference period
EcaTr eca_tr Tropical nights index per time period
EcaTx10p eca_tx10p Very cold days percent wrt 10th percentile of reference period
EcaTx90p eca_tx90p Very warm days percent wrt 90th percentile of reference period
No preview for this file type
......@@ -84,7 +84,7 @@ Vardup Miscellaneous
Rotuv Miscellaneous
Mastrfu Miscellaneous
Hi Miscellaneous
Tchill Miscellaneous
Wct Miscellaneous
EcaCdd ECA indices
EcaCfd ECA indices
EcaCsu ECA indices
......@@ -113,7 +113,7 @@ EcaRr1 ECA indices
EcaRx1day ECA indices
EcaRx5day ECA indices
EcaSdii ECA indices
EcaStrwind ECA indices
EcaStrwin ECA indices
EcaSu ECA indices
EcaTg10p ECA indices
EcaTg90p ECA indices
......
......@@ -146,18 +146,14 @@
\end{picture}
\begin{flushright}
\large\bf{Climate Data Operators \\ Version 1.0.6 \\ December 2006}
\large\bf{Climate Data Operators \\ Version 1.0.7 \\ January 2007}
\end{flushright}
\vfill
\Large\bf{Uwe Schulzweida, Luis Kornblueh}
\Large{\bf Uwe Schulzweida, Luis Kornblueh -- \sl MPI for Meteorology}
\Large{Max-Planck-Institute for Meteorology}
\Large\bf{Ralf Quast}
\Large{Brockmann Consult}
\Large{\bf Ralf Quast -- \sl Brockmann Consult}
\begin{picture}(16,1)
\linethickness{1.0mm}
......
......@@ -10,7 +10,7 @@
\put(0,0.0){\line(1,0){3.95}}
\end{picture}
\begin{flushright}
{\small{Climate Data Operators \\ Version 1.0.6 \\ December 2006}}
{\small{Climate Data Operators \\ Version 1.0.7 \\ January 2007}}
\end{flushright}
\vspace*{0mm}
......
@BeginModule
@Name = EcaCfd
@Title = Consecutive frost days
@Title = Consecutive frost days index per time period
@Section = ECA indices
@Class = ECA index
@Arguments = ifile ofile
......@@ -17,7 +17,7 @@ the last contributing time step in @file{ifile}.
@BeginOperator_eca_cfd
@Title = Consecutive frost days
@Title = Consecutive frost days index per time period
@EndOperator
......
@BeginModule
@Name = EcaHwdi
@Title = Heat wave duration index
@Title = Heat wave duration index wrt mean of reference period
@Section = ECA indices
@Class = ECA index
@Arguments = ifile1 ifile2 ofile
......@@ -22,7 +22,7 @@ the last contributing time step in @file{ifile1}.
@BeginOperator_eca_hwdi
@Title = Heat wave duration index
@Title = Heat wave duration index wrt mean of reference period
@Parameter = [nday] [T]
@EndOperator
......
@BeginModule
@NewPage
@Name = EcaHwfi
@Title = Warm-spell days index wrt 90th percentile of reference period
@Title = Warm spell days index wrt 90th percentile of reference period
@Section = ECA indices
@Class = ECA index
@Arguments = ifile1 ifile2 ofile
......@@ -23,7 +23,7 @@ the last contributing time step in @file{ifile1}.
@BeginOperator_eca_hwfi
@Title = Warm-spell days index wrt 90th percentile of reference period
@Title = Warm spell days index wrt 90th percentile of reference period
@Parameter = [nday]
@EndOperator
......
@BeginModule
@Name = EcaR20mm
@Title = Very heavy precipitation days
@Title = Very heavy precipitation days index per time period
@Section = ECA indices
@Class = ECA index
@Arguments = ifile ofile
......@@ -8,7 +8,7 @@
@BeginDescription
Let @file{ifile} be a time series of daily precipitation amounts RR,
then counted is the number of days where RR > 20 mm.
then counted is the number of days where RR is at least 20 mm.
The date information for a time step in @file{ofile} is the date of
the last contributing time step in @file{ifile}.
@EndDescription
......@@ -16,7 +16,7 @@ the last contributing time step in @file{ifile}.
@BeginOperator_eca_r20mm
@Title = Very heavy precipitation days
@Title = Very heavy precipitation days index per time period
@EndOperator
......
@BeginModule
@Name = EcaRx1day
@Title = Highest one-day precipitation amount per time period
@Title = Highest one day precipitation amount per time period
@Section = ECA indices
@Class = ECA index
@Arguments = ifile ofile
......@@ -18,7 +18,7 @@ the last contributing time step in @file{ifile}.
@BeginOperator_eca_rx1day
@Title = Highest one-day precipitation amount per time period
@Title = Highest one day precipitation amount per time period
@Parameter = [mode]
@EndOperator
......
......@@ -19,7 +19,7 @@ the last contributing time step in @file{ifile1}.
@BeginOperator_eca_tx10p
@Title = Very old days percent wrt 10th percentile of reference period
@Title = Very cold days percent wrt 10th percentile of reference period
@EndOperator
......
@BeginModule
@Name = Tchill
@Title = Windchill temperature
@Name = Wct
@Title = Windchill temperature (°C)
@Section = Miscellaneous
@Class = Arithmetic
@Arguments = ifile1 ifile2 ofile
@Operators = tchill
@Operators = wct
@BeginDescription
Let @file{ifile1} and @file{ifile2} be time series of temperature and
wind speed records, then a corresponding time series of resulting
windchil temperatures is written to @file{ofile}.
Let @file{ifile1} and @file{ifile2} be time series of temperature and wind
speed records, then a corresponding time series of resulting windchill
temperatures is written to @file{ofile}. Note that the temperature and
wind speed records must be given in units of °C and m/s, respectively.
@EndDescription
@EndModule
@BeginOperator_tchill
@Title = Windchill temperature
@BeginOperator_wct
@Title = Windchill temperature (°C)
@EndOperator
@BeginModule
@NewPage
@Name = Yseaspctl
@Title = Multi-year seasonal percentile values
@Section = Statistical values
......
This diff is collapsed.
......@@ -18,7 +18,7 @@
/*
This module contains the following operators:
Hi hi Compute the humidity index
Hi hi Compute the humidity index
*/
......@@ -30,16 +30,23 @@
#include "cdo_int.h"
#include "pstream.h"
#define TO_DEGREE_CELSIUS(x) ((x) - 273.15)
#define TO_KELVIN(x) ((x) + 273.15)
static const char HI_NAME[] = "hum_index";
static const char HI_LONGNAME[] = "Humindex describes empirically in units of temperature how the temperature and humidity influence the wellness of a human being. HI = T + 5/9 * (A - 10) with A = e * (6.112 * 10 ** ((7.5 * T)/(237.7 + T)) * R), T = air temperature in degree Celsius, R = relative humidity, e = vapour pressure. Humindex is only defined for temperatures of at least 26 degree Celsius and relative humidity of at least 40 percent.";
static const char HI_UNITS[] = "Celsius";
static const int FIRST_VAR = 0;
static double humidityIndex(double t, double e, double r, double missval)
{
const double c = TO_DEGREE_CELSIUS(t);
const double a = 0.01 * r * e * 6.112 * pow(10.0, 7.5 * c / (237.7 + c));
static const double tmin = 26.0;
static const double hmin = 0.40;
return TO_KELVIN(c + (5.0 / 9.0) * (a - 10.0));
if ( t < tmin || r < hmin )
return missval;
return t + (5.0 / 9.0) * ((0.01 * r * e * 6.112 * pow(10.0, (7.5 * t) / (237.7 + t))) - 10.0);
}
......@@ -88,12 +95,13 @@ static void farexpr(FIELD *field1, FIELD field2, FIELD field3, double (*expressi
void *Hi(void *argument)
{
static char func[] = "Hi";
enum {FILL_NONE, FILL_REC, FILL_TS};
int streamIDm, streamIDs, streamID1, streamID2, streamID3, streamID4;
int streamID1, streamID2, streamID3, streamID4;
int gridsize;
int nrecs, nrecs2, nrecs3, recID;
int tsID;
int varID, levelID;
int gridID, zaxisID;
int varID1, varID2, varID3, varID4;
int levelID1, levelID2, levelID3;
int vlistID1, vlistID2, vlistID3, vlistID4;
int taxisID1, taxisID2, taxisID3, taxisID4;
FIELD field1, field2, field3;
......@@ -108,9 +116,6 @@ void *Hi(void *argument)
streamID3 = streamOpenRead(cdoStreamName(2));
if ( streamID3 < 0 ) cdiError(streamID3, "Open failed on %s", cdoStreamName(2));
streamIDm = streamID1;
streamIDs = streamID2;
vlistID1 = streamInqVlist(streamID1);
vlistID2 = streamInqVlist(streamID2);
vlistID3 = streamInqVlist(streamID3);
......@@ -134,11 +139,22 @@ void *Hi(void *argument)
cdoPrint("Number of timesteps: file1 %d, file2 %d, file3 %d",
vlistNtsteps(vlistID1), vlistNtsteps(vlistID2), vlistNtsteps(vlistID3));
vlistID4 = vlistDuplicate(vlistID1);
vlistID4 = vlistCreate();
gridID = vlistInqVarGrid(vlistID1, FIRST_VAR);
zaxisID = vlistInqVarZaxis(vlistID1, FIRST_VAR);
varID4 = vlistDefVar(vlistID4, gridID, zaxisID, TIME_VARIABLE);
taxisID4 = taxisDuplicate(taxisID1);
taxisID4 = taxisCreate(TAXIS_RELATIVE);
taxisDefTunit(taxisID4, TUNIT_SECOND);
taxisDefCalendar(taxisID4, CALENDAR_PROLEPTIC);
taxisDefRdate(taxisID4, 19550101);
taxisDefRtime(taxisID4, 0);
vlistDefTaxis(vlistID4, taxisID4);
vlistDefVarName(vlistID4, varID4, HI_NAME);
vlistDefVarLongname(vlistID4, varID4, HI_LONGNAME);
vlistDefVarUnits(vlistID4, varID4, HI_UNITS);
streamID4 = streamOpenWrite(cdoStreamName(3), cdoFiletype());
if ( streamID4 < 0 ) cdiError(streamID4, "Open failed on %s", cdoStreamName(3));
......@@ -157,27 +173,33 @@ void *Hi(void *argument)
for ( recID = 0; recID < nrecs; recID++ )
{
streamInqRecord(streamID1, &varID, &levelID);
streamInqRecord(streamID1, &varID1, &levelID1);
streamReadRecord(streamID1, field1.ptr, &field1.nmiss);
streamInqRecord(streamID2, &varID, &levelID);
streamInqRecord(streamID2, &varID2, &levelID2);
streamReadRecord(streamID2, field2.ptr, &field2.nmiss);
streamInqRecord(streamID3, &varID, &levelID);
streamInqRecord(streamID3, &varID3, &levelID3);
streamReadRecord(streamID3, field3.ptr, &field3.nmiss);
field1.grid = vlistInqVarGrid(vlistID1, varID);
field1.missval = vlistInqVarMissval(vlistID1, varID);
field2.grid = vlistInqVarGrid(vlistID2, varID);
field2.missval = vlistInqVarMissval(vlistID2, varID);
field3.grid = vlistInqVarGrid(vlistID3, varID);
field3.missval = vlistInqVarMissval(vlistID3, varID);
if ( varID1 != varID2 || varID1 != varID3 || levelID1 != levelID2 || levelID1 != levelID3 )
cdoAbort("Input streams have different structure!");
if ( varID1 != FIRST_VAR )
continue;
field1.grid = vlistInqVarGrid(vlistID1, varID1);
field1.missval = vlistInqVarMissval(vlistID1, varID1);
field2.grid = vlistInqVarGrid(vlistID2, varID2);
field2.missval = vlistInqVarMissval(vlistID2, varID2);
field3.grid = vlistInqVarGrid(vlistID3, varID3);
field3.missval = vlistInqVarMissval(vlistID3, varID3);
farexpr(&field1, field2, field3, humidityIndex);
streamDefRecord(streamID4, varID, levelID);
streamDefRecord(streamID4, varID4, levelID1);
streamWriteRecord(streamID4, field1.ptr, field1.nmiss);
}
......
......@@ -183,7 +183,7 @@ cdo_SOURCES = Arith.c \
ecautil.c \
EcaIndices.c \
Hi.c \
Tchill.c \
Wct.c \
cdi.h
#
cdo_LDADD = $(LDFLAGS) -lm
......
......@@ -264,7 +264,7 @@ cdo_SOURCES = Arith.c \
ecautil.c \
EcaIndices.c \
Hi.c \
Tchill.c \
Wct.c \
cdi.h
#
......@@ -335,7 +335,7 @@ am_cdo_OBJECTS = Arith.$(OBJEXT) Arithc.$(OBJEXT) Arithdays.$(OBJEXT) \
vinterp.$(OBJEXT) zaxis.$(OBJEXT) pthread_debug.$(OBJEXT) \
color.$(OBJEXT) list.$(OBJEXT) percentiles.$(OBJEXT) \
nth_element.$(OBJEXT) ecacore.$(OBJEXT) ecautil.$(OBJEXT) \
EcaIndices.$(OBJEXT) Hi.$(OBJEXT) Tchill.$(OBJEXT)
EcaIndices.$(OBJEXT) Hi.$(OBJEXT) Wct.$(OBJEXT)
cdo_OBJECTS = $(am_cdo_OBJECTS)
cdo_DEPENDENCIES =
cdo_LDFLAGS =
......@@ -390,13 +390,13 @@ am__depfiles_maybe = depfiles
@AMDEP_TRUE@ ./$(DEPDIR)/Specinfo.Po ./$(DEPDIR)/Spectral.Po \
@AMDEP_TRUE@ ./$(DEPDIR)/Split.Po ./$(DEPDIR)/Splitrec.Po \
@AMDEP_TRUE@ ./$(DEPDIR)/Splittime.Po ./$(DEPDIR)/Splityear.Po \
@AMDEP_TRUE@ ./$(DEPDIR)/Subtrend.Po ./$(DEPDIR)/Tchill.Po \
@AMDEP_TRUE@ ./$(DEPDIR)/Templates.Po ./$(DEPDIR)/Test.Po \
@AMDEP_TRUE@ ./$(DEPDIR)/Timpctl.Po ./$(DEPDIR)/Timsort.Po \
@AMDEP_TRUE@ ./$(DEPDIR)/Timstat.Po ./$(DEPDIR)/Trend.Po \
@AMDEP_TRUE@ ./$(DEPDIR)/Trms.Po ./$(DEPDIR)/Vardup.Po \
@AMDEP_TRUE@ ./$(DEPDIR)/Vargen.Po ./$(DEPDIR)/Varrms.Po \
@AMDEP_TRUE@ ./$(DEPDIR)/Vertint.Po ./$(DEPDIR)/Vertstat.Po \
@AMDEP_TRUE@ ./$(DEPDIR)/Subtrend.Po ./$(DEPDIR)/Templates.Po \
@AMDEP_TRUE@ ./$(DEPDIR)/Test.Po ./$(DEPDIR)/Timpctl.Po \
@AMDEP_TRUE@ ./$(DEPDIR)/Timsort.Po ./$(DEPDIR)/Timstat.Po \
@AMDEP_TRUE@ ./$(DEPDIR)/Trend.Po ./$(DEPDIR)/Trms.Po \
@AMDEP_TRUE@ ./$(DEPDIR)/Vardup.Po ./$(DEPDIR)/Vargen.Po \
@AMDEP_TRUE@ ./$(DEPDIR)/Varrms.Po ./$(DEPDIR)/Vertint.Po \
@AMDEP_TRUE@ ./$(DEPDIR)/Vertstat.Po ./$(DEPDIR)/Wct.Po \
@AMDEP_TRUE@ ./$(DEPDIR)/Wind.Po ./$(DEPDIR)/Writegrid.Po \
@AMDEP_TRUE@ ./$(DEPDIR)/Writerandom.Po ./$(DEPDIR)/Ydaypctl.Po \
@AMDEP_TRUE@ ./$(DEPDIR)/Ydaystat.Po ./$(DEPDIR)/Ydrunpctl.Po \
......@@ -581,7 +581,6 @@ distclean-compile:
@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/Splittime.Po@am__quote@
@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/Splityear.Po@am__quote@
@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/Subtrend.Po@am__quote@
@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/Tchill.Po@am__quote@
@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/Templates.Po@am__quote@
@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/Test.Po@am__quote@
@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/Timpctl.Po@am__quote@
......@@ -594,6 +593,7 @@ distclean-compile:
@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/Varrms.Po@am__quote@
@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/Vertint.Po@am__quote@
@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/Vertstat.Po@am__quote@
@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/Wct.Po@am__quote@
@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/Wind.Po@am__quote@
@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/Writegrid.Po@am__quote@
@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/Writerandom.Po@am__quote@
......
......@@ -2,7 +2,7 @@
This file is part of CDO. CDO is a collection of Operators to
manipulate and analyse Climate model Data.
Copyright (C) 2006 Brockmann Consult
Copyright (C) 2003-2006 Brockmann Consult
See COPYING file for copying and redistribution conditions.
This program is free software; you can redistribute it and/or modify
......@@ -18,7 +18,7 @@
/*
This module contains the following operators:
Tchill tchill Compute the windchill temperature (K)
Wct wct Compute the windchill temperature (degree C)
*/
......@@ -30,14 +30,20 @@
#include "cdo_int.h"
#include "pstream.h"
#define TO_KELVIN(x) ((x) + 273.15)
static const char WCT_NAME[] = "wind_chill_temperature";
static const char WCT_LONGNAME[] = "Windchill temperature describes the fact that low temperatures are felt to be even lower in case of wind. It is based on the rate of heat loss from exposed skin caused by wind and cold. It is calculated according to the empirical formula: 33 + (T - 33) * (0.478 + 0.237 * (SQRT(ff*3.6) - 0.0124 * ff * 3.6) with T = air temperature in degree Celsius, ff = 10m wind speed in m/s. Windchill temperature is only defined for temperatures at or below 33 degree Celsius and wind speeds above 1.39 m/s. It is mainly used for freezing temperatures.";
static const char WCT_UNITS[] = "Celsius";
static const int FIRST_VAR = 0;
static double windchillTemperature(double t, double v, double missval)
{
static const double tmax = TO_KELVIN(33.0);
static const double tmax = 33.0;
static const double vmin = 1.39; /* minimum wind speed (m/s) */
return v < 0.0 || t > tmax ? missval : tmax + (0.478 + 0.237 * sqrt(v) - 0.0124 * v) * (t - tmax);
return v < vmin || t > tmax ? missval : tmax + (0.478 + 0.237 * sqrt(v) - 0.0124 * v) * (t - tmax);
}
......@@ -79,197 +85,100 @@ static void farexpr(FIELD *field1, FIELD field2, double (*expression)(double, do
}
void *Tchill(void *argument)
void *Wct(void *argument)
{
static char func[] = "Tchill";
enum {FILL_NONE, FILL_REC, FILL_TS};
int streamIDm, streamIDs, streamID1, streamID2, streamID3;
static char func[] = "Wct";
int streamID1, streamID2, streamID3;
int gridsize;
int nrecs, nrecs2, nvars = 0, nlev, recID;
int nrecs, nrecs2, recID;
int tsID;
int varID, levelID;
int offset;
int vlistIDm, vlistIDs, vlistID1, vlistID2, vlistID3;
int taxisIDm, taxisID1, taxisID2, taxisID3;
int filltype = FILL_NONE;
FIELD *fieldm, *fields, fieldtmp, field1, field2;
int **varnmiss = NULL;
double **vardata = NULL;
int gridID, zaxisID;
int varID1, varID2, varID3;
int levelID1, levelID2;
int vlistID1, vlistID2, vlistID3;
int taxisID1, taxisID2, taxisID3;
FIELD field1, field2;
cdoInitialize(argument);
cdoOperatorAdd("tchill", 0, 0, NULL);
cdoOperatorAdd("hi", 0, 0, NULL);
streamID1 = streamOpenRead(cdoStreamName(0));
if ( streamID1 < 0 ) cdiError(streamID1, "Open failed on %s", cdoStreamName(0));
streamID2 = streamOpenRead(cdoStreamName(1));
if ( streamID2 < 0 ) cdiError(streamID2, "Open failed on %s", cdoStreamName(1));
streamIDm = streamID1;
streamIDs = streamID2;
fieldm = &field1;
fields = &field2;
vlistID1 = streamInqVlist(streamID1);
vlistID2 = streamInqVlist(streamID2);
vlistIDm = vlistID1;
vlistIDs = vlistID2;
taxisID1 = vlistInqTaxis(vlistID1);
taxisID2 = vlistInqTaxis(vlistID2);
taxisIDm = taxisID1;
if ( vlistNrecs(vlistID1) != 1 && vlistNrecs(vlistID2) == 1 )
{
filltype = FILL_TS;
cdoPrint("Filling up stream2 >%s< by copying the first record.", cdoStreamName(1));
}
else if ( vlistNrecs(vlistID1) == 1 && vlistNrecs(vlistID2) != 1 )
{
filltype = FILL_TS;
cdoPrint("Filling up stream1 >%s< by copying the first record.", cdoStreamName(0));
streamIDm = streamID2;
streamIDs = streamID1;
vlistIDm = vlistID2;
vlistIDs = vlistID1;
taxisIDm = taxisID2;
fieldm = &field2;
fields = &field1;
}
if ( filltype == FILL_NONE )
vlistCompare(vlistID1, vlistID2, func_sft);
vlistCompare(vlistID1, vlistID2, func_sft);
nospec(vlistID1);