#include "tags.h"
#include "areatime.h"
#include "errorhandler.h"
#include "readfunc.h"
#include "readword.h"
#include "stock.h"
#include "stockprey.h"
#include "gadget.h"
#include "global.h"
Tags::Tags(CommentStream& infile, const char* givenname, const AreaClass* const Area,
const TimeClass* const TimeInfo, Keeper* const keeper, StockPtrVector stockvec)
: HasName(givenname) {
taggingstock = 0;
numtagtimesteps = 0;
char text[MaxStrLength];
strncpy(text, "", MaxStrLength);
ifstream subfile;
CommentStream subcomment(subfile);
keeper->addString("tags");
keeper->addString(givenname);
//Currently can only have one stock per tagging experiment
readWordAndValue(infile, "stock", text);
stocknames.resize(new char[strlen(text) + 1]);
strcpy(stocknames[0], text);
//Currently can only have one area per tagging experiment
readWordAndVariable(infile, "tagarea", tagarea);
tagarea = Area->getInnerArea(tagarea);
infile >> ws;
char c = infile.peek();
if ((c == 'e') || (c == 'E'))
readWordAndVariable(infile, "endyear", endyear);
else
endyear = TimeInfo->getLastYear();
int i, found = 0;
for (i = 0; i < stockvec.Size(); i++) {
if (strcasecmp(stockvec[i]->getName(), stocknames[0]) == 0) {
if (found == 0) {
taggingstock = stockvec[i];
LgrpDiv = new LengthGroupDivision(*(taggingstock->getLengthGroupDiv()));
if (LgrpDiv->Error())
handle.logMessage(LOGFAIL, "Error in tags - failed to create length group");
tagStocks.resize(taggingstock);
taggingstock->setTaggedStock();
}
found++;
}
}
if (found == 0)
handle.logMessage(LOGFAIL, "Error in tags - failed to match stock", stocknames[0]);
if (found > 1)
handle.logMessage(LOGFAIL, "Error in tags - repeated stock", stocknames[0]);
//Now read in the tagloss information
readWordAndVariable(infile, "tagloss", tagloss);
tagloss.Inform(keeper);
//read in the numbers format: tagid - length - number
readWordAndValue(infile, "numbers", text);
subfile.open(text, ios::in);
handle.checkIfFailure(subfile, text);
handle.Open(text);
readNumbers(subcomment, givenname, TimeInfo);
handle.Close();
subfile.close();
subfile.clear();
keeper->clearLast();
keeper->clearLast();
}
void Tags::readNumbers(CommentStream& infile, const char* tagname, const TimeClass* const TimeInfo) {
int year, step, count, reject;
int i, lenid, keepdata, timeid;
int numlen = LgrpDiv->numLengthGroups();
double tmplength, tmpnumber;
char tmpname[MaxStrLength];
strncpy(tmpname, "", MaxStrLength);
infile >> ws;
//Check the number of columns in the inputfile
if (countColumns(infile) != 5)
handle.logFileMessage(LOGFAIL, "wrong number of columns in inputfile - should be 5");
year = step = count = reject = 0;
while (!infile.eof()) {
keepdata = 1;
infile >> tmpname >> year >> step >> tmplength >> tmpnumber >> ws;
//only keep the data if tmpname matches tagname
if (strcasecmp(tagname, tmpname) != 0)
keepdata = 0;
//only keep the data if the length is valid
lenid = -1;
for (i = 0; i < numlen; i++)
if (isEqual(tmplength, LgrpDiv->minLength(i)))
lenid = i;
//OK the length doesnt match a minimum length so find the length group it is in
if ((lenid == -1) && (tmplength > LgrpDiv->minLength()) && (tmplength < LgrpDiv->maxLength())) {
for (i = 1; i < numlen; i++) {
if (tmplength < LgrpDiv->minLength(i)) {
lenid = i - 1;
break;
}
}
if (lenid == -1)
lenid = numlen - 1; //then this must be the last length group
}
if (lenid == -1)
keepdata = 0;
//only keep the data if the number is positive
if (tmpnumber < 0.0) {
handle.logMessage(LOGWARN, "Warning in tags - found negative number of tags", tmpnumber);
keepdata = 0;
}
timeid = -1;
if ((TimeInfo->isWithinPeriod(year, step)) && (keepdata == 1)) {
for (i = 0; i < Years.Size(); i++)
if ((Years[i] == year) && (Steps[i] == step))
timeid = i;
if (timeid == -1) {
Years.resize(1, year);
Steps.resize(1, step);
timeid = Years.Size() - 1;
NumberByLength.resize(new DoubleMatrix(1, numlen, 0.0));
}
} else
keepdata = 0;
if (keepdata == 1) {
count++;
(*NumberByLength[timeid])[0][lenid] += tmpnumber;
} else
reject++; //count number of rejected data points read from file
}
if (count == 0)
handle.logMessage(LOGWARN, "Warning in tags - found no data in the data file for", tagname);
if (reject != 0)
handle.logMessage(LOGMESSAGE, "Discarded invalid tags data - number of invalid entries", reject);
handle.logMessage(LOGMESSAGE, "Read tags data file - number of entries", count);
tagyear = 9999;
tagstep = 9999;
for (i = 0; i < Years.Size(); i++)
if ((Years[i] < tagyear) || (Years[i] == tagyear && Steps[i] < tagstep)) {
tagyear = Years[i];
tagstep = Steps[i];
}
timeid = -1;
for (i = 0; i < Years.Size(); i++)
if ((Years[i] == tagyear) && (Steps[i] == tagstep))
timeid = i;
if (timeid == -1)
handle.logMessage(LOGFAIL, "Error in tags - calculated invalid timestep");
numtagtimesteps = (TimeInfo->numSteps() * (endyear - tagyear)) - tagstep + 1;
if (endyear == TimeInfo->getLastYear())
numtagtimesteps += TimeInfo->getLastStep();
else
numtagtimesteps += TimeInfo->numSteps();
}
Tags::~Tags() {
int i;
for (i = 0; i < stocknames.Size(); i++)
delete[] stocknames[i];
for (i = 0; i < NumberByLength.Size(); i++)
delete NumberByLength[i];
for (i = 0; i < CI.Size(); i++)
delete CI[i];
while (AgeLengthStock.Nrow() > 0)
AgeLengthStock.Delete(0);
while (NumBeforeEating.Nrow() > 0)
NumBeforeEating.Delete(0);
delete LgrpDiv;
}
void Tags::Reset() {
int i;
while (AgeLengthStock.Nrow() > 0)
AgeLengthStock.Delete(0);
while (NumBeforeEating.Nrow() > 0)
NumBeforeEating.Delete(0);
while (CI.Size() > 0) {
delete CI[0];
CI.Delete(0);
}
for (i = 0; i < preyindex.Size(); i++)
preyindex[i] = -1;
for (i = 0; i < updated.Size(); i++)
updated[i] = 0;
}
void Tags::setStock(StockPtrVector& Stocks) {
int i, j, found;
StockPtrVector tmpStockVector;
char* stockname;
preyindex.resize(1, -1);
updated.resize(1, 0);
if (!(taggingstock->isInArea(tagarea)))
handle.logMessage(LOGFAIL, "Error in tags - stock isnt defined on tagging area");
if (taggingstock->doesMove()) {
tmpStockVector = taggingstock->getTransitionStocks();
for (i = 0; i < tmpStockVector.Size(); i++) {
transitionStocks.resize(tmpStockVector[i]);
preyindex.resize(1, -1);
updated.resize(1, 0);
tagStocks.resize(tmpStockVector[i]);
tmpStockVector[i]->setTaggedStock();
}
}
if (taggingstock->doesMature()) {
tmpStockVector = taggingstock->getMatureStocks();
for (i = 0; i < tmpStockVector.Size(); i++) {
matureStocks.resize(tmpStockVector[i]);
found = 0;
for (j = 0; j < transitionStocks.Size(); j++)
if (strcasecmp(transitionStocks[j]->getName(), tmpStockVector[i]->getName()) != 0)
found++;
if (found == 0) {
preyindex.resize(1, -1);
updated.resize(1, 0);
tagStocks.resize(tmpStockVector[i]);
tmpStockVector[i]->setTaggedStock();
}
}
}
if (taggingstock->doesStray()) {
tmpStockVector = taggingstock->getStrayStocks();
for (i = 0; i < tmpStockVector.Size(); i++) {
strayStocks.resize(tmpStockVector[i]);
found = 0;
for (j = 0; j < transitionStocks.Size(); j++)
if (strcasecmp(transitionStocks[j]->getName(), tmpStockVector[i]->getName()) != 0)
found++;
for (j = 0; j < matureStocks.Size(); j++)
if (strcasecmp(matureStocks[j]->getName(), tmpStockVector[i]->getName()) != 0)
found++;
if (found == 0) {
preyindex.resize(1, -1);
updated.resize(1, 0);
tagStocks.resize(tmpStockVector[i]);
tmpStockVector[i]->setTaggedStock();
}
}
}
for (i = 1; i < tagStocks.Size(); i++) {
stockname = new char[strlen(tagStocks[i]->getName()) + 1];
strcpy(stockname, tagStocks[i]->getName());
stocknames.resize(stockname);
}
}
//Must have set stocks according to stocknames using setStock before calling Update()
//Now we need to distribute the tagged fish to the same age/length groups as the tagged stock.
void Tags::Update(int timeid) {
int i, j;
PopInfoVector NumberInArea;
NumberInArea.resizeBlank(LgrpDiv->numLengthGroups());
const AgeBandMatrix* stockPopInArea;
const LengthGroupDivision* tmpLgrpDiv;
stockPopInArea = &(taggingstock->getCurrentALK(tagarea));
stockPopInArea->sumColumns(NumberInArea);
//Now we have total number of stock per length in tagarea, N(., l) (NumberInArea) and
//number of stock per age/length, N(a, l) (stockPopInArea) so we must initialise
//AgeLengthStock so that it can hold all information of number of tagged stock
//per area/age/length after endtime. We must make AgeBandMatrixPtrVector same size as
//the one in stock even though have only one area entry at the beginning
IntVector stockareas = taggingstock->getAreas();
int numareas = stockareas.Size();
int tagareaindex = -1;
j = 0;
while (j <= numareas && tagareaindex == -1) {
if (tagarea == stockareas[j])
tagareaindex = j;
j++;
}
if (tagareaindex == -1)
handle.logMessage(LOGFAIL, "Error in tags - invalid area for tagged stock");
int maxage = stockPopInArea->maxAge();
int minage = stockPopInArea->minAge();
int numberofagegroups = maxage - minage + 1;
int upperlgrp, minl, maxl, age, length, stockid;
double numfishinarea, numstockinarea;
IntVector lgrpsize(numberofagegroups, 0);
IntVector lgrpmin(numberofagegroups, 0);
for (i = 0; i < numberofagegroups; i++) {
lgrpmin[i]= stockPopInArea->minLength(i + minage);
upperlgrp = stockPopInArea->maxLength(i + minage);
lgrpsize[i] = upperlgrp - lgrpmin[i];
}
AgeLengthStock.resize(new AgeBandMatrixPtrVector(numareas, minage, lgrpmin, lgrpsize));
for (age = minage; age <= maxage; age++) {
minl = stockPopInArea->minLength(age);
maxl = stockPopInArea->maxLength(age);
for (length = minl; length < maxl; length++) {
numfishinarea = NumberInArea[length].N;
numstockinarea = (*stockPopInArea)[age][length].N;
if (numfishinarea > verysmall && numstockinarea > verysmall)
(*AgeLengthStock[0])[tagareaindex][age][length].N = (*NumberByLength[timeid])[0][length - minl] * numstockinarea / numfishinarea;
else
(*AgeLengthStock[0])[tagareaindex][age][length].N = 0.0;
}
}
taggingstock->addTags(AgeLengthStock[0], this, exp(-tagloss));
updated[0] = 1;
if (taggingstock->isEaten()) {
tmpLgrpDiv = taggingstock->getPrey()->getLengthGroupDiv();
lgrpmin.Reset();
lgrpsize.Reset();
lgrpmin.resize(numberofagegroups, 0);
lgrpsize.resize(numberofagegroups, tmpLgrpDiv->numLengthGroups());
NumBeforeEating.resize(new AgeBandMatrixPtrVector(numareas, minage, lgrpmin, lgrpsize));
CI.resize(new ConversionIndex(LgrpDiv, tmpLgrpDiv));
if (CI[CI.Size() - 1]->Error())
handle.logMessage(LOGFAIL, "Error in tags - error when checking length structure");
stockid = stockIndex(taggingstock->getName());
if (stockid < 0 || stockid >= preyindex.Size())
handle.logMessage(LOGFAIL, "Error in tags - invalid stock identifier");
preyindex[stockid] = NumBeforeEating.Nrow() - 1;
}
for (i = 1; i < tagStocks.Size(); i++) {
stockPopInArea = &tagStocks[i]->getCurrentALK(tagarea);
stockareas = tagStocks[i]->getAreas();
numareas = stockareas.Size();
maxage = stockPopInArea->maxAge();
minage = stockPopInArea->minAge();
numberofagegroups = maxage - minage + 1;
lgrpmin.Reset();
lgrpsize.Reset();
lgrpmin.resize(numberofagegroups, 0);
lgrpsize.resize(numberofagegroups, 0);
for (j = 0; j < numberofagegroups; j++) {
lgrpmin[j] = stockPopInArea->minLength(j + minage);
upperlgrp = stockPopInArea->maxLength(j + minage);
lgrpsize[j] = upperlgrp - lgrpmin[j];
}
AgeLengthStock.resize(new AgeBandMatrixPtrVector(numareas, minage, lgrpmin, lgrpsize));
if (tagStocks[i]->isEaten()) {
tmpLgrpDiv = tagStocks[i]->getPrey()->getLengthGroupDiv();
lgrpmin.Reset();
lgrpsize.Reset();
lgrpmin.resize(numberofagegroups, 0);
lgrpsize.resize(numberofagegroups, tmpLgrpDiv->numLengthGroups());
NumBeforeEating.resize(new AgeBandMatrixPtrVector(numareas, minage, lgrpmin, lgrpsize));
CI.resize(new ConversionIndex(LgrpDiv, tmpLgrpDiv));
if (CI[CI.Size() - 1]->Error())
handle.logMessage(LOGFAIL, "Error in tags - error when checking length structure");
stockid = stockIndex(tagStocks[i]->getName());
if (stockid < 0 || stockid >= preyindex.Size())
handle.logMessage(LOGFAIL, "Error in tags - invalid stock identifier");
preyindex[stockid] = NumBeforeEating.Nrow() - 1;
}
}
}
void Tags::updateTags(int year, int step) {
int i, timeid;
timeid = -1;
for (i = 0; i < Years.Size(); i++)
if (Years[i] == year && Steps[i] == step)
timeid = i;
if (timeid != -1) {
if (tagyear == year && tagstep == step)
this->Update(timeid);
else
this->addToTagStock(timeid);
}
}
void Tags::deleteStockTags() {
int i;
for (i = 0; i < tagStocks.Size(); i++) {
if (updated[i] == 1) {
tagStocks[i]->deleteTags(this->getName());
updated[i] = 2;
}
}
}
void Tags::updateMatureStock(const TimeClass* const TimeInfo) {
int i, id;
if (endyear <= TimeInfo->getYear())
handle.logMessage(LOGWARN, "Warning in tags - tagging experiment has finished");
else {
for (i = 0; i < matureStocks.Size(); i++) {
id = stockIndex(matureStocks[i]->getName());
if (id < 0 || id >= AgeLengthStock.Nrow())
handle.logMessage(LOGFAIL, "Error in tags - invalid stock identifier");
if (updated[id] == 0) {
matureStocks[i]->addTags(AgeLengthStock[id], this, exp(-tagloss));
updated[id] = 1;
}
}
}
}
void Tags::updateTransitionStock(const TimeClass* const TimeInfo) {
int i, id;
if (endyear <= TimeInfo->getYear())
handle.logMessage(LOGWARN, "Warning in tags - tagging experiment has finished");
else {
for (i = 0; i < transitionStocks.Size(); i++) {
id = stockIndex(transitionStocks[i]->getName());
if (id < 0 || id >= AgeLengthStock.Nrow())
handle.logMessage(LOGFAIL, "Error in tags - invalid stock identifier");
if (updated[id] == 0) {
transitionStocks[i]->addTags(AgeLengthStock[id], this, exp(-tagloss));
updated[id] = 1;
}
}
}
}
void Tags::updateStrayStock(const TimeClass* const TimeInfo) {
int i, id;
if (endyear <= TimeInfo->getYear())
handle.logMessage(LOGWARN, "Warning in tags - tagging experiment has finished");
else {
for (i = 0; i < strayStocks.Size(); i++) {
id = stockIndex(strayStocks[i]->getName());
if (id < 0 || id >= AgeLengthStock.Nrow())
handle.logMessage(LOGFAIL, "Error in tags - invalid stock identifier");
if (updated[id] == 0) {
strayStocks[i]->addTags(AgeLengthStock[id], this, exp(-tagloss));
updated[id] = 1;
}
}
}
}
int Tags::stockIndex(const char* stockname) {
int i;
for (i = 0; i < tagStocks.Size(); i++)
if (strcasecmp(stockname, tagStocks[i]->getName()) == 0)
return i;
return -1;
}
int Tags::areaIndex(const char* stockname, int area) {
int i, j;
for (i = 0; i < tagStocks.Size(); i++) {
if (strcasecmp(stockname, tagStocks[i]->getName()) == 0) {
IntVector stockareas = tagStocks[i]->getAreas();
for (j = 0; j < stockareas.Size(); j++)
if (stockareas[j] == area)
return j;
return -1;
}
}
return -1;
}
void Tags::storeConsumptionALK(int area, const char* stockname) {
int stockid, preyid, areaid;
stockid = stockIndex(stockname);
if (stockid < 0)
handle.logMessage(LOGFAIL, "Error in tags - invalid stock identifier");
preyid = preyindex[stockid];
if (preyid > NumBeforeEating.Nrow() || preyid < 0)
handle.logMessage(LOGFAIL, "Error in tags - invalid prey identifier");
areaid = areaIndex(stockname, area);
if (areaid < 0)
handle.logMessage(LOGFAIL, "Error in tags - invalid area identifier");
(*NumBeforeEating[preyid])[areaid].setToZero();
(*NumBeforeEating[preyid])[areaid].Add((*AgeLengthStock[stockid])[areaid], *CI[preyid]);
}
const AgeBandMatrix& Tags::getConsumptionALK(int area, const char* stockname) {
int stockid, preyid, areaid;
stockid = stockIndex(stockname);
if (stockid < 0)
handle.logMessage(LOGFAIL, "Error in tags - invalid stock identifier");
preyid = preyindex[stockid];
if (preyid > NumBeforeEating.Nrow() || preyid < 0)
handle.logMessage(LOGFAIL, "Error in tags - invalid prey identifier");
areaid = areaIndex(stockname, area);
if (areaid < 0)
handle.logMessage(LOGFAIL, "Error in tags - invalid area identifier");
return (*NumBeforeEating[preyid])[areaid];
}
int Tags::isWithinPeriod(int year, int step) {
if ((year > tagyear || (year == tagyear && step >= tagstep)) && (year <= endyear))
return 1;
return 0;
}
void Tags::addToTagStock(int timeid) {
int i, tagareaindex;
PopInfoVector NumberInArea;
NumberInArea.resizeBlank(LgrpDiv->numLengthGroups());
const AgeBandMatrix* stockPopInArea;
stockPopInArea = &(taggingstock->getCurrentALK(tagarea));
stockPopInArea->sumColumns(NumberInArea);
i = 0;
tagareaindex = -1;
IntVector stockareas = taggingstock->getAreas();
while (i <= stockareas.Size() && tagareaindex == -1) {
if (tagarea == stockareas[i])
tagareaindex = i;
i++;
}
if (tagareaindex == -1)
handle.logMessage(LOGFAIL, "Error in tags - invalid area for tagged stock");
int maxage = stockPopInArea->maxAge();
int minage = stockPopInArea->minAge();
int minl, maxl, age, length;
double numfishinarea, numstockinarea;
for (age = minage; age <= maxage; age++) {
minl = stockPopInArea->minLength(age);
maxl = stockPopInArea->maxLength(age);
for (length = minl; length < maxl; length++) {
numfishinarea = NumberInArea[length].N;
numstockinarea = (*stockPopInArea)[age][length].N;
if (numfishinarea > verysmall && numstockinarea > verysmall)
(*AgeLengthStock[0])[tagareaindex][age][length].N += (*NumberByLength[timeid])[0][length - minl] * numstockinarea / numfishinarea;
}
}
}
void Tags::Print(ofstream& outfile) const {
// outfile << "\nTag\nName" << sep << this->getName() << endl;
}