3a4,6
> #include <string>
> #include <valarray>
> #include <vector>
6c9
< #include <cmath> // fabs
---
> using namespace std;
8c11,12
< #include "FssimApp.h"
---
> #include "RegenStruct.h"
> #include "GroupDefStruct.h"
9a14,16
> #include "InvalidDataException.h"
> 
> 
12,17c19,20
<     int	func_agearea_index(int age, int NumAgeAreaElements, int period);
<     int	func_upperage_from_index(int index, int period);
<     int	func_group_membership(GROUP *group, int period, int when, char *zone, char *au, char *block, int age, int width, int ctype);
<     double	func_vph(TABLE *table, int age);
<     double	func_vph_recovered(TIMBER_CLASS *timber_class);
<     double	func_area_lt(double *agearea, /*int NumAgeAreaElements, */ int age, int period);
---
>     int	func_agearea_index(int age, int NumAgeAreaElements);
>     int	func_upperage_from_index(int index);
18a22,23
>     double func_area_lt(const valarray<double>* pAgeArea, int age);
> 
20,22c25
<     double	func_area_ge(double *agearea, int NumAgeAreaElements,  int age, int period/*, int width*/);
<     int	func_derive_integer(char *s);
<     double	func_derive_float(char *s);
---
>     double	func_area_ge(const valarray<double>* pAgeArea, int NumAgeAreaElements,  int age/*, int width*/);
24,25c27,28
<     ZONE 	*func_find_zone(char *label);
<     //ZONE 	*func_find_zone(int id);
---
>     int	func_derive_integer(char* s);
>     double	func_derive_float(char* s);
27,28c30
<     GROUP	*func_find_group(char *label);
<     GROUP	*func_find_group(int id);
---
>     bool func_compare_float(float f1, float f2);
30,31c32
<     TABLE	*func_find_table(char *label);
<     TABLE	*func_find_table(int id);
---
> } // end namespace tsb
33,42c34,39
<     AU* func_find_au(ZONE *zone, char *label);
<     //AU* func_find_au(ZONE *zone, int id);
< 
<     BLOCK	*func_find_block(char *label);
<     //BLOCK	*func_find_block(int id);
< 
<     int 	func_compare_float(float f1, float f2);
< 
<     // inlined
<     inline int func_compare_float(float f1, float f2)
---
> // inlined for speed
> inline bool tsb::func_compare_float(float f1, float f2)
> {
>     double tmp = fabs(f1-f2);
>     //tmp = static_cast<double>((f1-f2)); // MATH MISTAKE!!! (intentional)
>     if ( tmp < FssimUtil::m_kRoundError )
44,53c41
<         float tmp;
<         // cast from float to double
<         tmp = static_cast<double>(fabs(f1-f2));
<         //tmp = static_cast<double>((f1-f2)); // MATH MISTAKE!!! (intentional)
<         if ( tmp < FssimApp::m_kRoundError )
<         {
<             return(true);
<         } else {
<             return(false);
<         }
---
>         return true;
55,56c43
< 
<     inline BLOCK* func_find_block(char *label)
---
>     else
58,81c45
<         BLOCK	*block;
<         HASHREC	*hashrec;
<         FssimApp* pApp = FssimApp::getApp();
< 
<         pApp->m_blockList = (LIST_PTR *)List::reset(pApp->m_blockList);
<         hashrec = (HASHREC *)List::next(pApp->m_blockList);
<         while (hashrec != 0)
<         {
<             if (hashrec->id == static_cast<int>(label[0]))
<             {
<                 hashrec->candidates = (LIST_PTR *)List::reset(hashrec->candidates);
<                 block = (BLOCK *)List::next(hashrec->candidates);
<                 while (block != 0)
<                 {
<                     if (strcmp(block->label, label) == 0)
<                     {
<                         return block;
<                     }
<                     block = (BLOCK *)List::next(hashrec->candidates);
<                 }
<             }
<             hashrec = (HASHREC *)List::next(pApp->m_blockList);
<         }
<         return 0;
---
>         return false;
82a47
> }
84c49,52
< /*    BLOCK* func_find_block(int id)
---
> inline int tsb::func_derive_integer(char *s)
> {
>     int	ret;
>     for(int i = 0 ; i<(int)strlen(s) ; ++i)
86,92c54
<         BLOCK	*block;
<         HASHREC	*hashrec;
<         FssimApp* pApp = FssimApp::getApp();
< 
<         pApp->m_blockList = (LIST_PTR *)List::reset(pApp->m_blockList);
<         hashrec = (HASHREC *)List::next(pApp->m_blockList);
<         while (hashrec != 0)
---
>         if ( ( !isdigit(s[i]) ) && ( s[i]!='-' ) )
94,109c56
<             // error - the following line won't work without
<             // looking up label by block id (due to list hashing)
<             //if (hashrec->id == static_cast<int>(label[0]))
<             {
<                 hashrec->candidates = (LIST_PTR *)List::reset(hashrec->candidates);
<                 block = (BLOCK *)List::next(hashrec->candidates);
<                 while (block != 0)
<                 {
<                     if (block->id == id)
<                     {
<                         return block;
<                     }
<                     block = (BLOCK *)List::next(hashrec->candidates);
<                 }
<             }
<             hashrec = (HASHREC *)List::next(pApp->m_blockList);
---
>             return(FssimUtil::ERRORVAL);
111d57
<         return 0;
113c59,61
< */
---
>     ret = atoi(s);
>     return(ret);
> }
115c63,66
<     inline ZONE* func_find_zone(char *label)
---
> inline double tsb::func_derive_float(char *s)
> {
>     double	ret;
>     for(int i=0 ; i < static_cast<int>(strlen(s)) ; ++i)
117,121c68
<         ZONE* zone = 0;
<         FssimApp* pApp = FssimApp::getApp();
<         pApp->m_zoneList = (LIST_PTR *)List::reset(pApp->m_zoneList);
<         zone = (ZONE *)List::next(pApp->m_zoneList);
<         while (zone != 0)
---
>         if ((!isdigit(s[i])) && (s[i] != '.') && (s[i] != '-'))
123,127c70
<             if (strcmp(zone->label, label)==0)
<             {
<                 return zone;
<             }
<             zone = (ZONE *)List::next(pApp->m_zoneList);
---
>             return(FssimUtil::ERRORVAL);
129d71
<         return 0;
130a73,75
>     ret = atof(s);
>     return(ret);
> }
132,149c77,81
< /*
<     ZONE* func_find_zone(int id)
<     {
<         ZONE* zone = 0;
<         FssimApp* pApp = FssimApp::getApp();
<         pApp->m_zoneList = (LIST_PTR *)List::reset(pApp->m_zoneList);
<         zone = (ZONE *)List::next(pApp->m_zoneList);
<         while (zone != 0)
<         {
<             if (zone->id == id)
<             {
<                 return zone;
<             }
<             zone = (ZONE *)List::next(pApp->m_zoneList);
<         }
<         return 0;
<     }
< */
---
> inline int tsb::func_agearea_index(int age, int NumAgeAreaElements)
> {
>     int	index = 0;
>     FssimApp* pApp = FssimApp::getApp();
>     int periodLength = pApp->m_periodLength[pApp->m_currentPeriod];
151c83,86
<     inline GROUP* func_find_group(char *label)
---
>     // if age is negative or the last (0th) year of a pApp->m_currentPeriod,
>     // assign to agearea bin of preceding years
>     // e.g., for a 10-yr pApp->m_currentPeriod, yr 80 is assigned to the same agearea bin (index 9) as yrs 71 through 79
>     if ( (age % periodLength == 0) || (age < 0) )
153,166c88
<         GROUP* group;
<         group = 0;
<         FssimApp* pApp = FssimApp::getApp();
<         pApp->m_groupList = (LIST_PTR *)List::reset(pApp->m_groupList);
<         group = (GROUP *)List::next(pApp->m_groupList);
<         while (group != 0)
<         {
<             if (strcmp(group->label, label) == 0)
<             {
<                 return group;
<             }
<             group = (GROUP *)List::next(pApp->m_groupList);
<         }
<         return 0;
---
>         index = (age / periodLength) + 1 ;
168,169c90
< 
<     inline GROUP* func_find_group(int id)
---
>     else // if age is between the 1st and 2nd-last year of a period, inclusive
171,184c92
<         GROUP* group;
<         group = 0;
<         FssimApp* pApp = FssimApp::getApp();
<         pApp->m_groupList = (LIST_PTR *)List::reset(pApp->m_groupList);
<         group = (GROUP *)List::next(pApp->m_groupList);
<         while (group != 0)
<         {
<             if (group->id == id)
<             {
<                 return group;
<             }
<             group = (GROUP *)List::next(pApp->m_groupList);
<         }
<         return 0;
---
>         index = (age / periodLength) + 2 ;
186,201c94,98
< 
<     inline TABLE* func_find_table(char *label)
<     {
<         TABLE	*table;
<         table = 0;
<         FssimApp* pApp = FssimApp::getApp();
<         pApp->m_tableList = (LIST_PTR *)List::reset(pApp->m_tableList);
<         table = (TABLE *)List::next(pApp->m_tableList);
<         while (table != 0)
<         {
<             if (strcmp(table->label, label) == 0)
<             {
<                 return table;
<             }
<             table = (TABLE *)List::next(pApp->m_tableList);
<         }
---
>     if (index < 0)  // if age was originally more negative than the period length, assign
>                     // to same agearea bin as years -((3*periodLength)-1) through -(periodLength)
>                     // e.g., for a 10-year period, years -29 through -10 are assigned to bin index 0 (zero)
>                     // along with years -30 and more negative which produce negative indices
>     {               // in the above statement
204,205c101
< 
<     inline TABLE* func_find_table(int id)
---
>     if (index >= NumAgeAreaElements) // assign all really old timber classes to last bin of agearea array
207,220c103
<         TABLE	*table;
<         table = 0;
<         FssimApp* pApp = FssimApp::getApp();
<         pApp->m_tableList = (LIST_PTR *)List::reset(pApp->m_tableList);
<         table = (TABLE *)List::next(pApp->m_tableList);
<         while (table != 0)
<         {
<             if (table->id == id)
<             {
<                 return table;
<             }
<             table = (TABLE *)List::next(pApp->m_tableList);
<         }
<         return 0;
---
>         return (NumAgeAreaElements - 1);
221a105,106
>     return index;
> }
223c108,121
<     inline AU* func_find_au(ZONE *zone, char *label)
---
> inline int tsb::func_upperage_from_index(int index)
> {
> 	FssimApp* pApp = FssimApp::getApp();
>     int upperAge = (index - 1) * ( pApp->m_periodLength[pApp->m_currentPeriod] );
>     return upperAge;
> }
> 
> inline double tsb::func_area_lt(const valarray<double>* pAgeArea, /*int NumAgeAreaElements,*/ int age)
> {
>     FssimApp* pApp = FssimApp::getApp();
>     int period = pApp->m_currentPeriod;
> 
>     // returns area < age (LESS THAN)
>     if (pAgeArea == 0)
225,236c123
<         AU* au;
<         zone->aulist = (LIST_PTR *)List::reset(zone->aulist);
<         au = static_cast<AU*>(List::next(zone->aulist));
<         while (au != 0)
<         {
<             if (strcmp(au->label, label) == 0)
<             {
<                 return au;
<             }
<             au = static_cast<AU*>(List::next(zone->aulist));
<         }
<         return 0;
---
>     	return 0;
238,267c125
< /*
<     AU* func_find_au(ZONE *zone, int id)
<     {
<         AU* au;
<         zone->aulist = (LIST_PTR *)List::reset(zone->aulist);
<         au = static_cast<AU*>(List::next(zone->aulist));
<         while (au != 0)
<         {
<             if (au->id == id)
<             {
<                 return au;
<             }
<             au = static_cast<AU*>(List::next(zone->aulist));
<         }
<         return 0;
<     }
< */
<     inline int func_derive_integer(char *s)
<     {
<         int	ret;
<         for(int i = 0 ; i<(int)strlen(s) ; ++i)
<         {
<             if ( ( !isdigit(s[i]) ) && ( s[i]!='-' ) )
<             {
<                 return(FssimApp::ERRORVAL);
<             }
<         }
<         ret = atoi(s);
<         return(ret);
<     }
---
>     const valarray<double>& ageArea = *pAgeArea;
269c127
<     inline double func_derive_float(char *s)
---
>     if ( age <= -(pApp->m_periodLength[period]) )
271,280c129
<         double	ret;
<         for(int i=0 ; i < static_cast<int>(strlen(s)) ; ++i)
<         {
<             if ((!isdigit(s[i])) && (s[i] != '.') && (s[i] != '-'))
<             {
<                 return(FssimApp::ERRORVAL);
<             }
<         }
<         ret = atof(s);
<         return(ret);
---
>         return ageArea[0];
283,286c132,133
<     inline int func_agearea_index(int age, int NumAgeAreaElements, int period)
<     {
<         int	index = 0;
<         FssimApp* pApp = FssimApp::getApp();
---
>     float x = static_cast<float>(age - 1 + pApp->m_periodLength[period]) /
>               static_cast<float>(pApp->m_periodLength[period]);	// proportion
288,311c135,137
<         // if age is negative or the last (0th) year of a period,
<         // assign to agearea bin of preceding years
<         // e.g., for a 10-yr period, yr 80 is assigned to the same agearea bin (index 9) as yrs 71 through 79
<         if ( ((age % (pApp->m_periodLength[period])) == 0) || (age < 0) )
<         {
<             index = ( static_cast<int>( (age / (pApp->m_periodLength[period]) ) ) + 1 );
<         }
<         else // if age is between the 1st and 2nd-last year of a period, inclusive
<         {
<             index = ( static_cast<int>( (age / (pApp->m_periodLength[period]) ) ) + 2 );
<         }
<         if (index < 0)  // if age was originally more negative than the period length, assign
<                         // to same agearea bin as years -((3*periodLength)-1) through -(periodLength)
<                         // e.g., for a 10-year period, years -29 through -10 are assigned to bin index 0 (zero)
<                         // along with years -30 and more negative which produce negative indices
<         {               // in the above statement
<             index = 0;
<         }
<         if (index >= NumAgeAreaElements) // assign all really old timber classes to last bin of agearea array
<         {
<             index = NumAgeAreaElements - 1;
<         }
<         return index;
<     }
---
>     int index = static_cast<int>(x); 	// integer part
>     float y = x - index;		// float part
>     double result = 0.0;
313c139,141
<     inline int func_upperage_from_index(int index, int period)
---
>     // vectorize pAgeArea as valarray to replace for loop
> 
>     for(int i = 0 ; i <= index ; i++)
315,317c143
<     	FssimApp* pApp = FssimApp::getApp();
<         int upperAge = (index - 1) * ( pApp->m_periodLength[period] );
<         return upperAge;
---
>         result += ageArea[i];
320,322c146,152
<     inline double func_area_lt(double *agearea, /*int NumAgeAreaElements,*/ int age, int period)
<     {
<         FssimApp* pApp = FssimApp::getApp();
---
>     // turned off vectorization for debugging
>     // take slice of AgeArea valarray from (index + 1) to NumAgeAreaElements
>     /*
>     const int sliceStart = 0;
>     const int sliceSize = index + 1;
>     valarray<double> ageAreaSlice(0.0, sliceSize);
>     ageAreaSlice = (*pAgeArea)[slice(sliceStart, sliceSize, 1)];
324,328c154,155
<         /* returns area < age (LESS THAN) */
<         if (agearea == 0)
<         {
<         	return 0.0;
<         }
---
>     result += ageAreaSlice.sum(); // pAgeArea is a valarray
>     */
330,333c157
<         if ( age <= -(pApp->m_periodLength[period]) )
<         {
<             return agearea[0];
<         }
---
>     result += ageArea[index + 1] * y;
335,336c159,160
<         float x = static_cast<float>(age - 1 + pApp->m_periodLength[period]) /
<                   static_cast<float>(pApp->m_periodLength[period]);	// proportion
---
>     return result;
> }
338,347c162,165
<         int index = static_cast<int>(x); 	/* integer part */
<         float y = x - index;		/* float part */
<         double result = 0.0;
<         for(int i = 0 ; i <= index ; i++)
<         {
<             result += agearea[i];
<         }
<         result += agearea[index + 1] * y;
<         return result;
<     }
---
> // width added in v4.0.1 to account for width = 1 at time zero
> inline double tsb::func_area_ge(const valarray<double>* pAgeArea, int NumAgeAreaElements, int age)
> {
>     FssimApp* pApp = FssimApp::getApp();
349,350c167,170
<     // width added in v4.0.1 to account for width = 1 at time zero
<     inline double func_area_ge(double* agearea, int NumAgeAreaElements, int age, int period/*, int width*/)
---
>     int period = pApp->m_currentPeriod;
> 
>     // returns area >= age (GREATER/EQUAL)
>     if (pAgeArea == 0)
352c172,174
<         FssimApp* pApp = FssimApp::getApp();
---
>     	return 0;
>     }
>     const valarray<double>& ageArea = *pAgeArea;
354,358c176,180
<         // returns area >= age (GREATER/EQUAL)
<         if (agearea == 0)
<         {
<             return 0.0;
<         }
---
>     // calculate the percentage of the age area elements at the constraint age
>     // which will have partially reached the given age by the end of the period
>     int thisPeriodLength = pApp->m_periodLength[period];
>     float factor = static_cast<float>( (thisPeriodLength - (age % thisPeriodLength) + 1) % thisPeriodLength ) /
>                    static_cast<float>(thisPeriodLength);
360,364d181
<         // calculate the percentage of the age area elements at the constraint age
<         // which will have partially reached the given age by the end of the period
<         int thisPeriodLength = pApp->m_periodLength[period];
<         float factor = static_cast<float>( (thisPeriodLength - (age % thisPeriodLength) + 1) % thisPeriodLength ) /
<                        static_cast<float>(thisPeriodLength);
366d182
< 
368,370c184,185
<         cout << endl << setw(43) << "factor before adjustment = "
<              << setprecision(2) << factor
<              /*<< setw(20) << "width = " << setw(3) << width*/ << endl;
---
>     cout << endl << setw(43) << "factor before adjustment = "
>          << setprecision(2) << factor << endl;
373,376c188,191
<         if (factor == 0.0)
<         {
<             factor = 1.0;
<         }
---
>     if (factor == 0)
>     {
>         factor = 1.0;
>     }
379,381c194,196
<         // pApp->m_applyFactor added in v4.0.2 to remove recursive factoring
<         if (!pApp->m_applyFactor)
<         {
---
>     // pApp->m_applyFactor added in v4.0.2 to remove recursive factoring
>     if (!pApp->m_applyFactor)
>     {
383c198
<             cout << "apply factor is false" << endl;
---
>         cout << "apply factor is false" << endl;
385,386c200,202
<             factor = 1.0;
<         }
---
>         factor = 1.0;
>     }
> 
389,390c205,206
<         double areaGreaterOrEqual = 0.0;
<         int index = func_agearea_index(age, NumAgeAreaElements, period);
---
>     double areaGreaterOrEqual = 0.0;
>     int index = func_agearea_index(age, NumAgeAreaElements);
392,393c208,209
<         // take only a percentage of the age area element which overlaps the contraint age
<         areaGreaterOrEqual += agearea[index] * factor;
---
>     // take only a percentage of the age area element which overlaps the contraint age
>     areaGreaterOrEqual += ageArea[index] * factor;
396,398c212,214
<         cout << setw(20) << "area_ge(before) = " << setw(8) << setprecision(8) << areaGreaterOrEqual
<              << setw(20) << "age = " << setprecision(3) << setw(3) << age
<              << setw(20) << "NumEl = " << setw(3) << NumAgeAreaElements << endl;
---
>     cout << setw(20) << "area_ge(before) = " << setw(8) << setprecision(8) << areaGreaterOrEqual
>          << setw(20) << "age = " << setprecision(3) << setw(3) << age
>          << setw(20) << "NumEl = " << setw(3) << NumAgeAreaElements << endl;
401,405c217,223
<         // sum up all age area elements older than the constraint age
<         for( int i = index + 1 ; i < NumAgeAreaElements ; ++i)
<         {
<             areaGreaterOrEqual += agearea[i];
<         }
---
>     // turned off vectorization for debugging
>     // take slice of AgeArea valarray from (index + 1) to NumAgeAreaElements
>     /*
>     const int sliceStart = index + 1;
>     const int sliceSize = NumAgeAreaElements - sliceStart - 1;
>     valarray<double> ageAreaSlice(0.0, sliceSize);
>     ageAreaSlice = ageArea[slice(sliceStart, sliceSize, 1)];
407,412c225,227
< #ifdef _TSB_GEFACTOR_DEBUG
<         cout << setw(20) << "areaGreaterOrEqual = " << setw(8) << setprecision(8) << areaGreaterOrEqual
<              << setw(20) << "index = " << setprecision(3) << setw(3) << index
<              << setw(20) << "factor = " << setw(3) << factor
<              /*<< setw(20) << "width = " << setw(3) << width*/ << endl;
< #endif
---
>     // sum up all age area elements older than the constraint age
>     areaGreaterOrEqual += ageAreaSlice.sum(); // pAgeArea is a valarray
>     */
414,415c229
<         return areaGreaterOrEqual;
<     }
---
>     // vectorized pAgeArea as valarray to replace for loop
417c231
<     inline int func_group_membership(GROUP *group, int period, int when, char *zone, char *au, char *block, int age, int width, int ctype)
---
>     for( int i = index + 1 ; i < NumAgeAreaElements ; ++i)
419,507c233
<         GROUPDEF	*gd;
<         int		offset;
<         FssimApp* pApp = FssimApp::getApp();
< 
<         if ( when == FssimApp::BEFORE )
<         {
<             offset = pApp->m_periodLength[period];
<         } else {
<             offset = 0;
<         }
< 
<         /*
<         Iterate through all groupdefs records, looking for a match
<         If ages match; eight cases:   Zone      AU      Block
<         1)	-	-        -
<         1)	-	-        z
<         2)	-	y        -
<         2)	-	y        z
<         3)	x	-        -
<         3)	x	-        z
<         4)	x	y        -
<         4)	x	y        z
<         */
<         group->groupdefs = (LIST_PTR *)List::reset(group->groupdefs);
<         gd = (GROUPDEF *)List::next(group->groupdefs);
<         while (gd != 0)
<         {
<             if ((group->grouptype==FssimApp::STATIC_GROUP)&&(period==1)&&(when==FssimApp::BEFORE))
<             {
<             /*
<             FssimApp::STATIC_GROUP groups do not use look-ahead for membership;
<             assess membership at time 0 - beginning of planning horizon
<             before any aging takes place
<                 */
<                 if ((age>=gd->low_age)&&(age<=gd->high_age)||
<                     ((age-width+1)>=gd->low_age)&&((age-width+1)<=gd->high_age)||
<                     ((age-width+1)<=gd->low_age)&&(age>=gd->high_age))
<                 {
<                     if ((strcmp(gd->block, FssimApp::m_kNullString)==0)||(strcmp(gd->block,block)==0))
<                     {
<                         if ((strcmp(gd->zone,FssimApp::m_kNullString)==0) && (strcmp(gd->au,FssimApp::m_kNullString)==0))
<                         {
<                             if ((gd->class_type==FssimApp::DUMMYINPUT)||(ctype==gd->class_type)) return(true); /* case 1 */
<                         } else if ((strcmp(gd->zone,FssimApp::m_kNullString)==0) && (strcmp(gd->au,au)==0))
<                         {
<                             if ((gd->class_type==FssimApp::DUMMYINPUT)||(ctype==gd->class_type)) return(true); /* case 2 */
<                         } else if ((strcmp(gd->zone,zone)==0) && (strcmp(gd->au,FssimApp::m_kNullString)==0))
<                         {
<                             if ((gd->class_type==FssimApp::DUMMYINPUT)||(ctype==gd->class_type)) return(true); /* case 3 */
<                         } else if ((strcmp(gd->zone,zone)==0) && (strcmp(gd->au,au)==0))
<                         {
<                             if ((gd->class_type==FssimApp::DUMMYINPUT)||(ctype==gd->class_type)) return(true); /* case 4 */
<                         }
<                     }
<                 }
<             } else if (group->grouptype==FssimApp::FLOAT_GROUP)
<             {
<             /*
<             FssimApp::FLOAT_GROUP and REPORT type groups use look-ahead for membership;
<             if timber_class ages will age into group or
<             timber_class ages out of group during period, consider
<             timber_class to be member
<                 */
<                 if (((age+offset)<gd->low_age)||((age-width+1)>gd->high_age))
<                 {
<                     /* NON MEMBER */
<                 } else {
<                     /* MEMBER */
<                     if ((strcmp(gd->block, FssimApp::m_kNullString)==0)||(strcmp(gd->block,block)==0))
<                     {
<                         if ((strcmp(gd->zone,FssimApp::m_kNullString)==0) && (strcmp(gd->au,FssimApp::m_kNullString)==0))
<                         {
<                             if ((gd->class_type==FssimApp::DUMMYINPUT)||(ctype==gd->class_type)) return(true); /* case 1 */
<                         } else if ((strcmp(gd->zone,FssimApp::m_kNullString)==0) && (strcmp(gd->au,au)==0))
<                         {
<                             if ((gd->class_type==FssimApp::DUMMYINPUT)||(ctype==gd->class_type)) return(true); /* case 2 */
<                         } else if ((strcmp(gd->zone,zone)==0) && (strcmp(gd->au,FssimApp::m_kNullString)==0))
<                         {
<                             if ((gd->class_type==FssimApp::DUMMYINPUT)||(ctype==gd->class_type)) return(true); /* case 3 */
<                         } else if ((strcmp(gd->zone,zone)==0) && (strcmp(gd->au,au)==0))
<                         {
<                             if ((gd->class_type==FssimApp::DUMMYINPUT)||(ctype==gd->class_type)) return(true); /* case 4 */
<                         }
<                     }
<                 }
<             }
<             gd = (GROUPDEF *)List::next(group->groupdefs);
<         }
<         return(false);
---
>         areaGreaterOrEqual += ageArea[i];
510,512d235
<     inline double func_vph(TABLE *table, int age)
<     {
<         double result = 0.0;
514,548c237,240
<         // OAFs are applied to tables in data_implement_table_oafs()
<         // No need to apply OAFs here; straight table lookup
<         if (age < 0)
<         {
<             result = 0.0;
<         }
<         else if (age >= FssimApp::MAX_NUM_TABLE_VERTICES)
<         {
<             result = table->Volume[FssimApp::MAX_NUM_TABLE_VERTICES-1];
<         }
<         else
<         {
<             result = table->Volume[age];
<         }
<         return (result);
<     }
< 
<     inline double func_vph_recovered(TIMBER_CLASS *timber_class)
<     {
<         REGEN* regen = 0;
<         double	vph = 0.0;
<         double regenvph = 0.0;
<         double result = 0.0;
<         int	new_age = 0;
< 
<         if (timber_class == 0)
<         {
<             printf("\n\nERROR: func_vph_recovered() passed 0 timber_class parameter\n\n");
<             exit(1);
<         }
< 
<         vph = func_vph(timber_class->parent_au->table, timber_class->age);
< 
< #ifdef _TSB_HARVEST_DEBUG
<     cout << endl << "func_vph_recovered: vph: " << vph << endl;
---
> #ifdef _TSB_GEFACTOR_DEBUG
>     cout << setw(20) << "areaGreaterOrEqual = " << setw(8) << setprecision(8) << areaGreaterOrEqual
>          << setw(20) << "index = " << setprecision(3) << setw(3) << index
>          << setw(20) << "factor = " << setw(3) << factor << endl;
551,617c243
<         if (timber_class->parent_au->regenlist == 0)
<         {
<             regenvph = func_vph(timber_class->parent_au->table, 0);
< 
< #ifdef _TSB_HARVEST_DEBUG
<     	    cout << endl << "func_vph_recovered: regenvph: " << regenvph << endl;
< #endif
< 
<         }
<         else
<         {
<             regenvph = 0.0;
<             timber_class->parent_au->regenlist = (LIST_PTR *)List::reset(timber_class->parent_au->regenlist);
<             regen = (REGEN *)List::next(timber_class->parent_au->regenlist);
<             while (regen != 0)
<             {
<                 /* may regenerate to original (or any) age */
<                 if (regen->ToAge == -FssimApp::INFINITY)
<                 {
<                     new_age = timber_class->age + regen->Offset;
<                 }
<                 else
<                 {
<                     new_age = regen->ToAge + regen->Offset;
<                 }
< #ifdef _TSB_REGEN_DEBUG
<             	cout << endl << "func_vph_recovered: new age: " << new_age << endl;
<     			cout << "func_vph_recovered: regen percent : " << regen->percent << endl;
<     			cout << "func_vph_recovered: regen ToAge : " << regen->ToAge << endl;
<     			cout << "func_vph_recovered: regen Offset : " << regen->Offset << endl;
< #endif
<                 if (0 != regen->au && 0 != regen->au->table)
<                 {
<                     regenvph += func_vph(regen->au->table, new_age)*
<                                 (regen->percent/100.0);
<                     regen = (REGEN *)List::next(timber_class->parent_au->regenlist);
< 
< #ifdef _TSB_HARVEST_DEBUG
<     	        cout << endl << "func_vph_recovered (2): regenvph: " << regenvph << endl;
< #endif
< 
<                 }
<                 else
<                 {
<                     string message = "func_vph_recovered: regen->au or regen->au->table is null. Au label: ";
<                     message += timber_class->parent_au->label;
<                     throw message;
<                 }
<             }
<         }
<         result = vph - regenvph;
< #ifdef _TSB_FUNC_VPH_DEBUG
<         cout << endl << "func_vph_recovered (end): vph-regenvph: " << result << endl;
< #endif
<         if (result<0.0)
<         {
<             printf("\n\nERROR: func_vph_recovered() returned NEGATIVE m3/ha" );
<             printf("\n\t          Zone: %s", timber_class->parent_zone->label );
<             printf("\n\t            AU: %s", timber_class->parent_au->table->label );
<             printf("\n\t           Age: %d", timber_class->age );
<             printf("\n\t Pre-Regen VPH: %.2f", vph );
<             printf("\n\tPost-Regen VPH: %.2f", regenvph );
<             printf("\n\t Resultant VPH: %.2f\n\n", result );
<             exit(1);
<         }
<         return(result);
<     }
---
>     return areaGreaterOrEqual;