ALMaSS Rabbit ODdox  1.00
The rabbit model description following ODdox protocol
weather.cpp
Go to the documentation of this file.
1 //
2 // weather.cpp
3 //
4 /*
5 *******************************************************************************************************
6 Copyright (c) 2011, Christopher John Topping, Aarhus University
7 All rights reserved.
8 
9 Redistribution and use in source and binary forms, with or without modification, are permitted provided
10 that the following conditions are met:
11 
12 Redistributions of source code must retain the above copyright notice, this list of conditions and the
13 following disclaimer.
14 Redistributions in binary form must reproduce the above copyright notice, this list of conditions and
15 the following disclaimer in the documentation and/or other materials provided with the distribution.
16 
17 THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY EXPRESS OR
18 IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND
19 FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS
20 BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING,
21 BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR
22 BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
23 LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
24 SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
25 ********************************************************************************************************
26 */
27 
28 #define _CRT_SECURE_NO_DEPRECATE
29 
30 #include <cstdlib>
31 #include <cstdio>
32 
33 #include "ls.h"
34 
35 extern void FloatToDouble(double&,float);
36 
37 // First year of weather data used from datafile. Ignored if 0.
38 static CfgInt l_weather_starting_year("WEATHER_STARTING_YEAR",
39  CFG_CUSTOM, 0);
40 
41 class Weather *g_weather;
42 
43 const int WindDirections[12][100] = {
44 {0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0},
45 {0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0},
46 {0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0},
47 {0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0},
48 {0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0},
49 {0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0},
50 {0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0},
51 {0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0},
52 {0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0},
53 {0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0},
54 {0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0},
55 {0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0}
56 };
57 
58 void Weather::Tick( void ) {
60  int weatherday = (g_date->Date() + m_datemodulus) % m_datemodulus;
61  double snowtempthreshold = -1.0;
62  if (m_snowdepth > 0.0 && m_temp[ weatherday ] < snowtempthreshold) {
63  m_snowdepth -= 0.1; // We decay the snow depth by 0.1 cm per day when temperatures are sub zero
64  }
65  else if (m_snowdepth > 0.0) // If temperatures are above 0.0, we decay snow 1 cm per degree C
66  {
67  m_snowdepth -= m_temp[ weatherday ];
68  if (m_snowdepth < 0.0) m_snowdepth = 0.0;
69  }
70 
71  if (m_temp[ weatherday ] < snowtempthreshold && m_rain[ weatherday ] > 0.0) {
72  m_snowdepth += m_rain[ weatherday ]; // rain is in mm and snow in cm, so we get the conversion for free.
73  }
74  if (m_temp[ weatherday ] > snowtempthreshold && m_rain[ weatherday ] > 0.0) {
75  m_snowdepth -= m_rain[ weatherday ]; // rain is in mm and snow in cm, so we get the conversion for free.
76  if (m_snowdepth < 0.0) m_snowdepth = 0.0;
77  }
78  if (m_snowdepth > 0.0) {
79  m_snowtoday = true;
80  }
81  else {
82  m_snowtoday = false;
83  }
84 
85  m_insolation = c_insolation[ g_date->DayInYear() ];
86  m_temptoday = m_temp[ weatherday ];
87  m_raintoday = m_rain[ weatherday ];
88  m_windtoday = m_wind[ weatherday ];
89 
90  if (m_raintoday > 0.001)
91  m_rainingtoday = true;
92  else
93  m_rainingtoday = false;
100  double rainscale = 0;
101  for (int i = 0; i < 5; i++) {
103  rainscale += GetRain( (weatherday + m_datemodulus) - (i + 1) ) * pow( 0.5, i );
104  }
105  double temp = GetTempPeriod( weatherday + m_datemodulus - 1, 5 ) / 5.0;
106  if (temp > 10.0) rainscale *= 0.5;
107  if (temp > 20.0) rainscale *= 0.5;
108  if (temp > 30.0) rainscale *= 0.5;
109  m_humiditytoday = rainscale;
111  double chance = g_rand_uni();
112  if (chance < 0.1) m_winddirtoday = 2; // South
113  else if (chance < 0.3) m_winddirtoday = 1; // East
114  else if (chance < 0.5) m_winddirtoday = 0; // North
115  else m_winddirtoday = 3; // West
116 }
117 
118 
119 
120 double Weather::GetMeanTemp( long a_date, unsigned int a_period )
121 {
122  double sum = 0.0;
123 
124  for ( unsigned int i=0; i<a_period; i++ )
125  sum += GetTemp( a_date - i );
126 
127  return sum/(double)a_period;
128 }
129 
130 
131 
132 bool Weather::GetSnow( long a_date )
133 {
134  if ( a_date == g_date->Date() ) {
135  return m_snowtoday;
136  }
137 
138  int weatherday = a_date % m_datemodulus;
139  bool snow = false;
140 
141  if ( m_temp[ weatherday ] < 0 &&
142  m_rain[ weatherday ] > 1 &&
143  random(100) < 50 ) {
144  snow = true;
145  }
146 
147  // if ( ((dayinyear<90) || (dayinyear>330)) && (rand()%4==1) )
148  // snow = true;
149 
150  return snow;
151 }
152 
153 
154 
155 double Weather::GetDDDegs( long a_date )
156 {
157  double temp = m_temp[ a_date%m_datemodulus ];
158  if ( temp < 0.0 ) {
159  temp = 0.0;
160  }
161 
162  return temp;
163 }
164 
165 
166 
167 Weather::Weather( const char* a_weatherdatafile )
168 {
169  FILE *ifile;
170  int Magic, FileFormat, NoDays, Day, Month, Year, FirstYear=0, LastYear=0;
171  float fTemp, fRain, fWind;
172  double Temp, Rain, Wind;
173  //char degree[20];
174  //double DDegs = 0.0;
175  ifile = fopen(a_weatherdatafile, "r" );
176  if (!ifile) {
177  g_msg->Warn(WARN_FILE, "Weather::Weather(): Unable to open file",
178  a_weatherdatafile );
179  exit(1);
180  }
181  fscanf( ifile, "%d", &Magic );
182 
183  if ( Magic != -1 ) {
184  // Assume old file format.
185  NoDays = Magic;
186  FileFormat = 1;
187  } else {
188  fscanf( ifile, "%d %d", &FileFormat, &NoDays );
189  }
190 
191  /*
192  switch( FileFormat ) {
193  case 1:
194  break;
195  default:
196  g_msg->Warn(WARN_FILE, "Weather::Weather(): "
197  "Unknown file format specified on line 1: ",
198  a_weatherdatafile);
199  g_msg->Warn(WARN_FILE, "First line must be: ",
200  "-1 <file format specifier> <days worth of weather data>" );
201  exit(1);
202  }
203  */
204 
205  m_rain.resize( NoDays );
206  m_wind.resize( NoDays );
207  m_temp.resize( NoDays );
208 
209  bool storing = false; // Whether we are storing weather data.
210  unsigned int index = 0;
211 
212  g_date->SetLastYear( 0 );
213 
214  for ( int i=0; i<NoDays; i++) {
215  // Read in the next data line.
216  fscanf( ifile, "%d %d %d %g %g %g",
217  &Year, &Month, &Day, &fTemp, &fWind, &fRain );
218  FloatToDouble(Temp, fTemp);
219  FloatToDouble(Rain, fRain);
220  FloatToDouble(Wind, fWind);
221  if ( Month == 2 && Day == 29 ) {
222  // Skip leap days.
223  continue;
224  }
225 
226  if ( Month == 1 && Day == 1 && !storing &&
227  (
228  l_weather_starting_year.value() == 0 ||
229  l_weather_starting_year.value() == Year
230  )
231  ) {
232  // Commence storing of data from Jan. 1st of the first
233  // year requested.
234  storing = true;
235  g_date->SetFirstYear( Year );
236  FirstYear = Year;
237  g_date->Reset();
238  }
239 
240  if ( storing ) {
241  //if (Temp > 0.0 )
242  // DDegs += Temp;
243  m_rain[ index ] = Rain;
244  m_wind[ index ] = Wind;
245  m_temp[ index ] = Temp;
246  index++;
247  }
248 
249  if ( Month == 12 && Day == 31 && storing ) {
250  // Found yet another full year worth of weather data.
251  g_date->SetLastYear( Year );
252  LastYear = Year;
253  //sprintf( degree, "%f", DDegs );
254  //g_msg->Warn( WARN_FILE,
255  // "Sum of day degrees:",
256  // degree );
257  //DDegs = 0.0;
258  }
259  }
260 
261  // Did we find at least one full year worth of data?
262  if ( g_date->GetLastYear() == 0 ) {
263  // Nope...
264  g_msg->Warn(WARN_FILE, "Weather::Weather(): Weather data file did",
265  "not contain at least one year of data!" );
266  exit(1);
267  }
268  fclose( ifile );
269  m_datemodulus = (LastYear - FirstYear + 1)*365;
270  m_snowdepth = 0;
271  Tick();
272 }
273 
274 Weather::~Weather( void )
275 {
276 }
277 
278 double Weather::GetRainPeriod( long a_date, unsigned int a_period )
279 {
280  double sum = 0.0;
281 
282  for ( unsigned int i=0; i<a_period; i++ )
283  sum += GetRain( a_date - i );
284  return sum;
285 }
286 
287 double Weather::GetWindPeriod( long a_date, unsigned int a_period )
288 {
289  double sum = 0.0;
290 
291  for ( unsigned int i=0; i<a_period; i++ )
292  sum += GetWind( a_date - i );
293 
294  return sum;
295 }
296 
297 double Weather::GetTempPeriod( long a_date, unsigned int a_period )
298 {
305  double sum = 0.0;
306 
307  for ( unsigned int i=0; i<a_period; i++ )
308  sum += GetTemp( a_date - i );
309 
310  return sum;
311 }
Weather::m_temp
vector< double > m_temp
Definition: weather.h:390
Calendar::Date
long Date(void)
Definition: calendar.h:57
Weather::GetWind
double GetWind(void)
Definition: weather.h:426
Weather::m_winddirtoday
int m_winddirtoday
Definition: weather.h:395
Weather::m_humiditytoday
double m_humiditytoday
Definition: weather.h:399
g_rand_uni
boost::variate_generator< base_generator_type &, boost::uniform_real<> > g_rand_uni
Weather::GetRain
double GetRain(void)
Definition: weather.h:424
CfgInt
Integer configurator entry class.
Definition: configurator.h:85
Weather::~Weather
~Weather(void)
Definition: weather.cpp:274
Calendar::SetLastYear
void SetLastYear(int a_year)
Definition: calendar.h:77
WindDirections
const int WindDirections[12][100]
Definition: weather.cpp:43
Calendar::Reset
void Reset(void)
Definition: calendar.cpp:40
Weather::m_snowtoday
bool m_snowtoday
Definition: weather.h:396
Weather::m_rain
vector< double > m_rain
Definition: weather.h:387
g_weather
class Weather * g_weather
Definition: weather.cpp:41
Weather::m_datemodulus
long m_datemodulus
Definition: weather.h:404
Weather::m_snowdepth
double m_snowdepth
The snow depth in cm.
Definition: weather.h:407
Weather::GetTemp
double GetTemp(void)
Get the temperature today.
Definition: weather.h:419
g_msg
class MapErrorMsg * g_msg
Definition: maperrormsg.cpp:38
Weather::m_windtoday
double m_windtoday
Definition: weather.h:394
Weather::GetDDDegs
double GetDDDegs(long a_date)
Definition: weather.cpp:155
Calendar::DayInYear
int DayInYear(void)
Definition: calendar.h:58
Calendar::SetFirstYear
void SetFirstYear(int a_year)
Definition: calendar.h:76
Weather::GetMeanTemp
double GetMeanTemp(long a_date, unsigned int a_period)
Definition: weather.cpp:120
Weather::m_raintoday
double m_raintoday
Definition: weather.h:393
Weather::Tick
void Tick(void)
Definition: weather.cpp:58
Weather::m_rainingtoday
bool m_rainingtoday
Definition: weather.h:397
Weather::GetWindPeriod
double GetWindPeriod(long a_date, unsigned int a_period)
Definition: weather.cpp:287
Weather::GetTempPeriod
double GetTempPeriod(long a_date, unsigned int a_period)
Definition: weather.cpp:297
g_date
class Calendar * g_date
Definition: calendar.cpp:38
c_insolation
const double c_insolation[365]
Definition: weather.h:16
l_weather_starting_year
static CfgInt l_weather_starting_year("WEATHER_STARTING_YEAR", CFG_CUSTOM, 0)
CfgInt::value
int value(void)
Definition: configurator.h:92
MapErrorMsg::Warn
void Warn(MapErrorState a_level, std::string a_msg1, std::string a_msg2)
Definition: maperrormsg.cpp:56
Weather::m_insolation
double m_insolation
Definition: weather.h:398
Weather::GetSnow
bool GetSnow(void)
Definition: weather.h:429
Weather::GetRainPeriod
double GetRainPeriod(long a_date, unsigned int a_period)
Definition: weather.cpp:278
Weather::Weather
Weather(const char *a_weatherdatafile)
Definition: weather.cpp:167
FloatToDouble
void FloatToDouble(double &, float)
Weather::m_temptoday
double m_temptoday
Definition: weather.h:392
Weather::m_wind
vector< double > m_wind
Definition: weather.h:388
Calendar::GetLastYear
int GetLastYear(void)
Definition: calendar.h:66
Generated on Mon Dec 5 2016 09:53:46 for ALMaSS Rabbit ODdox by   doxygen 1.8.11