1 /***************************************************************************
2  *   Copyright (C) 2006 by BUI Quang Minh, Steffen Klaere, Arndt von Haeseler   *
3  *   minh.bui@univie.ac.at   *
4  *                                                                         *
5  *   This program is free software; you can redistribute it and/or modify  *
6  *   it under the terms of the GNU General Public License as published by  *
7  *   the Free Software Foundation; either version 2 of the License, or     *
8  *   (at your option) any later version.                                   *
9  *                                                                         *
10  *   This program is distributed in the hope that it will be useful,       *
11  *   but WITHOUT ANY WARRANTY; without even the implied warranty of        *
12  *   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the         *
13  *   GNU General Public License for more details.                          *
14  *                                                                         *
15  *   You should have received a copy of the GNU General Public License     *
16  *   along with this program; if not, write to the                         *
17  *   Free Software Foundation, Inc.,                                       *
18  *   59 Temple Place - Suite 330, Boston, MA  02111-1307, USA.             *
19  ***************************************************************************/
20 #include "mpdablock.h"
21 #include "pda/split.h"
22 #include "pda/splitgraph.h"
23 
MPdaBlock(SplitGraph * asgraph)24 MPdaBlock::MPdaBlock(SplitGraph *asgraph)
25  : NxsBlock()
26 {
27 	budget = -1;
28 	min_budget = -1;
29 	sub_size = 0;
30 	cost_constrained = false;
31 	id = "PDA";
32 	sgraph = asgraph;
33 }
34 
35 
Report(ostream & out)36 void MPdaBlock::Report(ostream &out)
37 {
38 	out << "Budget = " << budget << endl;
39 	out << "Taxa Costs = ";
40 	for (DoubleVector::iterator it = costs.begin(); it != costs.end(); it++)
41 		out << *it << " ";
42 	out << endl;
43 }
44 
Reset()45 void MPdaBlock::Reset()
46 {
47 	errormsg.clear();
48 	isEmpty			= true;
49 	isEnabled		= true;
50 	isUserSupplied	= false;
51 
52 }
53 
Read(NxsToken & token)54 void MPdaBlock::Read(NxsToken &token)
55 {
56 
57 	int ntax = sgraph->getNTaxa();
58 	if (ntax <= 0) {
59 		errormsg = "PDA Block should be preceded by Splits Block";
60 		throw NxsException(errormsg, token.GetFilePosition(), token.GetFileLine(), token.GetFileColumn());
61 
62 	}
63 
64 	//int nominal_ntax	= 0;
65 	//int nominal_nsplits = 0;
66 
67 	// This should be the semicolon after the block name
68 	//
69 	token.GetNextToken();
70 
71 	if (!token.Equals(";"))
72 	{
73 		errormsg = "Expecting ';' after PDA block name, but found ";
74 		errormsg += token.GetToken();
75 		errormsg += " instead";
76 		throw NxsException(errormsg, token.GetFilePosition(), token.GetFileLine(), token.GetFileColumn());
77 	}
78 
79 	for (;;)
80 	{
81 		token.GetNextToken();
82 
83 		if (token.Equals("PARAMETERS"))
84 		{
85 			// This should be the NTAX keyword
86 			//
87 			token.GetNextToken();
88 
89 			do {
90 
91 				if (token.Equals("BUDGET")) {
92 					// This should be the equals sign
93 					//
94 					token.GetNextToken();
95 
96 					if (!token.Equals("="))
97 					{
98 						errormsg = "Expecting '=', but found ";
99 						errormsg += token.GetToken();
100 						errormsg += " instead";
101 						throw NxsException(errormsg, token.GetFilePosition(), token.GetFileLine(), token.GetFileColumn());
102 					}
103 
104 					// This should be the integer budget
105 					//
106 					token.GetNextToken();
107 
108 					budget = convert_double(token.GetToken().c_str());
109 					if (budget <= 0)
110 					{
111 						errormsg = "BUDGET should be greater than zero (";
112 						errormsg += token.GetToken();
113 						errormsg += " was specified)";
114 						throw NxsException(errormsg, token.GetFilePosition(), token.GetFileLine(), token.GetFileColumn());
115 					}
116 				} else if (token.Equals("MIN BUDGET")) {
117 					// This should be the equals sign
118 					//
119 					token.GetNextToken();
120 
121 					if (!token.Equals("="))
122 					{
123 						errormsg = "Expecting '=', but found ";
124 						errormsg += token.GetToken();
125 						errormsg += " instead";
126 						throw NxsException(errormsg, token.GetFilePosition(), token.GetFileLine(), token.GetFileColumn());
127 					}
128 
129 					// This should be the integer budget
130 					//
131 					token.GetNextToken();
132 
133 					min_budget = convert_double(token.GetToken().c_str());
134 					if (budget < 0)
135 					{
136 						errormsg = "MIN_BUDGET should be greater than or equal to zero (";
137 						errormsg += token.GetToken();
138 						errormsg += " was specified)";
139 						throw NxsException(errormsg, token.GetFilePosition(), token.GetFileLine(), token.GetFileColumn());
140 					}
141 
142 				} else if (token.Equals("BUDGET CONSTRAINED")) {
143 					cost_constrained = true;
144 				} else if (token.Equals("K")) {
145 					// This should be the equals sign
146 					//
147 					token.GetNextToken();
148 
149 					if (!token.Equals("="))
150 					{
151 						errormsg = "Expecting '=', but found ";
152 						errormsg += token.GetToken();
153 						errormsg += " instead";
154 						throw NxsException(errormsg, token.GetFilePosition(), token.GetFileLine(), token.GetFileColumn());
155 					}
156 
157 					// This should be the integer budget
158 					//
159 					token.GetNextToken();
160 
161 					sub_size = atoi(token.GetToken().c_str());
162 					if (sub_size <= 1)
163 					{
164 						errormsg = "K should be greater than 1 (";
165 						errormsg += token.GetToken();
166 						errormsg += " was specified)";
167 						throw NxsException(errormsg, token.GetFilePosition(), token.GetFileLine(), token.GetFileColumn());
168 					}
169 				} else
170 				{
171 					errormsg = "Invalid PARAMETERS command: ";
172 					errormsg += token.GetToken();
173 					throw NxsException(errormsg, token.GetFilePosition(), token.GetFileLine(), token.GetFileColumn());
174 				}
175 
176 				token.GetNextToken();
177 
178 			} while (!token.AtEOF() && !token.Equals(";"));
179 
180 			if (!token.Equals(";"))
181 			{
182 				errormsg = "Expecting ';' to terminate PARAMETERS command, but found ";
183 				errormsg += token.GetToken();
184 				errormsg += " instead";
185 				throw NxsException(errormsg, token.GetFilePosition(), token.GetFileLine(), token.GetFileColumn());
186 			}
187 
188 		} // if (token.Equals("PARAMETERS"))
189 
190 		else if (token.Equals("TAXCOSTS")) {
191 
192 			costs.resize(ntax, -1);
193 			// This should be taxon name
194 			//
195 			token.GetNextToken();
196 
197 			do {
198 				int tax_id = -1;
199 
200 				try {
201 					tax_id = sgraph->getTaxa()->FindTaxon(token.GetToken());
202 				} catch (NxsTaxaBlock::NxsX_NoSuchTaxon) {
203 					tax_id = -1;
204 				}
205 
206 				if (tax_id < 0)
207 				{
208 					errormsg = "Taxon is not found (";
209 					errormsg += token.GetToken();
210 					errormsg += " was specified)";
211 					throw NxsException(errormsg, token.GetFilePosition(), token.GetFileLine(), token.GetFileColumn());
212 				}
213 
214 				// This should be the cost of taxon
215 				//
216 				token.GetNextToken();
217 
218 				int taxcost = convert_double(token.GetToken().c_str());
219 				if (taxcost < 0)
220 				{
221 					errormsg = "Taxon cost should be greater than or equal to zero (";
222 					errormsg += token.GetToken();
223 					errormsg += " was specified)";
224 					throw NxsException(errormsg, token.GetFilePosition(), token.GetFileLine(), token.GetFileColumn());
225 				}
226 				costs[tax_id] = taxcost;
227 
228 				token.GetNextToken();
229 			} while (!token.AtEOF() && !token.Equals(";"));
230 
231 			// This should be the terminating semicolon
232 			//
233 			//token.GetNextToken();
234 
235 			if (!token.Equals(";"))
236 			{
237 				errormsg = "Expecting ';' to terminate TAXCOSTS command, but found ";
238 				errormsg += token.GetToken();
239 				errormsg += " instead";
240 				throw NxsException(errormsg, token.GetFilePosition(), token.GetFileLine(), token.GetFileColumn());
241 			}
242 
243 			for (int i = 0; i < ntax; i++)
244 				if (costs[i] < 0) {
245 					costs[i] = 0;
246 					cout << "WARNING: taxon " << sgraph->getTaxa()->GetTaxonLabel(i)
247 						<< "has no cost! set to 0." << endl;
248 				}
249 		}	// if (token.Equals("TAXCOSTS"))
250 
251 
252 		else if (token.Equals("END") || token.Equals("ENDBLOCK"))
253 		{
254 			// Get the semicolon following END
255 			//
256 			token.GetNextToken();
257 
258 			if (!token.Equals(";"))
259 			{
260 				errormsg = "Expecting ';' to terminate the ENDBLOCK command, but found ";
261 				errormsg += token.GetToken();
262 				errormsg += " instead";
263 				throw NxsException(errormsg, token.GetFilePosition(), token.GetFileLine(), token.GetFileColumn());
264 			}
265 			break;
266 		}	// if (token.Equals("END") || token.Equals("ENDBLOCK"))
267 
268 		else
269 		{
270 			SkippingCommand(token.GetToken());
271 			do
272 			{
273 				token.GetNextToken();
274 			}
275 			while (!token.AtEOF() && !token.Equals(";"));
276 
277 			if (token.AtEOF())
278 			{
279 				errormsg = "Unexpected end of file encountered";
280 				throw NxsException(errormsg, token.GetFilePosition(), token.GetFileLine(), token.GetFileColumn());
281 			}
282 		}	// token not END, ENDBLOCK, COST
283 	}	// GetNextToken loop
284 
285 }
286 
readBudgetFile(Params & params)287 void MPdaBlock::readBudgetFile(Params &params) {
288 	ifstream in;
289 	in.exceptions(ios::failbit | ios::badbit);
290 	cout << "Reading budget information file " << params.budget_file << "..." << endl;
291 	NxsString taxname;
292 	int ntaxa = sgraph->getNTaxa() - params.is_rooted;
293 	int i;
294 
295 	try {
296 		costs.resize(ntaxa, -1);
297 		in.open(params.budget_file);
298 		in >> budget;
299 		if (budget < 0)
300 			throw "Negative total budget.";
301 		for (i = 0; i < ntaxa && !in.eof(); i++) {
302 			double taxcost;
303 			int tax_id = -1;
304 			taxname = "";
305 			in.exceptions(ios::badbit);
306 			in >> taxname;
307 			in.exceptions(ios::failbit | ios::badbit);
308 			if (taxname == "") break;
309 			in >> taxcost;
310 			if (taxcost < 0)
311 				throw "Negative taxa preservation cost.";
312 			tax_id = sgraph->getTaxa()->FindTaxon(taxname);
313 			costs[tax_id] = taxcost;
314 		}
315 		in.close();
316 	} catch (ios::failure) {
317 		outError(ERR_READ_INPUT);
318 	} catch (NxsTaxaBlock::NxsX_NoSuchTaxon) {
319 		outError(ERR_NO_TAXON, taxname);
320 	} catch (const char *str) {
321 		outError(str);
322 	} catch (...) {
323 		// anything else
324 		outError(ERR_READ_ANY);
325 	}
326 
327 	for (i = 0; i < ntaxa; i++)
328 		if (costs[i] < 0) {
329 			costs[i] = 0;
330 			cout << "WARNING: taxon " << sgraph->getTaxa()->GetTaxonLabel(i)
331 				<< "has no cost! set to 0." << endl;
332 		}
333 	cost_constrained = true;
334 }
335 
readBudgetAreaFile(Params & params)336 void MPdaBlock::readBudgetAreaFile(Params &params) {
337 	ifstream in;
338 	in.exceptions(ios::failbit | ios::badbit);
339 	cout << "Reading budget for areas information file " << params.budget_file << "..." << endl;
340 	string areaname;
341 	int nareas = sgraph->getNAreas();
342 	int i;
343 
344 	try {
345 		costs.resize(nareas, -1);
346 		in.open(params.budget_file);
347 		in >> budget;
348 		if (budget < 0)
349 			throw "Negative total budget.";
350 		for (i = 0; i < nareas && !in.eof(); i++) {
351 			double areacost;
352 			int area_id = -1;
353 			areaname = "";
354 			in.exceptions(ios::badbit);
355 			in >> areaname;
356 			in.exceptions(ios::failbit | ios::badbit);
357 			if (areaname == "") break;
358 			in >> areacost;
359 			if (areacost < 0)
360 				throw "Negative taxa preservation cost.";
361 			area_id = sgraph->getSetsBlock()->findArea(areaname);
362 			if (area_id < 0)
363 				outError(ERR_NO_AREA, areaname);
364 			costs[area_id] = areacost;
365 		}
366 		in.close();
367 	} catch (ios::failure) {
368 		outError(ERR_READ_INPUT);
369 	} catch (const char *str) {
370 		outError(str);
371 	} catch (...) {
372 		// anything else
373 		outError(ERR_READ_ANY);
374 	}
375 
376 	for (i = 0; i < nareas; i++)
377 		if (costs[i] < 0) {
378 			costs[i] = 0;
379 			cout << "WARNING: area " << sgraph->getSetsBlock()->getSet(i)->name
380 				<< "has no cost! set to 0." << endl;
381 		}
382 	cost_constrained = true;
383 }
384 
385 
SkippingCommand(NxsString commandName)386 void MPdaBlock::SkippingCommand(NxsString commandName) {
387 	cout << "   Skipping unknown command (" << commandName << ")..." << endl;
388 }
389 
390 
~MPdaBlock()391 MPdaBlock::~MPdaBlock()
392 {
393 }
394