#include "maturity.h"
#include "grower.h"
#include "agebandmatrix.h"
#include "gadget.h"
#include "doublematrix.h"
#include "matrix.h"
/* Update the agebandmatrix to reflect the calculated growth */
/* Lgrowth contains the ratio of each length group that grows */
/* by a certain number of length groups, and Wgrowth contains */
/* the weight increase for each entry in Lgrowth */
void AgeBandMatrix::Grow(const Matrix& Lgrowth, const Matrix& Wgrowth) {
int i, lgrp, grow, maxlgrp;
double num, wt, tmp;
// ofstream outfile;
// outfile.open("aaa", ios::out);
// _Lgrowth.Print(outfile);
// int a;
// cout << "1-----------------" << endl;
// cin >> a;
maxlgrp = Lgrowth.Nrow();
// double** lGrowth = ToMatrix(Lgrowth, maxlgrp, aux);
// double** wGrowth = ToMatrix(Wgrowth, maxlgrp, aux);
for (i = 0; i < nrow; i++) {
int maxCol = v[i]->maxCol();
int minCol = v[i]->minCol();
//the part that grows to or above the highest length group
num = 0.0;
wt = 0.0;
for (lgrp = maxCol - 1; lgrp >= maxCol - maxlgrp; lgrp--) {
for (grow = maxCol - lgrp - 1; grow < maxlgrp; grow++) {
tmp = Lgrowth[grow][lgrp] * (*v[i])[lgrp].N;
num += tmp;
wt += tmp * (Wgrowth[grow][lgrp] + (*v[i])[lgrp].W);
}
}
// cout << "wt1:" << wt << endl;
lgrp = maxCol - 1;
if (isZero(num) || (wt < verysmall)) {
(*v[i])[lgrp].setToZero();
} else {
(*v[i])[lgrp].W = wt / num;
(*v[i])[lgrp].N = num;
}
//the central diagonal part of the length division
for (lgrp = maxCol - 2; lgrp >= minCol + maxlgrp - 1; lgrp--) {
num = 0.0;
wt = 0.0;
for (grow = 0; grow < maxlgrp; grow++) {
tmp = Lgrowth[grow][lgrp - grow] * (*v[i])[lgrp - grow].N;
num += tmp;
wt += tmp * (Wgrowth[grow][lgrp - grow] + (*v[i])[lgrp - grow].W);
}
// cout << "wt2:" << wt << endl;
if (isZero(num) || (wt < verysmall)) {
(*v[i])[lgrp].setToZero();
} else {
(*v[i])[lgrp].W = wt / num;
(*v[i])[lgrp].N = num;
}
}
//the lowest part of the length division
for (lgrp = minCol + maxlgrp - 2; lgrp >= minCol; lgrp--) {
num = 0.0;
wt = 0.0;
for (grow = 0; grow <= lgrp - minCol; grow++) {
tmp = Lgrowth[grow][lgrp - grow] * (*v[i])[lgrp - grow].N;
num += tmp;
wt += tmp * (Wgrowth[grow][lgrp - grow] + (*v[i])[lgrp - grow].W);
}
if (isZero(num) || (wt < verysmall)) {
(*v[i])[lgrp].setToZero();
} else {
(*v[i])[lgrp].W = wt / num;
(*v[i])[lgrp].N = num;
}
}
// cout << "wt3:" << wt << endl;
}
}
/********************************************************************************/
/* JMB changed to deal with very small weights a bit better */
void AgeBandMatrix::Grow(const DoubleMatrix& Lgrowth, const DoubleMatrix& Wgrowth) {
int i, lgrp, grow, maxlgrp;
double num, wt, tmp;
maxlgrp = Lgrowth.Nrow();
// double** lGrowth = ToMatrix(Lgrowth, maxlgrp, aux);
// double** wGrowth = ToMatrix(Wgrowth, maxlgrp, aux);
for (i = 0; i < nrow; i++) {
int maxCol = v[i]->maxCol();
int minCol = v[i]->minCol();
//the part that grows to or above the highest length group
num = 0.0;
wt = 0.0;
for (lgrp = maxCol - 1; lgrp >= maxCol - maxlgrp; lgrp--) {
for (grow = maxCol - lgrp - 1; grow < maxlgrp; grow++) {
tmp = Lgrowth[grow][lgrp] * (*v[i])[lgrp].N;
num += tmp;
wt += tmp * (Wgrowth[grow][lgrp] + (*v[i])[lgrp].W);
}
}
lgrp = maxCol - 1;
if (isZero(num) || (wt < verysmall)) {
(*v[i])[lgrp].setToZero();
} else {
(*v[i])[lgrp].W = wt / num;
(*v[i])[lgrp].N = num;
}
//the central diagonal part of the length division
for (lgrp = maxCol - 2; lgrp >= minCol + maxlgrp - 1; lgrp--) {
num = 0.0;
wt = 0.0;
for (grow = 0; grow < maxlgrp; grow++) {
tmp = Lgrowth[grow][lgrp - grow] * (*v[i])[lgrp - grow].N;
num += tmp;
wt += tmp * (Wgrowth[grow][lgrp - grow] + (*v[i])[lgrp - grow].W);
}
if (isZero(num) || (wt < verysmall)) {
(*v[i])[lgrp].setToZero();
} else {
(*v[i])[lgrp].W = wt / num;
(*v[i])[lgrp].N = num;
}
}
//the lowest part of the length division
for (lgrp = minCol + maxlgrp - 2; lgrp >= minCol; lgrp--) {
num = 0.0;
wt = 0.0;
for (grow = 0; grow <= lgrp - minCol; grow++) {
tmp = Lgrowth[grow][lgrp - grow] * (*v[i])[lgrp - grow].N;
num += tmp;
wt += tmp * (Wgrowth[grow][lgrp - grow] + (*v[i])[lgrp - grow].W);
}
if (isZero(num) || (wt < verysmall)) {
(*v[i])[lgrp].setToZero();
} else {
(*v[i])[lgrp].W = wt / num;
(*v[i])[lgrp].N = num;
}
}
}
}
//Same program with certain num of fish made mature.
void AgeBandMatrix::Grow(const DoubleMatrix& Lgrowth, const DoubleMatrix& Wgrowth, Maturity* const Mat, int area) {
int i, lgrp, grow, maxlgrp, age;
double num, wt, matnum, tmp, ratio;
maxlgrp = Lgrowth.Nrow();
for (i = 0; i < nrow; i++) {
int maxCol = v[i]->maxCol();
int minCol = v[i]->minCol();
age = i + minage;
num = 0.0;
wt = 0.0;
matnum = 0.0;
for (lgrp = maxCol - 1; lgrp >= maxCol - maxlgrp; lgrp--) {
for (grow = maxCol - lgrp - 1; grow < maxlgrp; grow++) {
ratio = Mat->calcMaturation(age, lgrp, grow, (*v[i])[lgrp].W);
tmp = Lgrowth[grow][lgrp] * (*v[i])[lgrp].N;
matnum += (tmp * ratio);
num += tmp;
wt += tmp * (Wgrowth[grow][lgrp] + (*v[i])[lgrp].W);
}
}
lgrp = maxCol - 1;
if (isZero(num) || (wt < verysmall)) {
//no fish grow to this length cell
(*v[i])[lgrp].setToZero();
Mat->storeMatureStock(area, age, lgrp, 0.0, 0.0);
} else if (isZero(matnum)) {
//none of the fish that grow to this length cell mature
(*v[i])[lgrp].W = wt / num;
(*v[i])[lgrp].N = num;
Mat->storeMatureStock(area, age, lgrp, 0.0, 0.0);
} else if (isEqual(num, matnum) || (matnum > num)) {
//all the fish that grow to this length cell mature
(*v[i])[lgrp].setToZero();
Mat->storeMatureStock(area, age, lgrp, num, wt / num);
} else {
(*v[i])[lgrp].W = wt / num;
(*v[i])[lgrp].N = num - matnum;
Mat->storeMatureStock(area, age, lgrp, matnum, wt / num);
}
for (lgrp = maxCol - 2; lgrp >= minCol + maxlgrp - 1; lgrp--) {
num = 0.0;
wt = 0.0;
matnum = 0.0;
for (grow = 0; grow < maxlgrp; grow++) {
ratio = Mat->calcMaturation(age, lgrp, grow, (*v[i])[lgrp - grow].W);
tmp = Lgrowth[grow][lgrp - grow] * (*v[i])[lgrp - grow].N;
matnum += (tmp * ratio);
num += tmp;
wt += tmp * (Wgrowth[grow][lgrp - grow] + (*v[i])[lgrp - grow].W);
}
if (isZero(num) || (wt < verysmall)) {
//no fish grow to this length cell
(*v[i])[lgrp].setToZero();
Mat->storeMatureStock(area, age, lgrp, 0.0, 0.0);
} else if (isZero(matnum)) {
//none of the fish that grow to this length cell mature
(*v[i])[lgrp].W = wt / num;
(*v[i])[lgrp].N = num;
Mat->storeMatureStock(area, age, lgrp, 0.0, 0.0);
} else if (isEqual(num, matnum) || (matnum > num)) {
//all the fish that grow to this length cell mature
(*v[i])[lgrp].setToZero();
Mat->storeMatureStock(area, age, lgrp, num, wt / num);
} else {
(*v[i])[lgrp].W = wt / num;
(*v[i])[lgrp].N = num - matnum;
Mat->storeMatureStock(area, age, lgrp, matnum, wt / num);
}
}
for (lgrp = minCol + maxlgrp - 2; lgrp >= minCol; lgrp--) {
num = 0.0;
wt = 0.0;
matnum = 0.0;
for (grow = 0; grow <= lgrp - minCol; grow++) {
ratio = Mat->calcMaturation(age, lgrp, grow, (*v[i])[lgrp - grow].W);
tmp = Lgrowth[grow][lgrp - grow] * (*v[i])[lgrp - grow].N;
matnum += (tmp * ratio);
num += tmp;
wt += tmp * (Wgrowth[grow][lgrp - grow] + (*v[i])[lgrp - grow].W);
}
if (isZero(num) || (wt < verysmall)) {
//no fish grow to this length cell
(*v[i])[lgrp].setToZero();
Mat->storeMatureStock(area, age, lgrp, 0.0, 0.0);
} else if (isZero(matnum)) {
//none of the fish that grow to this length cell mature
(*v[i])[lgrp].W = wt / num;
(*v[i])[lgrp].N = num;
Mat->storeMatureStock(area, age, lgrp, 0.0, 0.0);
} else if (isEqual(num, matnum) || (matnum > num)) {
//all the fish that grow to this length cell mature
(*v[i])[lgrp].setToZero();
Mat->storeMatureStock(area, age, lgrp, num, wt / num);
} else {
(*v[i])[lgrp].W = wt / num;
(*v[i])[lgrp].N = num - matnum;
Mat->storeMatureStock(area, age, lgrp, matnum, wt / num);
}
}
}
}
//fleksibest formulation - weight read in from file (should be positive)
void AgeBandMatrix::Grow(const DoubleMatrix& Lgrowth, const DoubleVector& Weight) {
int i, lgrp, grow, maxlgrp;
double num;
maxlgrp = Lgrowth.Nrow();
for (i = 0; i < nrow; i++) {
int maxCol = v[i]->maxCol();
int minCol = v[i]->minCol();
num = 0.0;
for (lgrp = maxCol - 1; lgrp >= maxCol - maxlgrp; lgrp--)
for (grow = maxCol - lgrp - 1; grow < maxlgrp; grow++)
num += (Lgrowth[grow][lgrp] * (*v[i])[lgrp].N);
lgrp = maxCol - 1;
if (isZero(num)) {
(*v[i])[lgrp].setToZero();
} else {
(*v[i])[lgrp].N = num;
(*v[i])[lgrp].W = Weight[lgrp];
}
for (lgrp = maxCol - 2; lgrp >= minCol + maxlgrp - 1; lgrp--) {
num = 0.0;
for (grow = 0; grow < maxlgrp; grow++)
num += (Lgrowth[grow][lgrp - grow] * (*v[i])[lgrp - grow].N);
if (isZero(num)) {
(*v[i])[lgrp].setToZero();
} else {
(*v[i])[lgrp].N = num;
(*v[i])[lgrp].W = Weight[lgrp];
}
}
for (lgrp = minCol + maxlgrp - 2; lgrp >= minCol; lgrp--) {
num = 0.0;
for (grow = 0; grow <= lgrp - minCol; grow++)
num += (Lgrowth[grow][lgrp - grow] * (*v[i])[lgrp - grow].N);
if (isZero(num)) {
(*v[i])[lgrp].setToZero();
} else {
(*v[i])[lgrp].N = num;
(*v[i])[lgrp].W = Weight[lgrp];
}
}
}
}
//fleksibest formulation - weight read in from file (should be positive)
//Same program with certain num of fish made mature.
void AgeBandMatrix::Grow(const DoubleMatrix& Lgrowth, const DoubleVector& Weight, Maturity* const Mat, int area) {
int i, lgrp, grow, maxlgrp, age;
double num, matnum, tmp, ratio;
maxlgrp = Lgrowth.Nrow();
for (i = 0; i < nrow; i++) {
int maxCol = v[i]->maxCol();
int minCol = v[i]->minCol();
age = i + minage;
num = 0.0;
matnum = 0.0;
for (lgrp = maxCol - 1; lgrp >= maxCol - maxlgrp; lgrp--) {
for (grow = maxCol - lgrp - 1; grow < maxlgrp; grow++) {
ratio = Mat->calcMaturation(age, lgrp, grow, (*v[i])[lgrp].W);
tmp = Lgrowth[grow][lgrp] * (*v[i])[lgrp].N;
matnum += (tmp * ratio);
num += tmp;
}
}
lgrp = maxCol - 1;
if (isZero(num)) {
//no fish grow to this length cell
(*v[i])[lgrp].setToZero();
Mat->storeMatureStock(area, age, lgrp, 0.0, 0.0);
} else if (isZero(matnum)) {
//none of the fish that grow to this length cell mature
(*v[i])[lgrp].W = Weight[lgrp];
(*v[i])[lgrp].N = num;
Mat->storeMatureStock(area, age, lgrp, 0.0, 0.0);
} else if (isEqual(num, matnum) || (matnum > num)) {
//all the fish that grow to this length cell mature
(*v[i])[lgrp].setToZero();
Mat->storeMatureStock(area, age, lgrp, num, Weight[lgrp]);
} else {
(*v[i])[lgrp].W = Weight[lgrp];
(*v[i])[lgrp].N = num - matnum;
Mat->storeMatureStock(area, age, lgrp, matnum, Weight[lgrp]);
}
for (lgrp = maxCol - 2; lgrp >= minCol + maxlgrp - 1; lgrp--) {
num = 0.0;
matnum = 0.0;
for (grow = 0; grow < maxlgrp; grow++) {
ratio = Mat->calcMaturation(age, lgrp, grow, (*v[i])[lgrp - grow].W);
tmp = Lgrowth[grow][lgrp - grow] * (*v[i])[lgrp - grow].N;
matnum += (tmp * ratio);
num += tmp;
}
if (isZero(num)) {
//no fish grow to this length cell
(*v[i])[lgrp].setToZero();
Mat->storeMatureStock(area, age, lgrp, 0.0, 0.0);
} else if (isZero(matnum)) {
//none of the fish that grow to this length cell mature
(*v[i])[lgrp].W = Weight[lgrp];
(*v[i])[lgrp].N = num;
Mat->storeMatureStock(area, age, lgrp, 0.0, 0.0);
} else if (isEqual(num, matnum) || (matnum > num)) {
//all the fish that grow to this length cell mature
(*v[i])[lgrp].setToZero();
Mat->storeMatureStock(area, age, lgrp, num, Weight[lgrp]);
} else {
(*v[i])[lgrp].W = Weight[lgrp];
(*v[i])[lgrp].N = num - matnum;
Mat->storeMatureStock(area, age, lgrp, matnum, Weight[lgrp]);
}
}
for (lgrp = minCol + maxlgrp - 2; lgrp >= minCol; lgrp--) {
num = 0.0;
matnum = 0.0;
for (grow = 0; grow <= lgrp - minCol; grow++) {
ratio = Mat->calcMaturation(age, lgrp, grow, (*v[i])[lgrp - grow].W);
tmp = Lgrowth[grow][lgrp - grow] * (*v[i])[lgrp - grow].N;
matnum += (tmp * ratio);
num += tmp;
}
if (isZero(num)) {
//no fish grow to this length cell
(*v[i])[lgrp].setToZero();
Mat->storeMatureStock(area, age, lgrp, 0.0, 0.0);
} else if (isZero(matnum)) {
//none of the fish that grow to this length cell mature
(*v[i])[lgrp].W = Weight[lgrp];
(*v[i])[lgrp].N = num;
Mat->storeMatureStock(area, age, lgrp, 0.0, 0.0);
} else if (isEqual(num, matnum) || (matnum > num)) {
//all the fish that grow to this length cell mature
(*v[i])[lgrp].setToZero();
Mat->storeMatureStock(area, age, lgrp, num, Weight[lgrp]);
} else {
(*v[i])[lgrp].W = Weight[lgrp];
(*v[i])[lgrp].N = num - matnum;
Mat->storeMatureStock(area, age, lgrp, matnum, Weight[lgrp]);
}
}
}
}