/*
* Project: MoleCuilder
* Description: creates and alters molecular systems
* Copyright (C) 2010-2012 University of Bonn. All rights reserved.
*
*
* This file is part of MoleCuilder.
*
* MoleCuilder 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 2 of the License, or
* (at your option) any later version.
*
* MoleCuilder 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 MoleCuilder. If not, see .
*/
/*
* BoundaryLineSet.cpp
*
* Created on: Jul 29, 2010
* Author: heber
*/
// include config.h
#ifdef HAVE_CONFIG_H
#include
#endif
#include "CodePatterns/MemDebug.hpp"
#include "BoundaryLineSet.hpp"
#include
#include "BoundaryPointSet.hpp"
#include "BoundaryTriangleSet.hpp"
#include "Atom/TesselPoint.hpp"
#include "CodePatterns/Assert.hpp"
#include "CodePatterns/Info.hpp"
#include "CodePatterns/Log.hpp"
#include "CodePatterns/Verbose.hpp"
#include "tesselationhelpers.hpp"
#include "LinearAlgebra/Vector.hpp"
using namespace std;
/** Constructor of BoundaryLineSet.
*/
BoundaryLineSet::BoundaryLineSet() :
Nr(-1)
{
//Info FunctionInfo(__func__);
for (int i = 0; i < 2; i++)
endpoints[i] = NULL;
}
;
/** Constructor of BoundaryLineSet with two endpoints.
* Adds line automatically to each endpoints' LineMap
* \param *Point[2] array of two boundary points
* \param number number of the list
*/
BoundaryLineSet::BoundaryLineSet(BoundaryPointSet * const Point[2], const int number)
{
//Info FunctionInfo(__func__);
// set number
Nr = number;
// set endpoints in ascending order
SetEndpointsOrdered(endpoints, Point[0], Point[1]);
// add this line to the hash maps of both endpoints
Point[0]->AddLine(this); //Taken out, to check whether we can avoid unwanted double adding.
Point[1]->AddLine(this); //
// set skipped to false
skipped = false;
// clear triangles list
LOG(5, "DEBUG: New Line with endpoints " << *this << ".");
}
;
/** Constructor of BoundaryLineSet with two endpoints.
* Adds line automatically to each endpoints' LineMap
* \param *Point1 first boundary point
* \param *Point2 second boundary point
* \param number number of the list
*/
BoundaryLineSet::BoundaryLineSet(BoundaryPointSet * const Point1, BoundaryPointSet * const Point2, const int number) :
Nr(number),
skipped(false)
{
//Info FunctionInfo(__func__);
// set endpoints in ascending order
SetEndpointsOrdered(endpoints, Point1, Point2);
// add this line to the hash maps of both endpoints
Point1->AddLine(this); //Taken out, to check whether we can avoid unwanted double adding.
Point2->AddLine(this); //
// clear triangles list
LOG(5, "DEBUG: New Line with endpoints " << *this << ".");
}
;
/** Destructor for BoundaryLineSet.
* Removes itself from each endpoints' LineMap, calling RemoveTrianglePoint() when point not connected anymore.
* \note When removing lines from a class Tesselation, use RemoveTesselationLine()
*/
BoundaryLineSet::~BoundaryLineSet()
{
//Info FunctionInfo(__func__);
int Numbers[2];
// get other endpoint number of finding copies of same line
if (endpoints[1] != NULL)
Numbers[0] = endpoints[1]->Nr;
else
Numbers[0] = -1;
if (endpoints[0] != NULL)
Numbers[1] = endpoints[0]->Nr;
else
Numbers[1] = -1;
for (int i = 0; i < 2; i++) {
if (endpoints[i] != NULL) {
if (Numbers[i] != -1) { // as there may be multiple lines with same endpoints, we have to go through each and find in the endpoint's line list this line set
pair erasor = endpoints[i]->lines.equal_range(Numbers[i]);
for (LineMap::iterator Runner = erasor.first; Runner != erasor.second; Runner++)
if ((*Runner).second == this) {
//LOG(0, "Removing Line Nr. " << Nr << " in boundary point " << *endpoints[i] << ".");
endpoints[i]->lines.erase(Runner);
break;
}
} else { // there's just a single line left
if (endpoints[i]->lines.erase(Nr)) {
//LOG(0, "Removing Line Nr. " << Nr << " in boundary point " << *endpoints[i] << ".");
}
}
if (endpoints[i]->lines.empty()) {
//LOG(0, *endpoints[i] << " has no more lines it's attached to, erasing.");
if (endpoints[i] != NULL) {
delete (endpoints[i]);
endpoints[i] = NULL;
}
}
}
}
if (!triangles.empty())
ELOG(2, "Memory Leak! I " << *this << " am still connected to some triangles.");
}
;
/** Add triangle to TriangleMap of this boundary line.
* \param *triangle to add
*/
void BoundaryLineSet::AddTriangle(BoundaryTriangleSet * const triangle)
{
//Info FunctionInfo(__func__);
LOG(5, "DEBUG: Add " << triangle->Nr << " to line " << *this << ".");
triangles.insert(TrianglePair(triangle->Nr, triangle));
}
;
/** Checks whether we have a common endpoint with given \a *line.
* \param *line other line to test
* \return true - common endpoint present, false - not connected
*/
bool BoundaryLineSet::IsConnectedTo(const BoundaryLineSet * const line) const
{
//Info FunctionInfo(__func__);
if ((endpoints[0] == line->endpoints[0]) || (endpoints[1] == line->endpoints[0]) || (endpoints[0] == line->endpoints[1]) || (endpoints[1] == line->endpoints[1]))
return true;
else
return false;
}
;
/** Checks whether the adjacent triangles of a baseline are convex or not.
* We sum the two angles of each height vector with respect to the center of the baseline.
* If greater/equal M_PI than we are convex.
* \param *out output stream for debugging
* \return true - triangles are convex, false - concave or less than two triangles connected
*/
bool BoundaryLineSet::CheckConvexityCriterion() const
{
//Info FunctionInfo(__func__);
double angle = CalculateConvexity();
if (angle > -MYEPSILON) {
LOG(3, "ACCEPT: Angle is greater than pi: convex.");
return true;
} else {
LOG(3, "REJECT: Angle is less than pi: concave.");
return false;
}
}
/** Calculates the angle between two triangles with respect to their normal vector.
* We sum the two angles of each height vector with respect to the center of the baseline.
* \return angle > 0 then convex, if < 0 then concave
*/
double BoundaryLineSet::CalculateConvexity() const
{
//Info FunctionInfo(__func__);
Vector BaseLineCenter, BaseLineNormal, BaseLine, helper[2], NormalCheck;
// get the two triangles
if (triangles.size() != 2) {
ELOG(0, "Baseline " << *this << " is connected to less than two triangles, Tesselation incomplete!");
return true;
}
// check normal vectors
// have a normal vector on the base line pointing outwards
//LOG(0, "INFO: " << *this << " has vectors at " << *(endpoints[0]->node->node) << " and at " << *(endpoints[1]->node->node) << ".");
BaseLineCenter = (1./2.)*((endpoints[0]->node->getPosition()) + (endpoints[1]->node->getPosition()));
BaseLine = (endpoints[0]->node->getPosition()) - (endpoints[1]->node->getPosition());
//LOG(0, "INFO: Baseline is " << BaseLine << " and its center is at " << BaseLineCenter << ".");
BaseLineNormal.Zero();
NormalCheck.Zero();
double sign = -1.;
int i = 0;
class BoundaryPointSet *node = NULL;
for (TriangleMap::const_iterator runner = triangles.begin(); runner != triangles.end(); runner++) {
//LOG(0, "INFO: NormalVector of " << *(runner->second) << " is " << runner->second->NormalVector << ".");
NormalCheck += runner->second->NormalVector;
NormalCheck *= sign;
sign = -sign;
if (runner->second->NormalVector.NormSquared() > MYEPSILON)
BaseLineNormal = runner->second->NormalVector; // yes, copy second on top of first
else {
ELOG(0, "Triangle " << *runner->second << " has zero normal vector!");
}
node = runner->second->GetThirdEndpoint(this);
if (node != NULL) {
//LOG(0, "INFO: Third node for triangle " << *(runner->second) << " is " << *node << " at " << *(node->node->node) << ".");
helper[i] = (node->node->getPosition()) - BaseLineCenter;
helper[i].MakeNormalTo(BaseLine); // we want to compare the triangle's heights' angles!
//LOG(0, "INFO: Height vector with respect to baseline is " << helper[i] << ".");
i++;
} else {
ELOG(1, "I cannot find third node in triangle, something's wrong.");
return true;
}
}
//LOG(0, "INFO: BaselineNormal is " << BaseLineNormal << ".");
if (NormalCheck.NormSquared() < MYEPSILON) {
LOG(3, "ACCEPT: Normalvectors of both triangles are the same: convex.");
return true;
}
BaseLineNormal.Scale(-1.);
double angle = GetAngle(helper[0], helper[1], BaseLineNormal);
return (angle - M_PI);
}
/** Checks whether point is any of the two endpoints this line contains.
* \param *point point to test
* \return true - point is of the line, false - is not
*/
bool BoundaryLineSet::ContainsBoundaryPoint(const BoundaryPointSet * const point) const
{
//Info FunctionInfo(__func__);
for (int i = 0; i < 2; i++)
if (point == endpoints[i])
return true;
return false;
}
;
/** Returns other endpoint of the line.
* \param *point other endpoint
* \return NULL - if endpoint not contained in BoundaryLineSet::lines, or pointer to BoundaryPointSet otherwise
*/
class BoundaryPointSet *BoundaryLineSet::GetOtherEndpoint(const BoundaryPointSet * const point) const
{
//Info FunctionInfo(__func__);
if (endpoints[0] == point)
return endpoints[1];
else if (endpoints[1] == point)
return endpoints[0];
else
return NULL;
}
;
/** Returns other triangle of the line.
* \param *point other endpoint
* \return NULL - if triangle not contained in BoundaryLineSet::triangles, or pointer to BoundaryTriangleSet otherwise
*/
class BoundaryTriangleSet *BoundaryLineSet::GetOtherTriangle(const BoundaryTriangleSet * const triangle) const
{
//Info FunctionInfo(__func__);
if (triangles.size() == 2) {
for (TriangleMap::const_iterator TriangleRunner = triangles.begin(); TriangleRunner != triangles.end(); ++TriangleRunner)
if (TriangleRunner->second != triangle)
return TriangleRunner->second;
}
return NULL;
}
;
/** output operator for BoundaryLineSet.
* \param &ost output stream
* \param &a boundary line
*/
ostream & operator <<(ostream &ost, const BoundaryLineSet &a)
{
ost << "[" << a.Nr << "|" << a.endpoints[0]->node->getName() << "," << a.endpoints[1]->node->getName() << "]";
//ost << "[" << a.Nr << "|" << a.endpoints[0]->node->getName() << " at " << a.endpoints[0]->node->getPosition() << "," << a.endpoints[1]->node->getName() << " at " << a.endpoints[1]->node->getPosition() << "]";
return ost;
}
;