Changeset a567c3


Ignore:
Timestamp:
Jul 28, 2009, 10:33:53 AM (16 years ago)
Author:
Frederik Heber <heber@…>
Branches:
Action_Thermostats, Add_AtomRandomPerturbation, Add_FitFragmentPartialChargesAction, Add_RotateAroundBondAction, Add_SelectAtomByNameAction, Added_ParseSaveFragmentResults, AddingActions_SaveParseParticleParameters, Adding_Graph_to_ChangeBondActions, Adding_MD_integration_tests, Adding_ParticleName_to_Atom, Adding_StructOpt_integration_tests, AtomFragments, Automaking_mpqc_open, AutomationFragmentation_failures, Candidate_v1.5.4, Candidate_v1.6.0, Candidate_v1.6.1, ChangeBugEmailaddress, ChangingTestPorts, ChemicalSpaceEvaluator, CombiningParticlePotentialParsing, Combining_Subpackages, Debian_Package_split, Debian_package_split_molecuildergui_only, Disabling_MemDebug, Docu_Python_wait, EmpiricalPotential_contain_HomologyGraph, EmpiricalPotential_contain_HomologyGraph_documentation, Enable_parallel_make_install, Enhance_userguide, Enhanced_StructuralOptimization, Enhanced_StructuralOptimization_continued, Example_ManyWaysToTranslateAtom, Exclude_Hydrogens_annealWithBondGraph, FitPartialCharges_GlobalError, Fix_BoundInBox_CenterInBox_MoleculeActions, Fix_ChargeSampling_PBC, Fix_ChronosMutex, Fix_FitPartialCharges, Fix_FitPotential_needs_atomicnumbers, Fix_ForceAnnealing, Fix_IndependentFragmentGrids, Fix_ParseParticles, Fix_ParseParticles_split_forward_backward_Actions, Fix_PopActions, Fix_QtFragmentList_sorted_selection, Fix_Restrictedkeyset_FragmentMolecule, Fix_StatusMsg, Fix_StepWorldTime_single_argument, Fix_Verbose_Codepatterns, Fix_fitting_potentials, Fixes, ForceAnnealing_goodresults, ForceAnnealing_oldresults, ForceAnnealing_tocheck, ForceAnnealing_with_BondGraph, ForceAnnealing_with_BondGraph_continued, ForceAnnealing_with_BondGraph_continued_betteresults, ForceAnnealing_with_BondGraph_contraction-expansion, FragmentAction_writes_AtomFragments, FragmentMolecule_checks_bonddegrees, GeometryObjects, Gui_Fixes, Gui_displays_atomic_force_velocity, ImplicitCharges, IndependentFragmentGrids, IndependentFragmentGrids_IndividualZeroInstances, IndependentFragmentGrids_IntegrationTest, IndependentFragmentGrids_Sole_NN_Calculation, JobMarket_RobustOnKillsSegFaults, JobMarket_StableWorkerPool, JobMarket_unresolvable_hostname_fix, MoreRobust_FragmentAutomation, ODR_violation_mpqc_open, PartialCharges_OrthogonalSummation, PdbParser_setsAtomName, PythonUI_with_named_parameters, QtGui_reactivate_TimeChanged_changes, Recreated_GuiChecks, Rewrite_FitPartialCharges, RotateToPrincipalAxisSystem_UndoRedo, SaturateAtoms_findBestMatching, SaturateAtoms_singleDegree, StoppableMakroAction, Subpackage_CodePatterns, Subpackage_JobMarket, Subpackage_LinearAlgebra, Subpackage_levmar, Subpackage_mpqc_open, Subpackage_vmg, Switchable_LogView, ThirdParty_MPQC_rebuilt_buildsystem, TrajectoryDependenant_MaxOrder, TremoloParser_IncreasedPrecision, TremoloParser_MultipleTimesteps, TremoloParser_setsAtomName, Ubuntu_1604_changes, stable
Children:
d4d0dd
Parents:
fcc9b0
git-author:
Frederik Heber <heber@…> (07/28/09 10:31:26)
git-committer:
Frederik Heber <heber@…> (07/28/09 10:33:53)
Message:

