/*
* This file is part of rasdaman community.
*
* Rasdaman community is free software: you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation, either version 3 of the License, or
* (at your option) any later version.
*
* Rasdaman community is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with rasdaman community. If not, see .
*
* Copyright 2003, 2004, 2005, 2006, 2007, 2008, 2009 Peter Baumann /
rasdaman GmbH.
*
* For more information please see
* or contact Peter Baumann via .
/
/**
* SOURCE: interesttiling.cc
*
* MODULE: rasodmg
* CLASS: r_Interest_Tiling
*
* COMMENTS:
* None
*/
#include "rasodmg/interesttiling.hh"
#include "rasodmg/alignedtiling.hh"
#include "rasodmg/dirdecompose.hh"
#include "rasodmg/dirtiling.hh"
#include "raslib/rminit.hh"
#include
#include
#include
#include
#ifdef AIX
#include
#endif
// Uncoment the _VISUALIZE_2D_DECOMP_ line to generate ppm files the
// visualization of the domain decomposition done by the algoritm
// #define _VISUALIZE_2D_DECOMP_
// Uncoment the folling line for some debug printfs on this class
// #define _DEBUG_INTERESTTILING_
#ifdef _VISUALIZE_2D_DECOMP_
#include "tools/visualtiling2d.hh"
#endif
// This is a data structure for internal use within interesttiling.cc.
// It is defined here because of problems with ptrepository.
// -------------
// This structure consists of a r_Minterval and a counter. Also, it has
// some auxiliary constructors.
class Classified_Block
{
public:
int intersection_count; // Intersection count
r_Minterval block; // The actual block
// Default constructor
Classified_Block(int count = 0);
// Constructor with the actual block and the counter
Classified_Block(const r_Minterval& b, int count = 0);
// Same data structure (this operator is needed because of dlist)
bool operator==(const Classified_Block& other) const;
// Different data structure (this operator is needed because of dlist)
bool operator!=(const Classified_Block& other) const;
// Friend of std::ostream (this operator is needed because of dlist)
friend std::ostream& operator<<(std::ostream& os, const Classified_Block& block);
};
const char*
r_Interest_Tiling::description = "dimensions, areas of interest, tile size (in bytes) and tile size limit [NOLIMIT|REGROUP|SUBTILING|REGROUPSUBTILING] (ex: \"2;[0:9,0:9];[100:109,0:9];100;REGROUPSUBTILING\")";
const char* r_Interest_Tiling::tilesizelimit_name_nolimit ="NOLIMIT";
const char* r_Interest_Tiling::tilesizelimit_name_regroup ="REGROUP";
const char* r_Interest_Tiling::tilesizelimit_name_subtiling ="SUBTILING";
const char* r_Interest_Tiling::tilesizelimit_name_regroupandsubtiling ="REGROUPSUBTILING";
const char* r_Interest_Tiling::all_tilesizelimit_names[r_Interest_Tiling::NUMBER]={
tilesizelimit_name_nolimit,
tilesizelimit_name_regroup,
tilesizelimit_name_subtiling,
tilesizelimit_name_regroupandsubtiling
};
r_Interest_Tiling::Tilesize_Limit
r_Interest_Tiling::get_tilesize_limit_from_name(const char* name)
{
if(!name) {
RMInit::logOut << "r_Interest_Tiling::get_tilesize_limit_from_name(" << (name?name: "NULL") << ")" << endl;
return r_Interest_Tiling::NUMBER;
}
unsigned int i=r_Interest_Tiling::NUMBER;
for (i=0; i<(unsigned int)r_Interest_Tiling::NUMBER; i++)
{
if (strcasecmp(name, all_tilesizelimit_names[i]) == 0)
break;
}
return (r_Interest_Tiling::Tilesize_Limit)i;
}
const char*
r_Interest_Tiling::get_name_from_tilesize_limit(Tilesize_Limit tsl)
{
static const char* unknown="UNKNOWN";
unsigned int idx = (unsigned int)tsl;
if (idx >= (unsigned int)r_Interest_Tiling::NUMBER)
return unknown;
return all_tilesizelimit_names[idx];
}
r_Interest_Tiling::r_Interest_Tiling(const char* encoded) throw (r_Error)
: r_Dimension_Tiling(0,0)
{
if(!encoded)
{
RMInit::logOut << "r_Interest_Tiling::r_Interest_Tiling(" << (encoded?encoded: "NULL") << ")" << endl;
throw r_Error(TILINGPARAMETERNOTCORRECT);
}
std::vector vectInterestAreas;
r_Interest_Tiling::Tilesize_Limit tileSizeLimit;
r_Dimension tileD=0;
r_Bytes tileS=0, lenToConvert=0;
const char *pStart=NULL, *pRes=NULL, *pTemp=NULL, *pEnd=NULL;
char *pToConvert=NULL;
pStart=encoded;
pEnd=pStart+strlen(pStart);
pTemp=pStart;
pRes=strstr(pTemp,TCOLON);
if(!pRes)
{
RMInit::logOut << "r_Interest_Tiling::r_Interest_Tiling(" << encoded << "): Error decoding tile dimension from \"" << pTemp << "\"." << endl;
throw r_Error(TILINGPARAMETERNOTCORRECT);
}
//deal with dimension
lenToConvert=pRes-pTemp;
pToConvert=new char[lenToConvert+1];
memcpy(pToConvert,pTemp, lenToConvert);
pToConvert[lenToConvert]='\0';
tileD=strtol(pToConvert, (char**)NULL, DefaultBase);
if (!tileD)
{
RMInit::logOut << "r_Interest_Tiling::r_Interest_Tiling(" << encoded << "): Error decoding tile dimension from \"" << pToConvert << "\"." << endl;
delete[] pToConvert;
throw r_Error(TILINGPARAMETERNOTCORRECT);
}
if (tileD < 0)
{
RMInit::logOut << "r_Interest_Tiling::r_Interest_Tiling(" << encoded << "): Error decoding tile dimension from \"" << pToConvert << "\", is negative." << endl;
delete[] pToConvert;
throw r_Error(TILINGPARAMETERNOTCORRECT);
}
//skip COLON && free buffer
delete[] pToConvert;
if(pRes != (pEnd-1))
pRes++;
else
{
RMInit::logOut << "r_Interest_Tiling::r_Interest_Tiling(" << encoded << "): Error decoding interest areas, end of stream." << endl;
throw r_Error(TILINGPARAMETERNOTCORRECT);
}
//parse interest areas
pTemp=pRes;
pRes=strstr(pTemp,TCOLON);
if(!pRes)
{
RMInit::logOut << "r_Interest_Tiling::r_Interest_Tiling(" << encoded << "): Error decoding interest areas." << endl;
throw r_Error(TILINGPARAMETERNOTCORRECT);
}
while(pRes)
{
//is interest areas?
if(*pTemp != *LSQRBRA)
break;
//copy parsed interest area
lenToConvert=pRes-pTemp;
pToConvert=new char[lenToConvert+1];
memcpy(pToConvert, pTemp, lenToConvert);
pToConvert[lenToConvert]='\0';
//add to vector
try
{
r_Minterval a(pToConvert);
vectInterestAreas.push_back(a);
}
catch(r_Error& err)
{
RMInit::logOut << "r_Interest_Tiling::r_Interest_Tiling(" << encoded << "): Error decoding interest area from \"" << pToConvert << "\"." << endl;
RMInit::logOut << "Error " << err.get_errorno() << " : " << err.what() << endl;
delete [] pToConvert;
throw r_Error(TILINGPARAMETERNOTCORRECT);
}
//skip COLON
delete[] pToConvert;
if(pRes != (pEnd-1))
pRes++;
else
{
RMInit::logOut << "r_Interest_Tiling::r_Interest_Tiling(" << encoded << "): Error decoding interest areas, end of stream." << endl;
throw r_Error(TILINGPARAMETERNOTCORRECT);
}
//try next item
pTemp=pRes;
pRes=strstr(pTemp, TCOLON);
}
if(vectInterestAreas.empty())
{
RMInit::logOut << "r_Interest_Tiling::r_Interest_Tiling(" << encoded << "): Error decoding interest areas, no interest areas specified." << endl;
throw r_Error(TILINGPARAMETERNOTCORRECT);
}
//deal with tile size
lenToConvert=pRes-pTemp;
pToConvert=new char[lenToConvert+1];
memcpy(pToConvert, pTemp, lenToConvert);
pToConvert[lenToConvert]='\0';
tileS=strtol(pToConvert, (char**)NULL, DefaultBase);
if (!tileS)
{
RMInit::logOut << "r_Interest_Tiling::r_Interest_Tiling(" << encoded << "): Error decoding tile size from \"" << pToConvert << "\"." << endl;
delete[] pToConvert;
throw r_Error(TILINGPARAMETERNOTCORRECT);
}
if (tileS < 0)
{
RMInit::logOut << "r_Interest_Tiling::r_Interest_Tiling(" << encoded << "): Error decoding tile size from \"" << pToConvert << "\", is negative." << endl;
delete[] pToConvert;
throw r_Error(TILINGPARAMETERNOTCORRECT);
}
//skip COLON && free buffer
delete[] pToConvert;
if(pRes != (pEnd-1))
pRes++;
else
{
RMInit::logOut << "r_Interest_Tiling::r_Interest_Tiling(" << encoded << "): Error decoding tile size limit, end of stream." << endl;
throw r_Error(TILINGPARAMETERNOTCORRECT);
}
pTemp=pRes;
tileSizeLimit=r_Interest_Tiling::get_tilesize_limit_from_name(pTemp);
if(tileSizeLimit==r_Interest_Tiling::NUMBER)
{
RMInit::logOut << "r_Interest_Tiling::r_Interest_Tiling(" << encoded << "): Error decoding tile size limit from \"" << pTemp << "\"." << endl;
throw r_Error(TILINGPARAMETERNOTCORRECT);
}
iareas = vectInterestAreas;
ts_strat = tileSizeLimit;
dimension = tileD;
tile_size = tileS;
}
r_Interest_Tiling::r_Interest_Tiling(r_Dimension dim, const std::vector& interest_areas, r_Bytes ts, Tilesize_Limit strat) throw (r_Error)
: r_Dimension_Tiling(dim, ts),
iareas(interest_areas),
ts_strat(strat)
{
for (std::vector::iterator it = iareas.begin(); it != iareas.end(); it++)
if (it->dimension() != dimension)
{
RMInit::logOut << "r_Interest_Tiling::r_Interest_Tiling(" << dim << ", " << interest_areas << ", " << ts << ", " << strat << ") the interest area domain " << *it << " does not match the dimension of this tiling scheme (" << dimension << ")" << endl;
throw r_Edim_mismatch(dimension, it->dimension());
}
}
r_Interest_Tiling::~r_Interest_Tiling()
{
}
r_Tiling* r_Interest_Tiling::clone() const
{
r_Tiling* copy = new r_Interest_Tiling(dimension, iareas, tile_size, ts_strat);
return copy;
}
void r_Interest_Tiling::print_status(std::ostream& os) const
{
os << "r_Interest_Tiling[ ";
r_Dimension_Tiling::print_status(os);
os << " interest areas = " << iareas << ", tiling strategy = " << ts_strat << " ]";
}
r_Tiling_Scheme
r_Interest_Tiling::get_tiling_scheme() const
{
return r_InterestTiling;
}
static int r_Range_comp(const void *elem1, const void *elem2)
{
r_Range e1 = *((r_Range*) elem1);
r_Range e2 = *((r_Range*) elem2);
if (e1 == e2)
return 0;
if (e1 < e2)
return -1;
else
return +1;
}
std::vector*
r_Interest_Tiling::make_partition(const r_Minterval& domain) const
{
r_Dimension dim = domain.dimension();
int total = 2 * iareas.size();
// We need one decomp from each dimension
std::vector* part = new std::vector(dim);
// We have at most (number of interest areas + 2) intervals
r_Range* intervals = new r_Range[total + 2];
// Create iterator for interest areas
std::vector::const_iterator it = iareas.begin();
int total_iareas = iareas.size();
// For all dimensions
for (r_Dimension i = 0; i < dim; i++)
{
it = iareas.begin(); // Reset iterator
intervals[0] = domain[i].low(); // Input lower domain limit
intervals[total+1] = domain[i].high(); // Input higher domain limit
for (int j = 1; j < total + 1; j += 2, ++it) // For all possible intervals
{
if ((*it)[i].low()-1 <= domain[i].low()) // Input low iarea limit
intervals[j] = domain[i].low();
else
intervals[j] = (*it)[i].low()-1;
intervals[j+1] = (*it)[i].high(); // Input higher iarea limit
}
// Sort the table
qsort((void*) intervals, total+2, sizeof(r_Range), r_Range_comp);
// Create partition using the limits table
for (int k=0; k*
r_Interest_Tiling::group(std::vector& blocks, r_Bytes typelen, Blocks_Type btype) const
{
r_Bytes tilesize = get_tile_size();
int joins = 0;
bool group_blocks = true;
// The list of threated blocks
std::vector* treated = new std::vector;
// An iterator for the blocks list
std::vector::iterator blocks_it = blocks.begin();
// The current block
r_Minterval current_block;
// For all the blocks in list
while (!blocks.empty())
{
// Get first block from list
current_block = blocks.back();
blocks.pop_back();
//this is neccessary when the compiler optimizes out the .end() check
unsigned int numberOfLevels = blocks.size();
for (blocks_it = blocks.begin(); blocks_it != blocks.end(); blocks_it++)
{
if (numberOfLevels == 0)
{
RMInit::logOut << "r_Interest_Tiling::group() the for loop was incorrectly optimized. breaking the loop." << endl;
break;
}
r_Minterval aux = *blocks_it;
// In principle two blocks can't be merged
group_blocks = false;
// If they can be merged
if (current_block.is_mergeable(aux))
{
std::vector::iterator ia_it = blocks.begin();
switch (btype)
{
case BLOCKS_A:
group_blocks = true;
// Check if the two blocks belong exaclty to the same iareas
for (; ia_it != blocks.end(); ia_it++)
{
if (aux.intersects_with(*ia_it) != current_block.intersects_with(*ia_it))
{
group_blocks = false;
break;
}
}
break;
case BLOCKS_B:
for (; ia_it != blocks.end(); ia_it++) // For all iareas
{
// Find the one this block intersects
if (current_block.intersects_with(*ia_it))
{
// Check if the other area intersects it
if (!aux.intersects_with(*ia_it))
group_blocks = false;
else
group_blocks = true;
break;
}
}
if (!group_blocks)
break;
group_blocks = false;
case BLOCKS_C: // Falls in (this is, also applies to B);
// Only on this two strategies, tilesize should be looked at
if ((ts_strat == REGROUP) || (ts_strat == REGROUP_AND_SUBTILING))
{
// If the resulting size isn't larger than tilesize
if ((current_block.cell_count()+aux.cell_count()) * typelen
< get_tile_size())
group_blocks = true;
}
else
group_blocks = true;
break;
}
}
// take care of the iterator advance, if is possible
if (group_blocks)
{
// Merge them
blocks.erase(blocks_it); // Automatically advances to next block
// take care of the size of the blocks
numberOfLevels--;
current_block.closure_with(aux);
++joins;
}
else
numberOfLevels--;
}
// Update the treated list with the current block
treated->push_back(current_block);
}
std::vector* result = 0;
// If there was joins, the algoritm must be repeted
if (joins != 0)
{
result = group(*treated, typelen, btype);
delete treated;
}
else
{
result = treated;
}
// Return the result
return result;
}
std::vector*
r_Interest_Tiling::compute_tiles(const r_Minterval& domain, r_Bytes typelen) const throw (r_Error)
{
r_Dimension num_dims = domain.dimension(); // Dimensionality of dom
if (domain.dimension() != dimension)
{
RMInit::logOut << "r_Interest_Tiling::compute_tiles(" << domain << ", " << typelen << ") dimension (" << dimension << ") does not match dimension of object to tile (" << num_dims << ")" << endl;
throw r_Edim_mismatch(dimension, num_dims);
}
if (typelen > tile_size)
{
RMInit::logOut << "r_Interest_Tiling::compute_tiles(" << domain << ", " << typelen << ") tile size (" << tile_size << ") is smaller than type length (" << typelen << ")" << endl;
throw r_Error(TILESIZETOOSMALL);
}
#ifdef _VISUALIZE_2D_DECOMP_ // User wants a visual
static int count; // Number of decomps
++count; // Update num decomps
// of the 2D decomp.
Visual_Tiling_2D* vis;
if (domain.dimension() == 2)
{
// Create an object for visualization
char fname[80];
sprintf(fname, "2D_decomp_int_%d.ppm", count);
vis = new Visual_Tiling_2D(domain, fname);
}
#endif
// *** Main algoritm ***
// The result
std::vector* result = new std::vector;
// Create a partition for dir tiling
std::vector* part = make_partition(domain);
// Perform dirtiling
r_Dir_Tiling *dir_tiling= NULL;
std::vector* dir_domain=NULL;
try
{
dir_tiling=new r_Dir_Tiling(num_dims, *part, tile_size, r_Dir_Tiling::WITHOUT_SUBTILING);
// Compute all tiles (directional tiles)
dir_domain = dir_tiling->compute_tiles(domain, typelen);
delete dir_tiling;
}
catch(r_Error& err)
{
delete result;
delete part;
if(dir_tiling)
delete dir_tiling;
throw;
}
// Get an interator for interest areas
//std::vector::iterator interest_area;
// Create a list for holding the classifed blocks
std::vector part_domain;
// Finds how many intersections exist between a block an the interest areas
for (std::vector::iterator dir_block = dir_domain->begin(); dir_block != dir_domain->end(); dir_block++)
{
Classified_Block b(*dir_block, 0);
for (std::vector::const_iterator interest_area = iareas.begin(); interest_area != iareas.end(); interest_area++)
{
if (b.block.intersects_with(*interest_area))
++b.intersection_count;
}
part_domain.push_back(b);
}
// Lists used for grouping blocks
std::vector Out;
std::vector In_Unique;
std::vector In_Common;
// Divide blocks into lists according to their number of intersections
for (std::vector::iterator class_block = part_domain.begin(); class_block != part_domain.end(); class_block++)
{
switch ((*class_block).intersection_count)
{
case 0:
Out.push_back((*class_block).block);
break;
case 1:
In_Unique.push_back((*class_block).block);
break;
default:
In_Common.push_back((*class_block).block);
break;
}
}
// Group blocks
std::vector* Blocks_A = group(In_Common, typelen, BLOCKS_A);
std::vector* Blocks_B = group(In_Unique, typelen, BLOCKS_B);
std::vector* Blocks_C = group(Out, typelen, BLOCKS_C);
std::vector::iterator it_A = Blocks_A->begin();
std::vector::iterator it_B = Blocks_B->begin();
std::vector::iterator it_C = Blocks_C->begin();
// If may be necessary to perform sub-tiling
if ((ts_strat == SUB_TILING) || (ts_strat == REGROUP_AND_SUBTILING))
{
// Variable to hold result of sub-tiling
std::vector* subtiles = 0;
// We need an array to hold the 3 iterators of the resulting blocks
std::vector::iterator* blocks_it[3];
std::vector* blocks_vec[3];
// Put iterators on the list
blocks_it[0] = &it_A;
blocks_it[1] = &it_B;
blocks_it[2] = &it_C;
blocks_vec[0] = Blocks_A;
blocks_vec[1] = Blocks_B;
blocks_vec[2] = Blocks_C;
// For all the lists (Blocs_A, Blocks_B and Blocks_C)
for (int j=0; j<3; j++)
{
// Reset iterator to the begin
(*blocks_it[j]) = blocks_vec[j]->begin();
// Tile each block if necessary
for (; *blocks_it[j] != blocks_vec[j]->end(); (*blocks_it[j])++)
{
if ((**blocks_it[j]).cell_count()*typelen > get_tile_size())
{
// Create a specification of a regular n-dim cube grid
r_Minterval specs(num_dims);
for (r_Dimension i=0; i::iterator it = subtiles->begin(); it != subtiles->end(); it++)
{
result->push_back(*it);
#ifdef _VISUALIZE_2D_DECOMP_
if (domain.dimension() == 2)
(*vis) << (*it);
#endif
}
delete subtiles;
}
else // No subtiling needed
{
// Insert block as it is
result->push_back(**blocks_it[j]);
#ifdef _VISUALIZE_2D_DECOMP_
if (domain.dimension() == 2)
(*vis) << (**blocks_it[j]);
#endif
}
}
}
#ifdef _VISUALIZE_2D_DECOMP_
if (domain.dimension() == 2)
{
vis->set_pen(255, 255, 0);
interest_area.seek_begin();
while (interest_area.not_done())
{
(*vis) << (*interest_area);
++interest_area;
}
}
#endif
}
else
{
// The result is just the sum of all blocks
result->insert(result->end(), Blocks_A->begin(), Blocks_A->end());
result->insert(result->end(), Blocks_B->begin(), Blocks_B->end());
result->insert(result->end(), Blocks_C->begin(), Blocks_C->end());
// Visualization
#ifdef _VISUALIZE_2D_DECOMP_
if (domain.dimension() == 2)
{
vis->set_pen(255, 255, 0);
for (it_A.seek_begin(); it_A.not_done(); it_A++)
(*vis) << *it_A;
vis->set_pen(0, 255, 0);
for (it_B.seek_begin(); it_B.not_done(); it_B++)
(*vis) << *it_B;
vis->set_pen(0, 0, 255);
for (it_C.seek_begin(); it_C.not_done(); it_C++)
(*vis) << *it_C;
}
#endif
}
// *** End of the main algorithm ***
// *** Clean up ***
delete part;
delete dir_domain;
delete Blocks_A;
delete Blocks_B;
delete Blocks_C;
#ifdef _VISUALIZE_2D_DECOMP_
if (domain.dimension() == 2)
delete vis;
#endif
// Return result
return result;
}
std::ostream& operator<<(std::ostream& os, const Classified_Block block);
Classified_Block::Classified_Block(int count)
: intersection_count(count)
{
}
Classified_Block::Classified_Block(const r_Minterval& b, int count)
: intersection_count(count),
block(b)
{
}
bool Classified_Block::operator==(const Classified_Block& other) const
{
return (this->block == other.block)
&& (this->intersection_count == other.intersection_count);
}
bool Classified_Block::operator!=(const Classified_Block& other) const
{
return !(*this == other);
}
std::ostream& operator<<(std::ostream& os, const Classified_Block block)
{
os << "CBlock(" << block.intersection_count << "x: " << block.block << ")";
return os;
}