Annotation of /trunk/gadget/suitfunc.cc

Parent Directory | Revision Log

---

Revision 1 - (view) (download)

1 :	agomez	1	#include "suitfunc.h"
2 :			#include "errorhandler.h"
3 :			#include "gadget.h"
4 :			#include "global.h"
5 :
6 :			// ********************************************************
7 :			// Functions for base suitability function
8 :			// ********************************************************
9 :			void SuitFunc::readConstants(CommentStream& infile,
10 :			const TimeClass* const TimeInfo, Keeper* const keeper) {
11 :
12 :			coeff.read(infile, TimeInfo, keeper);
13 :			coeff.Update(TimeInfo);
14 :			}
15 :
16 :			void SuitFunc::setPredLength(double length) {
17 :			handle.logMessage(LOGWARN, "Warning in suitability - trying to set predator length for", this->getName());
18 :			}
19 :
20 :			void SuitFunc::setPreyLength(double length) {
21 :			handle.logMessage(LOGWARN, "Warning in suitability - trying to set prey length for", this->getName());
22 :			}
23 :
24 :			const ModelVariableVector& SuitFunc::getConstants() const {
25 :			return coeff;
26 :			}
27 :
28 :			void SuitFunc::updateConstants(const TimeClass* const TimeInfo) {
29 :			coeff.Update(TimeInfo);
30 :			}
31 :
32 :			int SuitFunc::didChange(const TimeClass* const TimeInfo) {
33 :			return coeff.didChange(TimeInfo);
34 :			}
35 :
36 :			// ********************************************************
37 :			// Functions for ExpSuitFuncA suitability function
38 :			// ********************************************************
39 :			ExpSuitFuncA::ExpSuitFuncA() : SuitFunc("ExponentialSuitFunc") {
40 :			coeff.setsize(4);
41 :			preyLength = -1.0;
42 :			predLength = -1.0;
43 :			}
44 :
45 :			double ExpSuitFuncA::calculate() {
46 :			double check = 0.0;
47 :
48 :			if (coeff[0] < 0.0 && coeff[1] < 0.0)
49 :			check = coeff[3] / (1.0 + exp(-(coeff[0] - (coeff[1] * preyLength) + (coeff[2] * predLength))));
50 :			else if (coeff[0] > 0.0 && coeff[1] > 0.0)
51 :			check = coeff[3] / (1.0 + exp(-(-coeff[0] + (coeff[1] * preyLength) + (coeff[2] * predLength))));
52 :			else
53 :			check = coeff[3] / (1.0 + exp(-(coeff[0] + (coeff[1] * preyLength) + (coeff[2] * predLength))));
54 :
55 :			if (check != check) { //check for NaN
56 :			handle.logMessageNaN(LOGWARN, "exponential suitability function");
57 :			return 0.0;
58 :			}
59 :
60 :			if (check < 0.0) {
61 :			handle.logMessage(LOGWARN, "Warning in suitability - function outside bounds", check);
62 :			return 0.0;
63 :			} else if (check > 1.0) {
64 :			handle.logMessage(LOGWARN, "Warning in suitability - function outside bounds", check);
65 :			return 1.0;
66 :			} else
67 :			return check;
68 :			}
69 :
70 :			// ********************************************************
71 :			// Functions for ConstSuitFunc suitability function
72 :			// ********************************************************
73 :			ConstSuitFunc::ConstSuitFunc() : SuitFunc("ConstantSuitFunc") {
74 :			coeff.setsize(1);
75 :			}
76 :
77 :			double ConstSuitFunc::calculate() {
78 :			if (coeff[0] < 0.0) {
79 :			handle.logMessage(LOGWARN, "Warning in suitability - function outside bounds", coeff[0]);
80 :			return 0.0;
81 :			} else if (coeff[0] > 1.0) {
82 :			handle.logMessage(LOGWARN, "Warning in suitability - function outside bounds", coeff[0]);
83 :			return 1.0;
84 :			} else
85 :			return coeff[0];
86 :			}
87 :
88 :			// ********************************************************
89 :			// Functions for AndersenSuitFunc suitability function
90 :			// ********************************************************
91 :			AndersenSuitFunc::AndersenSuitFunc() : SuitFunc("AndersenSuitFunc") {
92 :			coeff.setsize(5);
93 :			preyLength = -1.0;
94 :			predLength = -1.0;
95 :			}
96 :
97 :			double AndersenSuitFunc::calculate() {
98 :			double l = log(predLength / preyLength);
99 :			double e, q, check;
100 :
101 :			q = 0.0;
102 :			check = 0.0;
103 :			if (l > coeff[1])
104 :			q = coeff[3];
105 :			else
106 :			q = coeff[4];
107 :
108 :			if (isZero(q))
109 :			q = 1.0;
110 :			if (q < 0.0)
111 :			q = -q;
112 :
113 :			e = (l - coeff[1]) * (l - coeff[1]);
114 :			check = coeff[0] + coeff[2] * exp(-e / q);
115 :			if (check < 0.0) {
116 :			handle.logMessage(LOGWARN, "Warning in suitability - function outside bounds", check);
117 :			return 0.0;
118 :			} else if (check > 1.0) {
119 :			handle.logMessage(LOGWARN, "Warning in suitability - function outside bounds", check);
120 :			return 1.0;
121 :			} else
122 :			return check;
123 :			}
124 :
125 :			// ********************************************************
126 :			// Functions for ExpSuitFuncL50 suitability function
127 :			// ********************************************************
128 :			ExpSuitFuncL50::ExpSuitFuncL50() : SuitFunc("ExponentialL50SuitFunc") {
129 :			coeff.setsize(2);
130 :			preyLength = -1.0;
131 :			}
132 :
133 :			double ExpSuitFuncL50::calculate() {
134 :			double check = 1.0 / (1.0 + exp(-1.0 * coeff[0] * (preyLength - coeff[1])));
135 :
136 :			if (check != check) { //check for NaN
137 :			handle.logMessageNaN(LOGWARN, "exponential l50 suitability function");
138 :			return 0.0;
139 :			}
140 :
141 :			if (check < 0.0) {
142 :			handle.logMessage(LOGWARN, "Warning in suitability - function outside bounds", check);
143 :			return 0.0;
144 :			} else if (check > 1.0) {
145 :			handle.logMessage(LOGWARN, "Warning in suitability - function outside bounds", check);
146 :			return 1.0;
147 :			} else
148 :			return check;
149 :			}
150 :
151 :			// ********************************************************
152 :			// Functions for StraightSuitFunc suitability function
153 :			// ********************************************************
154 :			StraightSuitFunc::StraightSuitFunc() : SuitFunc("StraightLineSuitFunc") {
155 :			coeff.setsize(2);
156 :			preyLength = -1.0;
157 :			}
158 :
159 :			double StraightSuitFunc::calculate() {
160 :			double check = coeff[0] * preyLength + coeff[1];
161 :			if (check < 0.0) {
162 :			handle.logMessage(LOGWARN, "Warning in suitability - function outside bounds", check);
163 :			return 0.0;
164 :			} else if (check > 1.0) {
165 :			handle.logMessage(LOGWARN, "Warning in suitability - function outside bounds", check);
166 :			return 1.0;
167 :			} else
168 :			return check;
169 :			}
170 :
171 :			// ********************************************************
172 :			// Functions for InverseSuitFunc suitability function
173 :			// ********************************************************
174 :			InverseSuitFunc::InverseSuitFunc() : SuitFunc("InverseSuitFunc") {
175 :			coeff.setsize(2);
176 :			preyLength = -1.0;
177 :			}
178 :
179 :			double InverseSuitFunc::calculate() {
180 :			double check = 1.0 / (1.0 + exp(-1.0 * coeff[0] * (preyLength - coeff[1])));
181 :
182 :			if (check != check) { //check for NaN
183 :			handle.logMessageNaN(LOGWARN, "inverse suitability function");
184 :			return 0.0;
185 :			}
186 :
187 :			// make this value 1 - check to switch the direction of the slope
188 :			check = 1.0 - check;
189 :
190 :			if (check < 0.0) {
191 :			handle.logMessage(LOGWARN, "Warning in suitability - function outside bounds", check);
192 :			return 0.0;
193 :			} else if (check > 1.0) {
194 :			handle.logMessage(LOGWARN, "Warning in suitability - function outside bounds", check);
195 :			return 1.0;
196 :			} else
197 :			return check;
198 :			}
199 :
200 :			// ********************************************************
201 :			// Functions for StraightUnboundedSuitFunc suitability function
202 :			// ********************************************************
203 :			StraightUnboundedSuitFunc::StraightUnboundedSuitFunc() : SuitFunc("StraightLineUnboundedSuitFunc") {
204 :			coeff.setsize(2);
205 :			preyLength = -1.0;
206 :			}
207 :
208 :			double StraightUnboundedSuitFunc::calculate() {
209 :			double check = coeff[0] * preyLength + coeff[1];
210 :			if (check < 0.0) {
211 :			handle.logMessage(LOGWARN, "Warning in suitability - function outside bounds", check);
212 :			return 0.0;
213 :			} else
214 :			return check;
215 :			}
216 :
217 :			// ********************************************************
218 :			// Functions for Richards suitability function
219 :			// ********************************************************
220 :			RichardsSuitFunc::RichardsSuitFunc() : SuitFunc("RichardsSuitFunc") {
221 :			coeff.setsize(5);
222 :			preyLength = -1.0;
223 :			predLength = -1.0;
224 :			}
225 :
226 :			double RichardsSuitFunc::calculate() {
227 :
228 :			if (isZero(coeff[4])) {
229 :			handle.logMessage(LOGWARN, "Warning in suitability - divide by zero error");
230 :			return 1.0;
231 :			}
232 :
233 :			double check = 0.0;
234 :			if (coeff[0] < 0.0 && coeff[1] < 0.0)
235 :			check = pow(coeff[3] / (1.0 + exp(-(coeff[0] - coeff[1] * preyLength + coeff[2] * predLength))), (1.0 / coeff[4]));
236 :			else if (coeff[0] > 0.0 && coeff[1] > 0.0)
237 :			check = pow(coeff[3] / (1.0 + exp(-(-coeff[0] + coeff[1] * preyLength + coeff[2] * predLength))), (1.0 / coeff[4]));
238 :			else
239 :			check = pow(coeff[3] / (1.0 + exp(-(coeff[0] + coeff[1] * preyLength + coeff[2] * predLength))), (1.0 / coeff[4]));
240 :
241 :			if (check != check) { //check for NaN
242 :			handle.logMessageNaN(LOGWARN, "richards suitability function");
243 :			return 0.0;
244 :			}
245 :
246 :			if (check < 0.0) {
247 :			handle.logMessage(LOGWARN, "Warning in suitability - function outside bounds", check);
248 :			return 0.0;
249 :			} else if (check > 1.0) {
250 :			handle.logMessage(LOGWARN, "Warning in suitability - function outside bounds", check);
251 :			return 1.0;
252 :			} else
253 :			return check;
254 :			}
255 :
256 :			// ********************************************************
257 :			// Functions for GammaSuitFunc suitability function
258 :			// ********************************************************
259 :			GammaSuitFunc::GammaSuitFunc() : SuitFunc("GammaSuitFunc") {
260 :			coeff.setsize(3);
261 :			preyLength = -1.0;
262 :			}
263 :
264 :			double GammaSuitFunc::calculate() {
265 :
266 :			if (isZero(coeff[1]) || (isZero(coeff[2])) || (isEqual(coeff[0], 1.0))) {
267 :			handle.logMessage(LOGWARN, "Warning in suitability - divide by zero error");
268 :			return 1.0;
269 :			}
270 :
271 :			double check = exp(coeff[0] - 1.0 - (preyLength / (coeff[1] * coeff[2])));
272 :			check *= pow(preyLength / ((coeff[0] - 1.0) * coeff[1] * coeff[2]), (coeff[0] - 1.0));
273 :
274 :			if (check != check) { //check for NaN
275 :			handle.logMessageNaN(LOGWARN, "gamma suitability function");
276 :			return 0.0;
277 :			}
278 :
279 :			if (check < 0.0) {
280 :			handle.logMessage(LOGWARN, "Warning in suitability - function outside bounds", check);
281 :			return 0.0;
282 :			} else if (check > 1.0) {
283 :			handle.logMessage(LOGWARN, "Warning in suitability - function outside bounds", check);
284 :			return 1.0;
285 :			} else
286 :			return check;
287 :			}
288 :
289 :			// ********************************************************
290 :			// Functions for AndersenFleetSuitFunc suitability function
291 :			// ********************************************************
292 :			AndersenFleetSuitFunc::AndersenFleetSuitFunc() : SuitFunc("AndersenFleetSuitFunc") {
293 :			coeff.setsize(6);
294 :			preyLength = -1.0;
295 :			}
296 :
297 :			double AndersenFleetSuitFunc::calculate() {
298 :			double l = log(coeff[5] / preyLength);
299 :			double e, q, check;
300 :
301 :			q = 0.0;
302 :			check = 0.0;
303 :			if (l > coeff[1])
304 :			q = coeff[3];
305 :			else
306 :			q = coeff[4];
307 :
308 :			if (isZero(q))
309 :			q = 1.0;
310 :			if (q < 0.0)
311 :			q = -q;
312 :
313 :			e = (l - coeff[1]) * (l - coeff[1]);
314 :			check = coeff[0] + coeff[2] * exp(-e / q);
315 :			if (check < 0.0) {
316 :			handle.logMessage(LOGWARN, "Warning in suitability - function outside bounds", check);
317 :			return 0.0;
318 :			} else if (check > 1.0) {
319 :			handle.logMessage(LOGWARN, "Warning in suitability - function outside bounds", check);
320 :			return 1.0;
321 :			} else
322 :			return check;
323 :			}

---

root@forge.cesga.es	ViewVC Help
Powered by ViewVC 1.0.0