ConvexTesselation working again.

  • File was not written, as NULL was given in builder.cpp instead of argv[argptr]
  • Tesselation::TesselateOnBoundary() - lots of small errors with the normal vector of the triangle, the propagation direction check and so forth
  • DetermineCenterOfAll() has sign changed, this was adapted in boundary.cpp
  • heptan is NOT working yet, as too many boundary points are thrown away
Location:
src
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • src/boundary.cpp

    rfcc9b0 ra567c3  
    682682void write_tecplot_file(ofstream *out, ofstream *tecplot, class Tesselation *TesselStruct, class molecule *mol, int N)
    683683{
    684   if (tecplot != NULL)
    685     {
    686       *tecplot << "TITLE = \"3D CONVEX SHELL\"" << endl;
    687       *tecplot << "VARIABLES = \"X\" \"Y\" \"Z\"" << endl;
    688       *tecplot << "ZONE T=\"TRIANGLES" << N << "\", N="
    689           << TesselStruct->PointsOnBoundaryCount << ", E="
    690           << TesselStruct->TrianglesOnBoundaryCount
    691           << ", DATAPACKING=POINT, ZONETYPE=FETRIANGLE" << endl;
    692       int *LookupList = new int[mol->AtomCount];
    693       for (int i = 0; i < mol->AtomCount; i++)
    694         LookupList[i] = -1;
    695 
    696       // print atom coordinates
    697       *out << Verbose(2) << "The following triangles were created:";
    698       int Counter = 1;
    699       atom *Walker = NULL;
    700       for (PointMap::iterator target = TesselStruct->PointsOnBoundary.begin(); target
    701           != TesselStruct->PointsOnBoundary.end(); target++)
    702         {
    703           Walker = target->second->node;
    704           LookupList[Walker->nr] = Counter++;
    705           *tecplot << Walker->x.x[0] << " " << Walker->x.x[1] << " "
    706               << Walker->x.x[2] << " " << endl;
    707         }
    708       *tecplot << endl;
    709       // print connectivity
    710       for (TriangleMap::iterator runner =
    711           TesselStruct->TrianglesOnBoundary.begin(); runner
    712           != TesselStruct->TrianglesOnBoundary.end(); runner++)
    713         {
    714           *out << " " << runner->second->endpoints[0]->node->Name << "<->"
    715               << runner->second->endpoints[1]->node->Name << "<->"
    716               << runner->second->endpoints[2]->node->Name;
    717           *tecplot << LookupList[runner->second->endpoints[0]->node->nr] << " "
    718               << LookupList[runner->second->endpoints[1]->node->nr] << " "
    719               << LookupList[runner->second->endpoints[2]->node->nr] << endl;
    720         }
    721       delete[] (LookupList);
    722       *out << endl;
    723     }
     684  if (tecplot != NULL) {
     685    *tecplot << "TITLE = \"3D CONVEX SHELL\"" << endl;
     686    *tecplot << "VARIABLES = \"X\" \"Y\" \"Z\"" << endl;
     687    *tecplot << "ZONE T=\"TRIANGLES" << N << "\", N=" << TesselStruct->PointsOnBoundaryCount << ", E=" << TesselStruct->TrianglesOnBoundaryCount << ", DATAPACKING=POINT, ZONETYPE=FETRIANGLE" << endl;
     688    int *LookupList = new int[mol->AtomCount];
     689    for (int i = 0; i < mol->AtomCount; i++)
     690      LookupList[i] = -1;
     691
     692    // print atom coordinates
     693    *out << Verbose(2) << "The following triangles were created:";
     694    int Counter = 1;
     695    atom *Walker = NULL;
     696    for (PointMap::iterator target = TesselStruct->PointsOnBoundary.begin(); target != TesselStruct->PointsOnBoundary.end(); target++) {
     697      Walker = target->second->node;
     698      LookupList[Walker->nr] = Counter++;
     699      *tecplot << Walker->x.x[0] << " " << Walker->x.x[1] << " " << Walker->x.x[2] << " " << endl;
     700    }
     701    *tecplot << endl;
     702    // print connectivity
     703    for (TriangleMap::iterator runner = TesselStruct->TrianglesOnBoundary.begin(); runner != TesselStruct->TrianglesOnBoundary.end(); runner++) {
     704      *out << " " << runner->second->endpoints[0]->node->Name << "<->" << runner->second->endpoints[1]->node->Name << "<->" << runner->second->endpoints[2]->node->Name;
     705      *tecplot << LookupList[runner->second->endpoints[0]->node->nr] << " " << LookupList[runner->second->endpoints[1]->node->nr] << " " << LookupList[runner->second->endpoints[2]->node->nr] << endl;
     706    }
     707    delete[] (LookupList);
     708    *out << endl;
     709  }
    724710}
    725711
     
    791777  // 4. Store triangles in tecplot file
    792778  if (filename != NULL) {
    793     string OutputName(filename);
    794     OutputName.append(TecplotSuffix);
    795     ofstream *tecplot = new ofstream(OutputName.c_str());
    796     write_tecplot_file(out, tecplot, TesselStruct, mol, 0);
    797     tecplot->close();
    798     delete(tecplot);
    799     ofstream *rasterplot = new ofstream(OutputName.c_str());
    800     write_raster3d_file(out, rasterplot, TesselStruct, mol);
    801     rasterplot->close();
    802     delete(rasterplot);
     779    if (DoTecplotOutput) {
     780      string OutputName(filename);
     781      OutputName.append(TecplotSuffix);
     782      ofstream *tecplot = new ofstream(OutputName.c_str());
     783      write_tecplot_file(out, tecplot, TesselStruct, mol, 0);
     784      tecplot->close();
     785      delete(tecplot);
     786    }
     787    if (DoRaster3DOutput) {
     788      string OutputName(filename);
     789      OutputName.append(Raster3DSuffix);
     790      ofstream *rasterplot = new ofstream(OutputName.c_str());
     791      write_raster3d_file(out, rasterplot, TesselStruct, mol);
     792      rasterplot->close();
     793      delete(rasterplot);
     794    }
    803795  }
    804796
     
    11001092                  << A->second->node->Name << ","
    11011093                  << baseline->second.first->second->node->Name << ","
    1102                   << baseline->second.second->second->node->Name << " leave "
     1094                  << baseline->second.second->second->node->Name << " leaves "
    11031095                  << checker->second->node->Name << " outside the convex hull."
    11041096                  << endl;
     
    11851177 * -# if the lines contains to only one triangle
    11861178 * -# We search all points in the boundary
    1187  *    -# if the triangle with the baseline and the current point has the smallest of angles (comparison between normal vectors
    11881179 *    -# if the triangle is in forward direction of the baseline (at most 90 degrees angle between vector orthogonal to
    11891180 *       baseline in triangle plane pointing out of the triangle and normal vector of new triangle)
     1181 *    -# if the triangle with the baseline and the current point has the smallest of angles (comparison between normal vectors)
    11901182 *    -# then we have a new triangle, whose baselines we again add (or increase their TriangleCount)
    11911183 * \param *out output stream for debugging
     
    12071199    for (LineMap::iterator baseline = LinesOnBoundary.begin(); baseline != LinesOnBoundary.end(); baseline++)
    12081200      if (baseline->second->TrianglesCount == 1) {
    1209         *out << Verbose(2) << "Current baseline is between " << *(baseline->second) << "." << endl;
    12101201        // 5a. go through each boundary point if not _both_ edges between either endpoint of the current line and this point exist (and belong to 2 triangles)
    12111202        SmallestAngle = M_PI;
     1203
     1204        // get peak point with respect to this base line's only triangle
    12121205        BTS = baseline->second->triangles.begin()->second; // there is only one triangle so far
    1213         // get peak point with respect to this base line's only triangle
     1206        *out << Verbose(2) << "Current baseline is between " << *(baseline->second) << "." << endl;
    12141207        for (int i = 0; i < 3; i++)
    12151208          if ((BTS->endpoints[i] != baseline->second->endpoints[0]) && (BTS->endpoints[i] != baseline->second->endpoints[1]))
    12161209            peak = BTS->endpoints[i];
    12171210        *out << Verbose(3) << " and has peak " << *peak << "." << endl;
     1211
     1212        // prepare some auxiliary vectors
     1213        Vector BaseLineCenter, BaseLine;
     1214        BaseLineCenter.CopyVector(&baseline->second->endpoints[0]->node->x);
     1215        BaseLineCenter.AddVector(&baseline->second->endpoints[0]->node->x);
     1216        BaseLineCenter.Scale(1. / 2.); // points now to center of base line
     1217        BaseLine.CopyVector(&baseline->second->endpoints[0]->node->x);
     1218        BaseLine.SubtractVector(&baseline->second->endpoints[1]->node->x);
     1219
    12181220        // offset to center of triangle
    12191221        CenterVector.Zero();
     
    12221224        CenterVector.Scale(1. / 3.);
    12231225        *out << Verbose(4) << "CenterVector of base triangle is " << CenterVector << endl;
     1226
     1227        // normal vector of triangle
    12241228        NormalVector.CopyVector(MolCenter);
    12251229        NormalVector.SubtractVector(&CenterVector);
    1226         // normal vector of triangle
    12271230        BTS->GetNormalVector(NormalVector);
     1231        NormalVector.CopyVector(&BTS->NormalVector);
    12281232        *out << Verbose(4) << "NormalVector of base triangle is " << NormalVector << endl;
     1233
    12291234        // vector in propagation direction (out of triangle)
    12301235        // project center vector onto triangle plane (points from intersection plane-NormalVector to plane-CenterVector intersection)
    1231         TempVector.CopyVector(&baseline->second->endpoints[0]->node->x);
    1232         TempVector.SubtractVector(&baseline->second->endpoints[1]->node->x);
    1233         PropagationVector.MakeNormalVector(&TempVector, &NormalVector);
     1236        PropagationVector.MakeNormalVector(&BaseLine, &NormalVector);
    12341237        TempVector.CopyVector(&CenterVector);
    12351238        TempVector.SubtractVector(&baseline->second->endpoints[0]->node->x); // TempVector is vector on triangle plane pointing from one baseline egde towards center!
     
    12371240        if (PropagationVector.Projection(&TempVector) > 0) // make sure normal propagation vector points outward from baseline
    12381241          PropagationVector.Scale(-1.);
    1239         *out << Verbose(4) << "PropagationVector of base triangle is ";
    1240         PropagationVector.Output(out);
    1241         *out << endl;
     1242        *out << Verbose(4) << "PropagationVector of base triangle is " << PropagationVector << endl;
    12421243        winner = PointsOnBoundary.end();
    1243         for (PointMap::iterator target = PointsOnBoundary.begin(); target
    1244             != PointsOnBoundary.end(); target++)
    1245           if ((target->second != baseline->second->endpoints[0])
    1246               && (target->second != baseline->second->endpoints[1]))
    1247             { // don't take the same endpoints
    1248               *out << Verbose(3) << "Target point is " << *(target->second)
    1249                   << ":";
    1250               bool continueflag = true;
    1251 
    1252               VirtualNormalVector.CopyVector(
    1253                   &baseline->second->endpoints[0]->node->x);
    1254               VirtualNormalVector.AddVector(
    1255                   &baseline->second->endpoints[0]->node->x);
    1256               VirtualNormalVector.Scale(-1. / 2.); // points now to center of base line
    1257               VirtualNormalVector.AddVector(&target->second->node->x); // points from center of base line to target
    1258               TempAngle = VirtualNormalVector.Angle(&PropagationVector);
    1259               continueflag = continueflag && (TempAngle < (M_PI/2.)); // no bends bigger than Pi/2 (90 degrees)
    1260               if (!continueflag)
    1261                 {
    1262                   *out << Verbose(4)
    1263                       << "Angle between propagation direction and base line to "
    1264                       << *(target->second) << " is " << TempAngle
    1265                       << ", bad direction!" << endl;
    1266                   continue;
    1267                 }
    1268               else
    1269                 *out << Verbose(4)
    1270                     << "Angle between propagation direction and base line to "
    1271                     << *(target->second) << " is " << TempAngle
    1272                     << ", good direction!" << endl;
    1273               LineChecker[0] = baseline->second->endpoints[0]->lines.find(
    1274                   target->first);
    1275               LineChecker[1] = baseline->second->endpoints[1]->lines.find(
    1276                   target->first);
    1277               //            if (LineChecker[0] != baseline->second->endpoints[0]->lines.end())
    1278               //              *out << Verbose(4) << *(baseline->second->endpoints[0]) << " has line " << *(LineChecker[0]->second) << " to " << *(target->second) << " as endpoint with " << LineChecker[0]->second->TrianglesCount << " triangles." << endl;
    1279               //            else
    1280               //              *out << Verbose(4) << *(baseline->second->endpoints[0]) << " has no line to " << *(target->second) << " as endpoint." << endl;
    1281               //            if (LineChecker[1] != baseline->second->endpoints[1]->lines.end())
    1282               //              *out << Verbose(4) << *(baseline->second->endpoints[1]) << " has line " << *(LineChecker[1]->second) << " to " << *(target->second) << " as endpoint with " << LineChecker[1]->second->TrianglesCount << " triangles." << endl;
    1283               //            else
    1284               //              *out << Verbose(4) << *(baseline->second->endpoints[1]) << " has no line to " << *(target->second) << " as endpoint." << endl;
    1285               // check first endpoint (if any connecting line goes to target or at least not more than 1)
    1286               continueflag = continueflag && (((LineChecker[0]
    1287                   == baseline->second->endpoints[0]->lines.end())
    1288                   || (LineChecker[0]->second->TrianglesCount == 1)));
    1289               if (!continueflag)
    1290                 {
    1291                   *out << Verbose(4) << *(baseline->second->endpoints[0])
    1292                       << " has line " << *(LineChecker[0]->second)
    1293                       << " to " << *(target->second)
    1294                       << " as endpoint with "
    1295                       << LineChecker[0]->second->TrianglesCount
    1296                       << " triangles." << endl;
    1297                   continue;
    1298                 }
    1299               // check second endpoint (if any connecting line goes to target or at least not more than 1)
    1300               continueflag = continueflag && (((LineChecker[1]
    1301                   == baseline->second->endpoints[1]->lines.end())
    1302                   || (LineChecker[1]->second->TrianglesCount == 1)));
    1303               if (!continueflag)
    1304                 {
    1305                   *out << Verbose(4) << *(baseline->second->endpoints[1])
    1306                       << " has line " << *(LineChecker[1]->second)
    1307                       << " to " << *(target->second)
    1308                       << " as endpoint with "
    1309                       << LineChecker[1]->second->TrianglesCount
    1310                       << " triangles." << endl;
    1311                   continue;
    1312                 }
    1313               // check whether the envisaged triangle does not already exist (if both lines exist and have same endpoint)
    1314               continueflag = continueflag && (!(((LineChecker[0]
    1315                   != baseline->second->endpoints[0]->lines.end())
    1316                   && (LineChecker[1]
    1317                       != baseline->second->endpoints[1]->lines.end())
    1318                   && (GetCommonEndpoint(LineChecker[0]->second,
    1319                       LineChecker[1]->second) == peak))));
    1320               if (!continueflag)
    1321                 {
    1322                   *out << Verbose(4) << "Current target is peak!" << endl;
    1323                   continue;
    1324                 }
    1325               // in case NOT both were found
    1326               if (continueflag)
    1327                 { // create virtually this triangle, get its normal vector, calculate angle
    1328                   flag = true;
    1329                   VirtualNormalVector.MakeNormalVector(
    1330                       &baseline->second->endpoints[0]->node->x,
    1331                       &baseline->second->endpoints[1]->node->x,
    1332                       &target->second->node->x);
    1333                   // make it always point inward
    1334                   if (baseline->second->endpoints[0]->node->x.Projection(
    1335                       &VirtualNormalVector) > 0)
    1336                     VirtualNormalVector.Scale(-1.);
    1337                   // calculate angle
    1338                   TempAngle = NormalVector.Angle(&VirtualNormalVector);
    1339                   *out << Verbose(4) << "NormalVector is ";
    1340                   VirtualNormalVector.Output(out);
    1341                   *out << " and the angle is " << TempAngle << "." << endl;
    1342                   if (SmallestAngle > TempAngle)
    1343                     { // set to new possible winner
    1344                       SmallestAngle = TempAngle;
    1345                       winner = target;
    1346                     }
    1347                 }
     1244
     1245        // loop over all points and calculate angle between normal vector of new and present triangle
     1246        for (PointMap::iterator target = PointsOnBoundary.begin(); target != PointsOnBoundary.end(); target++) {
     1247          if ((target->second != baseline->second->endpoints[0]) && (target->second != baseline->second->endpoints[1])) { // don't take the same endpoints
     1248            *out << Verbose(3) << "Target point is " << *(target->second) << ":" << endl;
     1249
     1250            // first check direction, so that triangles don't intersect
     1251            VirtualNormalVector.CopyVector(&target->second->node->x);
     1252            VirtualNormalVector.SubtractVector(&BaseLineCenter); // points from center of base line to target
     1253            VirtualNormalVector.ProjectOntoPlane(&NormalVector);
     1254            TempAngle = VirtualNormalVector.Angle(&PropagationVector);
     1255            if (TempAngle > (M_PI/2.)) { // no bends bigger than Pi/2 (90 degrees)
     1256              *out << Verbose(4) << "Angle on triangle plane between propagation direction and base line to " << *(target->second) << " is " << TempAngle << ", bad direction!" << endl;
     1257              continue;
     1258            } else
     1259              *out << Verbose(4) << "Angle on triangle plane between propagation direction and base line to " << *(target->second) << " is " << TempAngle << ", good direction!" << endl;
     1260
     1261            // check first and second endpoint (if any connecting line goes to target has at least not more than 1 triangle)
     1262            LineChecker[0] = baseline->second->endpoints[0]->lines.find(target->first);
     1263            LineChecker[1] = baseline->second->endpoints[1]->lines.find(target->first);
     1264            if (((LineChecker[0] != baseline->second->endpoints[0]->lines.end()) && (LineChecker[0]->second->TrianglesCount == 2))) {
     1265              *out << Verbose(4) << *(baseline->second->endpoints[0]) << " has line " << *(LineChecker[0]->second) << " to " << *(target->second) << " as endpoint with " << LineChecker[0]->second->TrianglesCount << " triangles." << endl;
     1266              continue;
    13481267            }
     1268            if (((LineChecker[1] != baseline->second->endpoints[1]->lines.end()) && (LineChecker[1]->second->TrianglesCount == 2))) {
     1269              *out << Verbose(4) << *(baseline->second->endpoints[1]) << " has line " << *(LineChecker[1]->second) << " to " << *(target->second) << " as endpoint with " << LineChecker[1]->second->TrianglesCount << " triangles." << endl;
     1270              continue;
     1271            }
     1272
     1273            // check whether the envisaged triangle does not already exist (if both lines exist and have same endpoint)
     1274            if ((((LineChecker[0] != baseline->second->endpoints[0]->lines.end()) && (LineChecker[1] != baseline->second->endpoints[1]->lines.end()) && (GetCommonEndpoint(LineChecker[0]->second, LineChecker[1]->second) == peak)))) {
     1275              *out << Verbose(4) << "Current target is peak!" << endl;
     1276              continue;
     1277            }
     1278
     1279            // in case NOT both were found, create virtually this triangle, get its normal vector, calculate angle
     1280            flag = true;
     1281            VirtualNormalVector.MakeNormalVector(&baseline->second->endpoints[0]->node->x, &baseline->second->endpoints[1]->node->x, &target->second->node->x);
     1282            TempVector.CopyVector(&baseline->second->endpoints[0]->node->x);
     1283            TempVector.AddVector(&baseline->second->endpoints[1]->node->x);
     1284            TempVector.AddVector(&target->second->node->x);
     1285            TempVector.Scale(1./3.);
     1286            TempVector.SubtractVector(MolCenter);
     1287            // make it always point outward
     1288            if (VirtualNormalVector.Projection(&TempVector) < 0)
     1289              VirtualNormalVector.Scale(-1.);
     1290            // calculate angle
     1291            TempAngle = NormalVector.Angle(&VirtualNormalVector);
     1292            *out << Verbose(4) << "NormalVector is " << VirtualNormalVector << " and the angle is " << TempAngle << "." << endl;
     1293            if (SmallestAngle > TempAngle) { // set to new possible winner
     1294              SmallestAngle = TempAngle;
     1295              winner = target;
     1296            }
     1297          }
     1298        } // end of loop over all boundary points
     1299
    13491300        // 5b. The point of the above whose triangle has the greatest angle with the triangle the current line belongs to (it only belongs to one, remember!): New triangle
    1350         if (winner != PointsOnBoundary.end())
    1351           {
    1352             *out << Verbose(2) << "Winning target point is "
    1353                 << *(winner->second) << " with angle " << SmallestAngle
    1354                 << "." << endl;
    1355             // create the lins of not yet present
    1356             BLS[0] = baseline->second;
    1357             // 5c. add lines to the line set if those were new (not yet part of a triangle), delete lines that belong to two triangles)
    1358             LineChecker[0] = baseline->second->endpoints[0]->lines.find(
    1359                 winner->first);
    1360             LineChecker[1] = baseline->second->endpoints[1]->lines.find(
    1361                 winner->first);
    1362             if (LineChecker[0]
    1363                 == baseline->second->endpoints[0]->lines.end())
    1364               { // create
    1365                 BPS[0] = baseline->second->endpoints[0];
    1366                 BPS[1] = winner->second;
    1367                 BLS[1] = new class BoundaryLineSet(BPS,
    1368                     LinesOnBoundaryCount);
    1369                 LinesOnBoundary.insert(LinePair(LinesOnBoundaryCount,
    1370                     BLS[1]));
    1371                 LinesOnBoundaryCount++;
    1372               }
    1373             else
    1374               BLS[1] = LineChecker[0]->second;
    1375             if (LineChecker[1]
    1376                 == baseline->second->endpoints[1]->lines.end())
    1377               { // create
    1378                 BPS[0] = baseline->second->endpoints[1];
    1379                 BPS[1] = winner->second;
    1380                 BLS[2] = new class BoundaryLineSet(BPS,
    1381                     LinesOnBoundaryCount);
    1382                 LinesOnBoundary.insert(LinePair(LinesOnBoundaryCount,
    1383                     BLS[2]));
    1384                 LinesOnBoundaryCount++;
    1385               }
    1386             else
    1387               BLS[2] = LineChecker[1]->second;
    1388             BTS = new class BoundaryTriangleSet(BLS,
    1389                 TrianglesOnBoundaryCount);
    1390             TrianglesOnBoundary.insert(TrianglePair(
    1391                 TrianglesOnBoundaryCount, BTS));
    1392             TrianglesOnBoundaryCount++;
    1393           }
    1394         else
    1395           {
    1396             *out << Verbose(1)
    1397                 << "I could not determine a winner for this baseline "
    1398                 << *(baseline->second) << "." << endl;
    1399           }
     1301        if (winner != PointsOnBoundary.end()) {
     1302          *out << Verbose(2) << "Winning target point is " << *(winner->second) << " with angle " << SmallestAngle << "." << endl;
     1303          // create the lins of not yet present
     1304          BLS[0] = baseline->second;
     1305          // 5c. add lines to the line set if those were new (not yet part of a triangle), delete lines that belong to two triangles)
     1306          LineChecker[0] = baseline->second->endpoints[0]->lines.find(winner->first);
     1307          LineChecker[1] = baseline->second->endpoints[1]->lines.find(winner->first);
     1308          if (LineChecker[0] == baseline->second->endpoints[0]->lines.end()) { // create
     1309            BPS[0] = baseline->second->endpoints[0];
     1310            BPS[1] = winner->second;
     1311            BLS[1] = new class BoundaryLineSet(BPS, LinesOnBoundaryCount);
     1312            LinesOnBoundary.insert(LinePair(LinesOnBoundaryCount, BLS[1]));
     1313            LinesOnBoundaryCount++;
     1314          } else
     1315            BLS[1] = LineChecker[0]->second;
     1316          if (LineChecker[1] == baseline->second->endpoints[1]->lines.end()) { // create
     1317            BPS[0] = baseline->second->endpoints[1];
     1318            BPS[1] = winner->second;
     1319            BLS[2] = new class BoundaryLineSet(BPS, LinesOnBoundaryCount);
     1320            LinesOnBoundary.insert(LinePair(LinesOnBoundaryCount, BLS[2]));
     1321            LinesOnBoundaryCount++;
     1322          } else
     1323            BLS[2] = LineChecker[1]->second;
     1324          BTS = new class BoundaryTriangleSet(BLS, TrianglesOnBoundaryCount);
     1325          TrianglesOnBoundary.insert(TrianglePair(TrianglesOnBoundaryCount, BTS));
     1326          TrianglesOnBoundaryCount++;
     1327        } else {
     1328          *out << Verbose(1) << "I could not determine a winner for this baseline " << *(baseline->second) << "." << endl;
     1329        }
    14001330
    14011331        // 5d. If the set of lines is not yet empty, go to 5. and continue
    1402       }
    1403       else
    1404         *out << Verbose(2) << "Baseline candidate " << *(baseline->second)
    1405             << " has a triangle count of "
    1406             << baseline->second->TrianglesCount << "." << endl;
     1332      } else
     1333        *out << Verbose(2) << "Baseline candidate " << *(baseline->second) << " has a triangle count of " << baseline->second->TrianglesCount << "." << endl;
    14071334  } while (flag);
    14081335  delete(MolCenter);
     
    27952722    }
    27962723  }
    2797   if (1) { //failflag) {
     2724  if (filename != 0) {
    27982725    *out << Verbose(1) << "Writing final tecplot file\n";
    27992726    if (DoTecplotOutput) {
  • src/builder.cpp

    rfcc9b0 ra567c3  
    18131813                LinkedCell LCList(mol, 10.);
    18141814                class Tesselation *TesselStruct = NULL;
    1815                 Find_convex_border((ofstream *)&cout, mol, TesselStruct, &LCList, NULL);
     1815                Find_convex_border((ofstream *)&cout, mol, TesselStruct, &LCList, argv[argptr]);
    18161816                double clustervolume = VolumeOfConvexEnvelope((ofstream *)&cout, TesselStruct, &configuration);
    18171817                delete(TesselStruct);
Note: See TracChangeset for help on using the changeset viewer.