Log In | Get Help   
Home My Page Projects Code Snippets Project Openings Mareframe
Summary Activity Forums Tracker Lists Tasks Docs Surveys News SCM Files
[mareframe] Diff of /trunk/gadget/grow.cc
[mareframe] / trunk / gadget / grow.cc Repository:
ViewVC logotype

Diff of /trunk/gadget/grow.cc

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 1, Mon Feb 10 17:09:07 2014 UTC revision 2, Wed Apr 29 12:55:30 2015 UTC
# Line 2  Line 2 
2  #include "grower.h"  #include "grower.h"
3  #include "agebandmatrix.h"  #include "agebandmatrix.h"
4  #include "gadget.h"  #include "gadget.h"
5    #include "doublematrix.h"
6    #include "matrix.h"
7    
8  /* Update the agebandmatrix to reflect the calculated growth  */  /* Update the agebandmatrix to reflect the calculated growth  */
9  /* Lgrowth contains the ratio of each length group that grows */  /* Lgrowth contains the ratio of each length group that grows */
10  /* by a certain number of length groups, and Wgrowth contains */  /* by a certain number of length groups, and Wgrowth contains */
11  /* the weight increase for each entry in Lgrowth              */  /* the weight increase for each entry in Lgrowth              */
12    
13    
14      void AgeBandMatrix::Grow(const Matrix& Lgrowth, const Matrix& Wgrowth) {
15        int i, lgrp, grow, maxlgrp;
16        double num, wt, tmp;
17    //    ofstream outfile;
18    //      outfile.open("aaa", ios::out);
19    //    _Lgrowth.Print(outfile);
20    //    int a;
21    //    cout << "1-----------------" << endl;
22    //    cin >> a;
23    
24    
25        maxlgrp = Lgrowth.Nrow();
26    
27    
28      //   double** lGrowth = ToMatrix(Lgrowth, maxlgrp, aux);
29      //   double** wGrowth = ToMatrix(Wgrowth, maxlgrp, aux);
30    
31    
32        for (i = 0; i < nrow; i++) {
33            int maxCol = v[i]->maxCol();
34            int minCol = v[i]->minCol();
35          //the part that grows to or above the highest length group
36          num = 0.0;
37          wt = 0.0;
38    
39          for (lgrp = maxCol - 1; lgrp >= maxCol - maxlgrp; lgrp--) {
40            for (grow = maxCol - lgrp - 1; grow < maxlgrp; grow++) {
41              tmp = Lgrowth[grow][lgrp] * (*v[i])[lgrp].N;
42              num += tmp;
43              wt += tmp * (Wgrowth[grow][lgrp] + (*v[i])[lgrp].W);
44            }
45          }
46    //      cout << "wt1:" << wt << endl;
47    
48          lgrp = maxCol - 1;
49          if (isZero(num) || (wt < verysmall)) {
50            (*v[i])[lgrp].setToZero();
51          } else {
52            (*v[i])[lgrp].W = wt / num;
53            (*v[i])[lgrp].N = num;
54          }
55    
56    
57          //the central diagonal part of the length division
58          for (lgrp = maxCol - 2; lgrp >= minCol + maxlgrp - 1; lgrp--) {
59            num = 0.0;
60            wt = 0.0;
61            for (grow = 0; grow < maxlgrp; grow++) {
62              tmp = Lgrowth[grow][lgrp - grow] * (*v[i])[lgrp - grow].N;
63              num += tmp;
64              wt += tmp * (Wgrowth[grow][lgrp - grow] + (*v[i])[lgrp - grow].W);
65            }
66    //        cout << "wt2:" << wt << endl;
67            if (isZero(num) || (wt < verysmall)) {
68              (*v[i])[lgrp].setToZero();
69            } else {
70              (*v[i])[lgrp].W = wt / num;
71              (*v[i])[lgrp].N = num;
72            }
73          }
74    
75    
76          //the lowest part of the length division
77          for (lgrp = minCol + maxlgrp - 2; lgrp >= minCol; lgrp--) {
78            num = 0.0;
79            wt = 0.0;
80            for (grow = 0; grow <= lgrp - minCol; grow++) {
81              tmp = Lgrowth[grow][lgrp - grow] * (*v[i])[lgrp - grow].N;
82              num += tmp;
83              wt += tmp * (Wgrowth[grow][lgrp - grow] + (*v[i])[lgrp - grow].W);
84            }
85    
86            if (isZero(num) || (wt < verysmall)) {
87              (*v[i])[lgrp].setToZero();
88            } else {
89              (*v[i])[lgrp].W = wt / num;
90              (*v[i])[lgrp].N = num;
91            }
92          }
93    //      cout << "wt3:" << wt << endl;
94        }
95      }
96    
97    
98    
99    
100    
101    
102    
103    
104    
105      /********************************************************************************/
106    
107  /* JMB changed to deal with very small weights a bit better   */  /* JMB changed to deal with very small weights a bit better   */
108  void AgeBandMatrix::Grow(const DoubleMatrix& Lgrowth, const DoubleMatrix& Wgrowth) {  void AgeBandMatrix::Grow(const DoubleMatrix& Lgrowth, const DoubleMatrix& Wgrowth) {
109    int i, lgrp, grow, maxlgrp;    int i, lgrp, grow, maxlgrp;
110    double num, wt, tmp;    double num, wt, tmp;
111    
112    maxlgrp = Lgrowth.Nrow();    maxlgrp = Lgrowth.Nrow();
113    
114    //   double** lGrowth = ToMatrix(Lgrowth, maxlgrp, aux);
115    //   double** wGrowth = ToMatrix(Wgrowth, maxlgrp, aux);
116    
117    
118    for (i = 0; i < nrow; i++) {    for (i = 0; i < nrow; i++) {
119            int maxCol = v[i]->maxCol();
120            int minCol = v[i]->minCol();
121      //the part that grows to or above the highest length group      //the part that grows to or above the highest length group
122      num = 0.0;      num = 0.0;
123      wt = 0.0;      wt = 0.0;
124      for (lgrp = v[i]->maxCol() - 1; lgrp >= v[i]->maxCol() - maxlgrp; lgrp--) {  
125        for (grow = v[i]->maxCol() - lgrp - 1; grow < maxlgrp; grow++) {      for (lgrp = maxCol - 1; lgrp >= maxCol - maxlgrp; lgrp--) {
126          for (grow = maxCol - lgrp - 1; grow < maxlgrp; grow++) {
127          tmp = Lgrowth[grow][lgrp] * (*v[i])[lgrp].N;          tmp = Lgrowth[grow][lgrp] * (*v[i])[lgrp].N;
128          num += tmp;          num += tmp;
129          wt += tmp * (Wgrowth[grow][lgrp] + (*v[i])[lgrp].W);          wt += tmp * (Wgrowth[grow][lgrp] + (*v[i])[lgrp].W);
130        }        }
131      }      }
132    
133      lgrp = v[i]->maxCol() - 1;      lgrp = maxCol - 1;
134      if (isZero(num) || (wt < verysmall)) {      if (isZero(num) || (wt < verysmall)) {
135        (*v[i])[lgrp].setToZero();        (*v[i])[lgrp].setToZero();
136      } else {      } else {
# Line 34  Line 138 
138        (*v[i])[lgrp].N = num;        (*v[i])[lgrp].N = num;
139      }      }
140    
141    
142      //the central diagonal part of the length division      //the central diagonal part of the length division
143      for (lgrp = v[i]->maxCol() - 2; lgrp >= v[i]->minCol() + maxlgrp - 1; lgrp--) {      for (lgrp = maxCol - 2; lgrp >= minCol + maxlgrp - 1; lgrp--) {
144        num = 0.0;        num = 0.0;
145        wt = 0.0;        wt = 0.0;
146        for (grow = 0; grow < maxlgrp; grow++) {        for (grow = 0; grow < maxlgrp; grow++) {
# Line 53  Line 158 
158      }      }
159    
160      //the lowest part of the length division      //the lowest part of the length division
161      for (lgrp = v[i]->minCol() + maxlgrp - 2; lgrp >= v[i]->minCol(); lgrp--) {      for (lgrp = minCol + maxlgrp - 2; lgrp >= minCol; lgrp--) {
162        num = 0.0;        num = 0.0;
163        wt = 0.0;        wt = 0.0;
164        for (grow = 0; grow <= lgrp - v[i]->minCol(); grow++) {        for (grow = 0; grow <= lgrp - minCol; grow++) {
165          tmp = Lgrowth[grow][lgrp - grow] * (*v[i])[lgrp - grow].N;          tmp = Lgrowth[grow][lgrp - grow] * (*v[i])[lgrp - grow].N;
166          num += tmp;          num += tmp;
167          wt += tmp * (Wgrowth[grow][lgrp - grow] + (*v[i])[lgrp - grow].W);          wt += tmp * (Wgrowth[grow][lgrp - grow] + (*v[i])[lgrp - grow].W);
# Line 78  Line 183 
183    int i, lgrp, grow, maxlgrp, age;    int i, lgrp, grow, maxlgrp, age;
184    double num, wt, matnum, tmp, ratio;    double num, wt, matnum, tmp, ratio;
185    
186    
187    maxlgrp = Lgrowth.Nrow();    maxlgrp = Lgrowth.Nrow();
188    for (i = 0; i < nrow; i++) {    for (i = 0; i < nrow; i++) {
189              int maxCol = v[i]->maxCol();
190                int minCol = v[i]->minCol();
191      age = i + minage;      age = i + minage;
192      num = 0.0;      num = 0.0;
193      wt = 0.0;      wt = 0.0;
194      matnum = 0.0;      matnum = 0.0;
195      for (lgrp = v[i]->maxCol() - 1; lgrp >= v[i]->maxCol() - maxlgrp; lgrp--) {      for (lgrp = maxCol - 1; lgrp >= maxCol - maxlgrp; lgrp--) {
196        for (grow = v[i]->maxCol() - lgrp - 1; grow < maxlgrp; grow++) {        for (grow = maxCol - lgrp - 1; grow < maxlgrp; grow++) {
197          ratio = Mat->calcMaturation(age, lgrp, grow, (*v[i])[lgrp].W);          ratio = Mat->calcMaturation(age, lgrp, grow, (*v[i])[lgrp].W);
198          tmp = Lgrowth[grow][lgrp] * (*v[i])[lgrp].N;          tmp = Lgrowth[grow][lgrp] * (*v[i])[lgrp].N;
199          matnum += (tmp * ratio);          matnum += (tmp * ratio);
# Line 94  Line 202 
202        }        }
203      }      }
204    
205      lgrp = v[i]->maxCol() - 1;      lgrp = maxCol - 1;
206      if (isZero(num) || (wt < verysmall)) {      if (isZero(num) || (wt < verysmall)) {
207        //no fish grow to this length cell        //no fish grow to this length cell
208        (*v[i])[lgrp].setToZero();        (*v[i])[lgrp].setToZero();
# Line 114  Line 222 
222        Mat->storeMatureStock(area, age, lgrp, matnum, wt / num);        Mat->storeMatureStock(area, age, lgrp, matnum, wt / num);
223      }      }
224    
225      for (lgrp = v[i]->maxCol() - 2; lgrp >= v[i]->minCol() + maxlgrp - 1; lgrp--) {  
226        for (lgrp = maxCol - 2; lgrp >= minCol + maxlgrp - 1; lgrp--) {
227        num = 0.0;        num = 0.0;
228        wt = 0.0;        wt = 0.0;
229        matnum = 0.0;        matnum = 0.0;
# Line 146  Line 255 
255        }        }
256      }      }
257    
258      for (lgrp = v[i]->minCol() + maxlgrp - 2; lgrp >= v[i]->minCol(); lgrp--) {      for (lgrp = minCol + maxlgrp - 2; lgrp >= minCol; lgrp--) {
259        num = 0.0;        num = 0.0;
260        wt = 0.0;        wt = 0.0;
261        matnum = 0.0;        matnum = 0.0;
262        for (grow = 0; grow <= lgrp - v[i]->minCol(); grow++) {        for (grow = 0; grow <= lgrp - minCol; grow++) {
263          ratio = Mat->calcMaturation(age, lgrp, grow, (*v[i])[lgrp - grow].W);          ratio = Mat->calcMaturation(age, lgrp, grow, (*v[i])[lgrp - grow].W);
264          tmp = Lgrowth[grow][lgrp - grow] * (*v[i])[lgrp - grow].N;          tmp = Lgrowth[grow][lgrp - grow] * (*v[i])[lgrp - grow].N;
265          matnum += (tmp * ratio);          matnum += (tmp * ratio);
# Line 186  Line 295 
295    double num;    double num;
296    
297    maxlgrp = Lgrowth.Nrow();    maxlgrp = Lgrowth.Nrow();
298    
299    for (i = 0; i < nrow; i++) {    for (i = 0; i < nrow; i++) {
300              int maxCol = v[i]->maxCol();
301                int minCol = v[i]->minCol();
302      num = 0.0;      num = 0.0;
303      for (lgrp = v[i]->maxCol() - 1; lgrp >= v[i]->maxCol() - maxlgrp; lgrp--)      for (lgrp = maxCol - 1; lgrp >= maxCol - maxlgrp; lgrp--)
304        for (grow = v[i]->maxCol() - lgrp - 1; grow < maxlgrp; grow++)        for (grow = maxCol - lgrp - 1; grow < maxlgrp; grow++)
305          num += (Lgrowth[grow][lgrp] * (*v[i])[lgrp].N);          num += (Lgrowth[grow][lgrp] * (*v[i])[lgrp].N);
306    
307      lgrp = v[i]->maxCol() - 1;      lgrp = maxCol - 1;
308      if (isZero(num)) {      if (isZero(num)) {
309        (*v[i])[lgrp].setToZero();        (*v[i])[lgrp].setToZero();
310      } else {      } else {
# Line 200  Line 312 
312        (*v[i])[lgrp].W = Weight[lgrp];        (*v[i])[lgrp].W = Weight[lgrp];
313      }      }
314    
315      for (lgrp = v[i]->maxCol() - 2; lgrp >= v[i]->minCol() + maxlgrp - 1; lgrp--) {      for (lgrp = maxCol - 2; lgrp >= minCol + maxlgrp - 1; lgrp--) {
316        num = 0.0;        num = 0.0;
317        for (grow = 0; grow < maxlgrp; grow++)        for (grow = 0; grow < maxlgrp; grow++)
318          num += (Lgrowth[grow][lgrp - grow] * (*v[i])[lgrp - grow].N);          num += (Lgrowth[grow][lgrp - grow] * (*v[i])[lgrp - grow].N);
# Line 213  Line 325 
325        }        }
326      }      }
327    
328      for (lgrp = v[i]->minCol() + maxlgrp - 2; lgrp >= v[i]->minCol(); lgrp--) {      for (lgrp = minCol + maxlgrp - 2; lgrp >= minCol; lgrp--) {
329        num = 0.0;        num = 0.0;
330        for (grow = 0; grow <= lgrp - v[i]->minCol(); grow++)        for (grow = 0; grow <= lgrp - minCol; grow++)
331          num += (Lgrowth[grow][lgrp - grow] * (*v[i])[lgrp - grow].N);          num += (Lgrowth[grow][lgrp - grow] * (*v[i])[lgrp - grow].N);
332    
333        if (isZero(num)) {        if (isZero(num)) {
# Line 236  Line 348 
348    double num, matnum, tmp, ratio;    double num, matnum, tmp, ratio;
349    
350    maxlgrp = Lgrowth.Nrow();    maxlgrp = Lgrowth.Nrow();
351    
352    for (i = 0; i < nrow; i++) {    for (i = 0; i < nrow; i++) {
353              int maxCol = v[i]->maxCol();
354                int minCol = v[i]->minCol();
355      age = i + minage;      age = i + minage;
356      num = 0.0;      num = 0.0;
357      matnum = 0.0;      matnum = 0.0;
358      for (lgrp = v[i]->maxCol() - 1; lgrp >= v[i]->maxCol() - maxlgrp; lgrp--) {      for (lgrp = maxCol - 1; lgrp >= maxCol - maxlgrp; lgrp--) {
359        for (grow = v[i]->maxCol() - lgrp - 1; grow < maxlgrp; grow++) {        for (grow = maxCol - lgrp - 1; grow < maxlgrp; grow++) {
360          ratio = Mat->calcMaturation(age, lgrp, grow, (*v[i])[lgrp].W);          ratio = Mat->calcMaturation(age, lgrp, grow, (*v[i])[lgrp].W);
361          tmp = Lgrowth[grow][lgrp] * (*v[i])[lgrp].N;          tmp = Lgrowth[grow][lgrp] * (*v[i])[lgrp].N;
362          matnum += (tmp * ratio);          matnum += (tmp * ratio);
# Line 249  Line 364 
364        }        }
365      }      }
366    
367      lgrp = v[i]->maxCol() - 1;      lgrp = maxCol - 1;
368      if (isZero(num)) {      if (isZero(num)) {
369        //no fish grow to this length cell        //no fish grow to this length cell
370        (*v[i])[lgrp].setToZero();        (*v[i])[lgrp].setToZero();
# Line 269  Line 384 
384        Mat->storeMatureStock(area, age, lgrp, matnum, Weight[lgrp]);        Mat->storeMatureStock(area, age, lgrp, matnum, Weight[lgrp]);
385      }      }
386    
387      for (lgrp = v[i]->maxCol() - 2; lgrp >= v[i]->minCol() + maxlgrp - 1; lgrp--) {      for (lgrp = maxCol - 2; lgrp >= minCol + maxlgrp - 1; lgrp--) {
388        num = 0.0;        num = 0.0;
389        matnum = 0.0;        matnum = 0.0;
390        for (grow = 0; grow < maxlgrp; grow++) {        for (grow = 0; grow < maxlgrp; grow++) {
# Line 299  Line 414 
414        }        }
415      }      }
416    
417      for (lgrp = v[i]->minCol() + maxlgrp - 2; lgrp >= v[i]->minCol(); lgrp--) {      for (lgrp = minCol + maxlgrp - 2; lgrp >= minCol; lgrp--) {
418        num = 0.0;        num = 0.0;
419        matnum = 0.0;        matnum = 0.0;
420        for (grow = 0; grow <= lgrp - v[i]->minCol(); grow++) {        for (grow = 0; grow <= lgrp - minCol; grow++) {
421          ratio = Mat->calcMaturation(age, lgrp, grow, (*v[i])[lgrp - grow].W);          ratio = Mat->calcMaturation(age, lgrp, grow, (*v[i])[lgrp - grow].W);
422          tmp = Lgrowth[grow][lgrp - grow] * (*v[i])[lgrp - grow].N;          tmp = Lgrowth[grow][lgrp - grow] * (*v[i])[lgrp - grow].N;
423          matnum += (tmp * ratio);          matnum += (tmp * ratio);
Legend:
Removed from v.1  
changed lines
  Added in v.2
root@forge.cesga.es
 ViewVC Help
Powered by ViewVC 1.0.0  

Powered By FusionForge