use Point/Trig/Chart-Id in stl-meshing, more to come

This commit is contained in:
Joachim Schöberl 2019-09-20 22:21:39 +02:00
parent 8381ce58ba
commit 61c9e669c3
12 changed files with 270 additions and 227 deletions

View File

@ -759,6 +759,7 @@ namespace netgen
{
for (int i = 1; i <= chartboundpoints.Size(); i++)
{
pindex.Append(-1);
plainpoints.Append (chartboundpoints.Get(i));
locpoints.Append (chartboundpoints3d.Get(i));
legalpoints.Append (0);
@ -1986,6 +1987,20 @@ namespace netgen
namespace netgen
{
void glrender (int wait)
{ ; }
{ ;
/*
if (multithread.drawing)
{
// vssurfacemeshing.Render();
// Render ();
if (wait || multithread.testmode)
{
multithread.pause = 1;
}
while (multithread.pause);
}
*/
}
}
#endif

View File

@ -175,7 +175,7 @@ namespace netgen
constexpr operator int () const { return i; }
PointIndex operator++ (int) { PointIndex hi(*this); i++; return hi; }
PointIndex operator-- (int) { PointIndex hi(*this); i--; return hi; }
PointIndex operator++ () { i++; return *this; }
PointIndex & operator++ () { i++; return *this; }
PointIndex operator-- () { i--; return *this; }
void Invalidate() { i = PointIndex::BASE-1; }
bool IsValid() const { return i != PointIndex::BASE-1; }

View File

@ -64,7 +64,7 @@ void STLMeshing (STLGeometry & geom,
STLGeometry :: ~STLGeometry()
{
for (auto p : atlas) delete p;
// for (auto p : atlas) delete p;
delete edgedata;
delete ref;
}
@ -1389,7 +1389,7 @@ void STLGeometry :: DestroyDirtyTrigs()
{
for (k = i+1; k <= GetNT(); k++)
{
trias.Elem(k-1) = trias.Get(k);
trias[k-1] = trias[k];
// readtrias: not longer permanent, JS
// readtrias.Elem(k-1) = readtrias.Get(k);
}
@ -1495,8 +1495,8 @@ void STLGeometry :: PrintSelectInfo()
//int p = GetTriangle(trig).PNum(GetNodeOfSelTrig());
PrintMessage(1,"touch triangle ", GetSelectTrig()
, ", local node ", GetNodeOfSelTrig()
, " (=", GetTriangle(GetSelectTrig()).PNum(GetNodeOfSelTrig()), ")");
, ", local node ", GetNodeOfSelTrig()
, " (=", int(GetTriangle(GetSelectTrig()).PNum(GetNodeOfSelTrig())), ")");
if (AtlasMade() && GetSelectTrig() >= 1 && GetSelectTrig() <= GetNT())
{
PrintMessage(1," chartnum=",GetChartNr(GetSelectTrig()));
@ -1533,11 +1533,11 @@ void STLGeometry :: ShowSelectedTrigCoords()
if (st >= 1 && st <= GetNT())
{
PrintMessage(1, "coordinates of selected trig ", st, ":");
PrintMessage(1, " p1 = ", GetTriangle(st).PNum(1), " = ",
PrintMessage(1, " p1 = ", int(GetTriangle(st).PNum(1)), " = ",
Point3d (GetPoint(GetTriangle(st).PNum(1))));
PrintMessage(1, " p2 = ", GetTriangle(st).PNum(2), " = ",
PrintMessage(1, " p2 = ", int(GetTriangle(st).PNum(2)), " = ",
Point3d (GetPoint(GetTriangle(st).PNum(2))));
PrintMessage(1, " p3 = ", GetTriangle(st).PNum(3), " = ",
PrintMessage(1, " p3 = ", int(GetTriangle(st).PNum(3)), " = ",
Point3d (GetPoint(GetTriangle(st).PNum(3))));
}
}
@ -3111,10 +3111,10 @@ void STLGeometry :: BuildSmoothEdges ()
int STLGeometry :: IsSmoothEdge (int pi1, int pi2) const
bool STLGeometry :: IsSmoothEdge (int pi1, int pi2) const
{
if (!smoothedges)
return 0;
return false;
INDEX_2 i2(pi1, pi2);
i2.Sort();
return smoothedges->Used (i2);
@ -3138,7 +3138,7 @@ int IsInArray(int n, const NgArray<int>& ia)
void STLGeometry :: AddConeAndSpiralEdges(const STLParameters& stlparam)
{
PrintMessage(5,"have now ", GetNE(), " edges with yellow angle = ", stlparam.yangle, " degree");
PrintFnStart("AddConeAndSpiralEdges");
int i,j,k,n;

View File

@ -23,6 +23,11 @@
#include <meshing.hpp>
#include "stltopology.hpp"
#include "stltool.hpp"
#include "stlline.hpp"
namespace netgen
{
@ -44,10 +49,6 @@ namespace netgen
#include "stltopology.hpp"
#include "stltool.hpp"
#include "stlline.hpp"
@ -138,9 +139,9 @@ namespace netgen
//spiralpoints:
NgArray<int> spiralpoints;
//
NgArray<STLChart*> atlas;
Array<unique_ptr<STLChart>, ChartId> atlas;
//marks all already charted trigs with chartnumber
NgArray<int> chartmark;
NgArray<ChartId> chartmark;
//outerchartspertrig, ascending sorted
TABLE<int> outerchartspertrig;
@ -367,8 +368,8 @@ namespace netgen
void AddConeAndSpiralEdges(const STLParameters& stlparam);
void AddFaceEdges(); //each face should have at least one starting edge (outherwise it won't be meshed)
void GetDirtyChartTrigs(int chartnum, STLChart& chart, const NgArray<int>& outercharttrigs,
NgArray<int>& chartpointchecked, NgArray<int>& dirtytrigs);
void GetDirtyChartTrigs(int chartnum, STLChart& chart, const NgArray<ChartId>& outercharttrigs,
NgArray<ChartId>& chartpointchecked, NgArray<int>& dirtytrigs);
void ClearSpiralPoints();
void SetSpiralPoint(int pn) {spiralpoints.Elem(pn) = 1;};
@ -378,7 +379,7 @@ namespace netgen
// smooth edges: sharp geometric edges not declared as edges
void BuildSmoothEdges ();
int IsSmoothEdge (int pi1, int pi2) const;
bool IsSmoothEdge (int pi1, int pi2) const;
//make charts with regions of a max. angle
@ -394,16 +395,16 @@ namespace netgen
//get chart number of a trig or 0 if unmarked
int GetChartNr(int i) const;
int GetMarker(int i) const
ChartId GetMarker(int i) const
{ return chartmark.Get(i); }
void SetMarker(int nr, int m);
int GetNOCharts() const;
void SetMarker(int nr, ChartId m);
int GetNOCharts() const { return atlas.Size(); }
//get a chart from atlas
const STLChart& GetChart(int nr) const;
STLChart& GetChart(int nr) {return *(atlas.Get(nr));};
const STLChart& GetChart(ChartId nr) const { return *atlas[nr];};
STLChart & GetChart(ChartId nr) { return *atlas[nr];};
int AtlasMade() const;
void GetInnerChartLimes(NgArray<twoint>& limes, int chartnum);
void GetInnerChartLimes(NgArray<twoint>& limes, ChartId chartnum);
//FOR MESHING
int GetMeshChartNr () { return meshchart; }
@ -452,7 +453,7 @@ namespace netgen
DLL_HEADER void RestrictLocalH(class Mesh & mesh, double gh, const STLParameters& stlparam);
void RestrictLocalHCurv(class Mesh & mesh, double gh, const STLParameters& stlparam);
void RestrictHChartDistOneChart(int chartnum, NgArray<int>& acttrigs, class Mesh & mesh,
void RestrictHChartDistOneChart(ChartId chartnum, NgArray<int>& acttrigs, class Mesh & mesh,
double gh, double fact, double minh, const STLParameters& stlparam);
friend class MeshingSTLSurface;

View File

@ -69,10 +69,10 @@ void STLGeometry :: MakeAtlas(Mesh & mesh, const MeshingParameters& mparam, cons
double sinchartangle = sin(chartangle);
double sinouterchartangle = sin(outerchartangle);
NgArray<int> outermark(GetNT()); //marks all trigs form actual outer region
NgArray<int> outertested(GetNT()); //marks tested trigs for outer region
NgArray<int> pointstochart(GetNP()); //point in chart becomes chartnum
NgArray<int> innerpointstochart(GetNP()); //point in chart becomes chartnum
NgArray<ChartId> outermark(GetNT()); //marks all trigs form actual outer region
NgArray<ChartId> outertested(GetNT()); //marks tested trigs for outer region
NgArray<ChartId> pointstochart(GetNP()); //point in chart becomes chartnum
NgArray<ChartId> innerpointstochart(GetNP()); //point in chart becomes chartnum
NgArray<int> chartpoints; //point in chart becomes chartnum
NgArray<int> innerchartpoints;
NgArray<Point<3>> innerchartpts;
@ -86,38 +86,32 @@ void STLGeometry :: MakeAtlas(Mesh & mesh, const MeshingParameters& mparam, cons
//int chartboundarydivisions = 10;
markedsegs.SetSize(0); //for testing!!!
NgArray<int> chartpointchecked(GetNP()); //for dirty-chart-trigs
NgArray<ChartId> chartpointchecked(GetNP()); //for dirty-chart-trigs
pointstochart.SetSize(GetNP());
innerpointstochart.SetSize(GetNP());
chartmark.SetSize(GetNT());
for (int i = 1; i <= GetNP(); i++)
{
innerpointstochart.Elem(i) = 0;
pointstochart.Elem(i) = 0;
chartpointchecked.Elem(i) = 0;
}
innerpointstochart = ChartId::INVALID;
pointstochart = ChartId::INVALID;
chartpointchecked = ChartId::INVALID;
double eps = 1e-12 * Dist (boundingbox.PMin(), boundingbox.PMax());
int spiralcheckon = stldoctor.spiralcheck;
if (!spiralcheckon) {PrintWarning("++++++++++++\nspiral deactivated by user!!!!\n+++++++++++++++"); }
chartmark = 0;
outermark = 0;
outertested = 0;
chartmark = ChartId::INVALID;
outermark = ChartId::INVALID;
outertested = ChartId::INVALID;
double atlasarea = Area();
double workedarea = 0;
double showinc = 100.*5000./(double)GetNT();
double nextshow = 0;
// Point<3> startp;
int lastunmarked = 1;
PrintMessage(5,"one dot per 5000 triangles: ");
int markedtrigcnt = 0;
size_t markedtrigcnt = 0;
while (markedtrigcnt < GetNT())
{
if (multithread.terminate) { PopStatus(); return; }
@ -128,11 +122,9 @@ void STLGeometry :: MakeAtlas(Mesh & mesh, const MeshingParameters& mparam, cons
SetThreadPercent(100.0 * workedarea / atlasarea);
STLChart * chart = new STLChart(this, stlparam);
atlas.Append(chart);
// *testout << "Chart " << atlas.Size() << endl;
atlas.Append (make_unique<STLChart> (this, stlparam));
STLChart & chart = *atlas.Last();
//find unmarked trig
int prelastunmarked = lastunmarked;
@ -149,29 +141,29 @@ void STLGeometry :: MakeAtlas(Mesh & mesh, const MeshingParameters& mparam, cons
innerchartpoints.SetSize(0);
innerchartpts.SetSize(0);
chartbound.Clear();
chartbound.SetChart(chart);
chartbound.SetChart(&chart);
chartbound.BuildSearchTree(); // different !!!
if (!found) { PrintSysError("Make Atlas, no starttrig found"); return; }
if (!found) { throw Exception("Make Atlas, no starttrig found"); }
//find surrounding trigs
// int starttrig = j;
int starttrig = lastunmarked;
Point<3> startp = GetPoint(GetTriangle(starttrig).PNum(1));
Point<3> startp = GetPoint(GetTriangle(starttrig)[0]);
int accepted;
int chartnum = GetNOCharts();
bool accepted;
ChartId chartnum = GetNOCharts();
Vec<3> sn = GetTriangle(starttrig).Normal();
chart->SetNormal (startp, sn);
chart.SetNormal (startp, sn);
// *testout << "first trig " << starttrig << ", n = " << sn << endl;
SetMarker(starttrig, chartnum);
markedtrigcnt++;
chart->AddChartTrig(starttrig);
chart.AddChartTrig(starttrig);
chartbound.AddTriangle(GetTriangle(starttrig));
workedarea += GetTriangle(starttrig).Area(points);
@ -195,11 +187,11 @@ void STLGeometry :: MakeAtlas(Mesh & mesh, const MeshingParameters& mparam, cons
{
changed = false;
oldstartic2 = oldstartic;
oldstartic = chart->GetNT();
oldstartic = chart.GetNT();
// for (ic = oldstartic2; ic <= chart->GetNT(); ic++)
for (int ic = oldstartic2; ic <= oldstartic; ic++)
{
int i = chart->GetTrig(ic);
int i = chart.GetTrig(ic);
if (GetMarker(i) == chartnum)
{
for (int j = 1; j <= NONeighbourTrigs(i); j++)
@ -215,7 +207,7 @@ void STLGeometry :: MakeAtlas(Mesh & mesh, const MeshingParameters& mparam, cons
if ( (n2 * sn) >= coschartangle )
{
// *testout << "good angle " << endl;
accepted = 1;
accepted = true;
/*
//alter spiralentest, schnell, aber ungenau
for (k = 1; k <= 3; k++)
@ -261,7 +253,7 @@ void STLGeometry :: MakeAtlas(Mesh & mesh, const MeshingParameters& mparam, cons
Vec<3> n3 = GetTriangle(nnt).Normal();
if ( (n3 * sn) >= coschartangle &&
IsSmoothEdge (nnp1, nnp2) )
accepted = 1;
accepted = true;
}
if (!accepted)
break;
@ -275,7 +267,7 @@ void STLGeometry :: MakeAtlas(Mesh & mesh, const MeshingParameters& mparam, cons
changed = true;
markedtrigcnt++;
workedarea += GetTriangle(nt).Area(points);
chart->AddChartTrig(nt);
chart.AddChartTrig(nt);
chartbound.AddTriangle(GetTriangle(nt));
@ -320,11 +312,11 @@ void STLGeometry :: MakeAtlas(Mesh & mesh, const MeshingParameters& mparam, cons
{
changed = false;
oldstartic2 = oldstartic;
oldstartic = chart->GetNT();
oldstartic = chart.GetNT();
for (int ic = oldstartic2; ic <= oldstartic; ic++)
{
int i = chart->GetTrig(ic);
int i = chart.GetTrig(ic);
if (outermark.Get(i) != chartnum) continue;
for (int j = 1; j <= NONeighbourTrigs(i); j++)
@ -350,7 +342,7 @@ void STLGeometry :: MakeAtlas(Mesh & mesh, const MeshingParameters& mparam, cons
//abfragen, ob noch im tolerierten Winkel
if ( (n2 * sn) >= cosouterchartangle )
{
accepted = 1;
accepted = true;
// NgProfiler::StartTimer (timer4);
bool isdirtytrig = false;
@ -392,7 +384,7 @@ void STLGeometry :: MakeAtlas(Mesh & mesh, const MeshingParameters& mparam, cons
Vec<3> n3 = GetTriangle(nnt).Normal();
if ( (n3 * sn) >= cosouterchartangle &&
IsSmoothEdge (nnp1, nnp2) )
accepted = 1;
accepted = true;
// NgProfiler::StopTimer (timer4e);
}
// NgProfiler::StopTimer (timer4b);
@ -411,12 +403,12 @@ void STLGeometry :: MakeAtlas(Mesh & mesh, const MeshingParameters& mparam, cons
if (accepted)
{
// NgProfiler::StartTimer (timer5a);
accepted = 0;
accepted = false;
for (int k = 1; k <= 3; k++)
if (innerpointstochart.Get(ntrig.PNum(k)) == chartnum)
{
accepted = 1;
accepted = true;
break;
}
@ -445,7 +437,7 @@ void STLGeometry :: MakeAtlas(Mesh & mesh, const MeshingParameters& mparam, cons
double tdist = Dist2(pt, innerchartpts[l]);
if (tdist < 4 * h2)
{
accepted = 1;
accepted = true;
break;
}
}
@ -464,7 +456,7 @@ void STLGeometry :: MakeAtlas(Mesh & mesh, const MeshingParameters& mparam, cons
if (GetMarker(nt) != chartnum)
{
chartbound.AddTriangle(GetTriangle(nt));
chart->AddOuterTrig(nt);
chart.AddOuterTrig(nt);
for (int k = 1; k <= 3; k++)
{
if (pointstochart.Get(GetTriangle(nt).PNum(k))
@ -486,26 +478,26 @@ void STLGeometry :: MakeAtlas(Mesh & mesh, const MeshingParameters& mparam, cons
// NgProfiler::StartTimer (timere);
// NgProfiler::StartTimer (timere1);
//end of while loop for outer chart
GetDirtyChartTrigs(chartnum, *chart, outermark, chartpointchecked, dirtycharttrigs);
GetDirtyChartTrigs(chartnum, chart, outermark, chartpointchecked, dirtycharttrigs);
//dirtycharttrigs are local (chart) point numbers!!!!!!!!!!!!!!!!
if (dirtycharttrigs.Size() != 0 &&
(dirtycharttrigs.Size() != chart->GetNChartT() || dirtycharttrigs.Size() != 1))
(dirtycharttrigs.Size() != chart.GetNChartT() || dirtycharttrigs.Size() != 1))
{
if (dirtycharttrigs.Size() == chart->GetNChartT() && dirtycharttrigs.Size() != 1)
if (dirtycharttrigs.Size() == chart.GetNChartT() && dirtycharttrigs.Size() != 1)
{
//if all trigs would be eliminated -> leave 1 trig!
dirtycharttrigs.SetSize(dirtycharttrigs.Size() - 1);
}
for (int k = 1; k <= dirtycharttrigs.Size(); k++)
{
int tn = chart->GetChartTrig(dirtycharttrigs.Get(k));
int tn = chart.GetChartTrig(dirtycharttrigs.Get(k));
outermark.Elem(tn) = 0; //not necessary, for later use
SetMarker(tn, 0);
markedtrigcnt--;
workedarea -= GetTriangle(tn).Area(points);
}
chart->MoveToOuterChart(dirtycharttrigs);
chart.MoveToOuterChart(dirtycharttrigs);
lastunmarked = 1;
lastunmarked = prelastunmarked;
}
@ -537,7 +529,8 @@ void STLGeometry :: MakeAtlas(Mesh & mesh, const MeshingParameters& mparam, cons
int cnttrias = 0;
outerchartspertrig.SetSize(GetNT());
for (int i = 1; i <= atlas.Size(); i++)
// for (int i = 1; i <= atlas.Size(); i++)
for (ChartId i : atlas.Range())
{
for (int j = 1; j <= GetChart(i).GetNT(); j++)
{
@ -633,23 +626,10 @@ int STLGeometry :: GetMarker(int i) const
return chartmark.Get(i);
}
*/
void STLGeometry :: SetMarker(int nr, int m)
void STLGeometry :: SetMarker(int nr, ChartId m)
{
chartmark.Elem(nr) = m;
}
int STLGeometry :: GetNOCharts() const
{
return atlas.Size();
}
const STLChart& STLGeometry :: GetChart(int nr) const
{
if (nr > atlas.Size())
{
PrintSysError("GetChart(", nr, ") not possible!!!");
nr = 1;
}
return *(atlas.Get(nr));
}
int STLGeometry :: AtlasMade() const
{
@ -671,7 +651,7 @@ int AddIfNotExists(NgArray<int>& list, int x)
}
*/
void STLGeometry :: GetInnerChartLimes(NgArray<twoint>& limes, int chartnum)
void STLGeometry :: GetInnerChartLimes(NgArray<twoint>& limes, ChartId chartnum)
{
int j, k;
@ -732,8 +712,8 @@ void STLGeometry :: GetInnerChartLimes(NgArray<twoint>& limes, int chartnum)
void STLGeometry :: GetDirtyChartTrigs(int chartnum, STLChart& chart,
const NgArray<int>& outercharttrigs,
NgArray<int>& chartpointchecked,
const NgArray<ChartId>& outercharttrigs,
NgArray<ChartId>& chartpointchecked,
NgArray<int>& dirtytrigs)
{
dirtytrigs.SetSize(0);

View File

@ -1116,7 +1116,7 @@ void STLGeometry :: RestrictLocalH(class Mesh & mesh, double gh, const STLParame
}
}
void STLGeometry :: RestrictHChartDistOneChart(int chartnum, NgArray<int>& acttrigs,
void STLGeometry :: RestrictHChartDistOneChart(ChartId chartnum, NgArray<int>& acttrigs,
class Mesh & mesh, double gh, double fact, double minh,
const STLParameters& stlparam)
{

View File

@ -581,59 +581,57 @@ int STLLine :: GetRightTrig(int nr) const
return righttrigs.Get(nr);
};
double STLLine :: GetSegLen(const NgArray<Point<3> >& ap, int nr) const
double STLLine :: GetSegLen(const Array<Point<3>,STLPointId>& ap, int nr) const
{
return Dist(ap.Get(PNum(nr)),ap.Get(PNum(nr+1)));
return Dist(ap[PNum(nr)],ap[PNum(nr+1)]);
}
double STLLine :: GetLength(const NgArray<Point<3> >& ap) const
double STLLine :: GetLength(const Array<Point<3>,STLPointId>& ap) const
{
double len = 0;
for (int i = 2; i <= pts.Size(); i++)
{
len += (ap.Get(pts.Get(i)) - ap.Get(pts.Get(i-1))).Length();
}
len += (ap[pts.Get(i)] - ap[pts.Get(i-1)]).Length();
return len;
}
void STLLine :: GetBoundingBox (const NgArray<Point<3> > & ap, Box<3> & box) const
void STLLine :: GetBoundingBox (const Array<Point<3>,STLPointId> & ap, Box<3> & box) const
{
box.Set (ap.Get (pts[0]));
box.Set (ap[pts[0]]);
for (int i = 1; i < pts.Size(); i++)
box.Add (ap.Get(pts[i]));
box.Add (ap[pts[i]]);
}
Point<3> STLLine ::
GetPointInDist(const NgArray<Point<3> >& ap, double dist, int& index) const
GetPointInDist(const Array<Point<3>,STLPointId>& ap, double dist, int& index) const
{
if (dist <= 0)
{
index = 1;
return ap.Get(StartP());
return ap[StartP()];
}
double len = 0;
int i;
for (i = 1; i < pts.Size(); i++)
{
double seglen = Dist (ap.Get(pts.Get(i)),
ap.Get(pts.Get(i+1)));
double seglen = Dist (ap[pts.Get(i)],
ap[pts.Get(i+1)]);
if (len + seglen > dist)
{
index = i;
double relval = (dist - len) / (seglen + 1e-16);
Vec3d v (ap.Get(pts.Get(i)), ap.Get(pts.Get(i+1)));
return ap.Get(pts.Get(i)) + relval * v;
Vec3d v (ap[pts.Get(i)], ap[pts.Get(i+1)]);
return ap[pts.Get(i)] + relval * v;
}
len += seglen;
}
index = pts.Size() - 1;
return ap.Get(EndP());
return ap[EndP()];
}
@ -644,7 +642,7 @@ double GetH(const Point3d& p, double x)
return stlgh;//+0.5)*(x+0.5);
}
*/
STLLine* STLLine :: Mesh(const NgArray<Point<3> >& ap,
STLLine* STLLine :: Mesh(const Array<Point<3>,STLPointId>& ap,
NgArray<Point3d>& mp, double ghi,
class Mesh& mesh) const
{
@ -720,7 +718,7 @@ STLLine* STLLine :: Mesh(const NgArray<Point<3> >& ap,
int j = 1;
p = ap.Get(StartP());
p = ap[StartP()];
int pn = AddPointIfNotExists(mp, p, 1e-10*diam);
int segn = 1;
@ -773,7 +771,7 @@ STLLine* STLLine :: Mesh(const NgArray<Point<3> >& ap,
NgProfiler::StartTimer (timer3);
p = ap.Get(EndP());
p = ap[EndP()];
pn = AddPointIfNotExists(mp, p, 1e-10*diam);
segn = GetNS();
line->AddPoint(pn);

View File

@ -9,6 +9,8 @@
/* Date: 20. Nov. 99 */
/**************************************************************************/
namespace netgen {
class STLGeometry;
class STLTopology;
@ -158,11 +160,11 @@ public:
int NP() const {return pts.Size();}
int GetNS() const;
void GetSeg(int nr, int& p1, int& p2) const;
double GetSegLen(const NgArray<Point<3> >& ap, int nr) const;
double GetSegLen(const Array<Point<3>,STLPointId>& ap, int nr) const;
int GetLeftTrig(int nr) const;
int GetRightTrig(int nr) const;
double GetDist(int nr) const { return dists.Get(nr);};
void GetBoundingBox (const NgArray<Point<3> > & ap, Box<3> & box) const;
void GetBoundingBox (const Array<Point<3>,STLPointId> & ap, Box<3> & box) const;
void AddLeftTrig(int nr) {lefttrigs.Append(nr);}
void AddRightTrig(int nr) {righttrigs.Append(nr);}
@ -170,14 +172,14 @@ public:
int StartP() const {return pts.Get(1);}
int EndP() const {return pts.Get(pts.Size());}
double GetLength(const NgArray<Point<3> >& ap) const;
double GetLength(const Array<Point<3>,STLPointId>& ap) const;
//suche punkt in entfernung (in linienkoordinaten) dist
//in index ist letzter punkt VOR dist (d.h. max pts.Size()-1)
Point<3> GetPointInDist(const NgArray<Point<3> >& ap, double dist, int& index) const;
Point<3> GetPointInDist(const Array<Point<3>,STLPointId>& ap, double dist, int& index) const;
//return a meshed polyline
STLLine* Mesh(const NgArray<Point<3> >& ap,
STLLine* Mesh(const Array<Point<3>,STLPointId>& ap,
NgArray<Point3d>& mp, double ghi,
class Mesh& mesh) const;
@ -185,4 +187,7 @@ public:
int ShouldSplit() const {return split;}
};
} // namespace netgen
#endif

View File

@ -282,7 +282,7 @@ STLReadTriangle :: STLReadTriangle (const Point<3> * apts,
STLTriangle :: STLTriangle(const int * apts)
STLTriangle :: STLTriangle(const STLPointId * apts)
{
pts[0] = apts[0];
pts[1] = apts[1];
@ -349,11 +349,11 @@ int STLTriangle :: GetNeighbourPointsAndOpposite(const STLTriangle& t, int& p1,
return 0;
}
Vec<3> STLTriangle :: GeomNormal(const NgArray<Point<3> >& ap) const
Vec<3> STLTriangle :: GeomNormal(const Array<Point<3>,STLPointId>& ap) const
{
const Point<3> & p1 = ap.Get(PNum(1));
const Point<3> & p2 = ap.Get(PNum(2));
const Point<3> & p3 = ap.Get(PNum(3));
const Point<3> & p1 = ap[PNum(1)];
const Point<3> & p2 = ap[PNum(2)];
const Point<3> & p3 = ap[PNum(3)];
return Cross(p2-p1, p3-p1);
}
@ -382,13 +382,13 @@ void STLTriangle :: ChangeOrientation()
double STLTriangle :: Area(const NgArray<Point<3> >& ap) const
double STLTriangle :: Area(const Array<Point<3>,STLPointId>& ap) const
{
return 0.5 * Cross(ap.Get(PNum(2))-ap.Get(PNum(1)),
ap.Get(PNum(3))-ap.Get(PNum(1))).Length();
return 0.5 * Cross(ap[PNum(2)]-ap[PNum(1)],
ap[PNum(3)]-ap[PNum(1)]).Length();
}
double STLTriangle :: MinHeight(const NgArray<Point<3> >& ap) const
double STLTriangle :: MinHeight(const Array<Point<3>,STLPointId>& ap) const
{
double ml = MaxLength(ap);
if (ml != 0) {return 2.*Area(ap)/ml;}
@ -396,19 +396,19 @@ double STLTriangle :: MinHeight(const NgArray<Point<3> >& ap) const
return 0;
}
double STLTriangle :: MaxLength(const NgArray<Point<3> >& ap) const
double STLTriangle :: MaxLength(const Array<Point<3>,STLPointId>& ap) const
{
return max3(Dist(ap.Get(PNum(1)),ap.Get(PNum(2))),
Dist(ap.Get(PNum(2)),ap.Get(PNum(3))),
Dist(ap.Get(PNum(3)),ap.Get(PNum(1))));
return max3(Dist(ap[PNum(1)],ap[PNum(2)]),
Dist(ap[PNum(2)],ap[PNum(3)]),
Dist(ap[PNum(3)],ap[PNum(1)]));
}
void STLTriangle :: ProjectInPlain(const NgArray<Point<3> >& ap,
void STLTriangle :: ProjectInPlain(const Array<Point<3>,STLPointId>& ap,
const Vec<3> & n, Point<3> & pp) const
{
const Point<3> & p1 = ap.Get(PNum(1));
const Point<3> & p2 = ap.Get(PNum(2));
const Point<3> & p3 = ap.Get(PNum(3));
const Point<3> & p1 = ap[PNum(1)];
const Point<3> & p2 = ap[PNum(2)];
const Point<3> & p3 = ap[PNum(3)];
Vec<3> v1 = p2 - p1;
Vec<3> v2 = p3 - p1;
@ -430,13 +430,13 @@ void STLTriangle :: ProjectInPlain(const NgArray<Point<3> >& ap,
}
int STLTriangle :: ProjectInPlain (const NgArray<Point<3> >& ap,
int STLTriangle :: ProjectInPlain (const Array<Point<3>,STLPointId>& ap,
const Vec<3> & nproj,
Point<3> & pp, Vec<3> & lam) const
{
const Point<3> & p1 = ap.Get(PNum(1));
const Point<3> & p2 = ap.Get(PNum(2));
const Point<3> & p3 = ap.Get(PNum(3));
const Point<3> & p1 = ap[PNum(1)];
const Point<3> & p2 = ap[PNum(2)];
const Point<3> & p3 = ap[PNum(3)];
Vec<3> v1 = p2-p1;
Vec<3> v2 = p3-p1;
@ -468,12 +468,12 @@ int STLTriangle :: ProjectInPlain (const NgArray<Point<3> >& ap,
void STLTriangle :: ProjectInPlain(const NgArray<Point<3> >& ap,
void STLTriangle :: ProjectInPlain(const Array<Point<3>,STLPointId>& ap,
Point<3> & pp) const
{
const Point<3> & p1 = ap.Get(PNum(1));
const Point<3> & p2 = ap.Get(PNum(2));
const Point<3> & p3 = ap.Get(PNum(3));
const Point<3> & p1 = ap[PNum(1)];
const Point<3> & p2 = ap[PNum(2)];
const Point<3> & p3 = ap[PNum(3)];
Vec<3> v1 = p2 - p1;
Vec<3> v2 = p3 - p1;
@ -488,12 +488,12 @@ void STLTriangle :: ProjectInPlain(const NgArray<Point<3> >& ap,
pp = pp + (nfact) * nt;
}
int STLTriangle :: PointInside(const NgArray<Point<3> > & ap,
int STLTriangle :: PointInside(const Array<Point<3>,STLPointId> & ap,
const Point<3> & pp) const
{
const Point<3> & p1 = ap.Get(PNum(1));
const Point<3> & p2 = ap.Get(PNum(2));
const Point<3> & p3 = ap.Get(PNum(3));
const Point<3> & p1 = ap[PNum(1)];
const Point<3> & p2 = ap[PNum(2)];
const Point<3> & p3 = ap[PNum(3)];
Vec<3> v1 = p2 - p1;
Vec<3> v2 = p3 - p1;
@ -532,7 +532,7 @@ int STLTriangle :: PointInside(const NgArray<Point<3> > & ap,
return 0;
}
double STLTriangle :: GetNearestPoint(const NgArray<Point<3> >& ap,
double STLTriangle :: GetNearestPoint(const Array<Point<3>,STLPointId>& ap,
Point<3> & p3d) const
{
Point<3> p = p3d;
@ -548,7 +548,7 @@ double STLTriangle :: GetNearestPoint(const NgArray<Point<3> >& ap,
for (int j = 1; j <= 3; j++)
{
p = p3d;
dist = GetDistFromLine(ap.Get(PNum(j)), ap.Get(PNumMod(j+1)), p);
dist = GetDistFromLine(ap[PNum(j)], ap[PNumMod(j+1)], p);
if (dist < nearest)
{
nearest = dist;
@ -690,15 +690,15 @@ void STLChart :: AddOuterTrig(int i)
{searchtree->Insert (pmin, pmax, i);}
}
int STLChart :: IsInWholeChart(int nr) const
bool STLChart :: IsInWholeChart(int nr) const
{
for (int i = 1; i <= charttrigs.Size(); i++)
if (charttrigs.Get(i) == nr) return 1;
if (charttrigs.Get(i) == nr) return true;
for (int i = 1; i <= outertrigs.Size(); i++)
if (outertrigs.Get(i) == nr) return 1;
if (outertrigs.Get(i) == nr) return true;
return 0;
return false;
}
void STLChart :: GetTrianglesInBox (const Point3d & pmin,
@ -1040,8 +1040,8 @@ void STLBoundary ::AddTriangle(const STLTriangle & t)
// NgProfiler::StopTimer (timer_new);
}
int STLBoundary :: TestSeg(const Point<3>& p1, const Point<3> & p2, const Vec<3> & sn,
double sinchartangle, int divisions, NgArray<Point<3> >& points, double eps)
bool STLBoundary :: TestSeg(const Point<3>& p1, const Point<3> & p2, const Vec<3> & sn,
double sinchartangle, int divisions, Array<Point<3>,STLPointId>& points, double eps)
{
if (usechartnormal)
return TestSegChartNV (p1, p2, sn);
@ -1276,7 +1276,7 @@ void STLBoundary :: DeleteSearchTree()
}
// checks, whether 2d projection intersects
int STLBoundary :: TestSegChartNV(const Point3d & p1, const Point3d& p2,
bool STLBoundary :: TestSegChartNV(const Point3d & p1, const Point3d& p2,
const Vec3d& sn)
{
// static int timerquick = NgProfiler::CreateTimer ("TestSegChartNV-searchtree");

View File

@ -1,7 +1,6 @@
#ifndef FILE_STLTOOL
#define FILE_STLTOOL
//#include "gprim/gprim.hh"
/**************************************************************************/
@ -11,7 +10,7 @@
/* Date: 20. Nov. 99 */
/**************************************************************************/
namespace netgen {
// use one normal vector for whole chart
extern int usechartnormal;
@ -40,6 +39,32 @@ typedef NgArray <int> * ArrayINTPTR;
class STLGeometry;
class STLParameters;
// typedef int ChartId
class ChartId
{
int i;
public:
class t_invalid { public: constexpr t_invalid() = default; };
static constexpr t_invalid INVALID{};
ChartId() { }
constexpr ChartId(t_invalid inv) : i(0) { ; }
constexpr ChartId(int ai) : i(ai) { }
operator int() const { return i; }
ChartId operator++ (int) { ChartId hi(*this); i++; return hi; }
ChartId & operator++ () { i++; return *this; }
};
}
namespace ngcore
{
template<>
constexpr netgen::ChartId IndexBASE<netgen::ChartId> () { return netgen::ChartId(1); }
}
namespace netgen {
class STLChart
{
private:
@ -60,7 +85,7 @@ public:
void AddChartTrig(int i);
void AddOuterTrig(int i);
int IsInWholeChart(int nr) const;
bool IsInWholeChart(int nr) const;
int GetChartTrig(int i) const {return charttrigs.Get(i);}
int GetOuterTrig(int i) const {return outertrigs.Get(i);}
@ -118,13 +143,13 @@ class STLBoundarySeg
// Point<2> p2dmin, p2dmax;
double rad;
int i1, i2;
STLPointId i1, i2;
int smoothedge;
public:
STLBoundarySeg () { ; }
STLBoundarySeg (int ai1, int ai2, const NgArray<Point<3> > & points,
STLBoundarySeg (STLPointId ai1, STLPointId ai2, const Array<Point<3>,STLPointId> & points,
const STLChart * chart)
: p1(points.Get(ai1)), p2(points.Get(ai2)),
: p1(points[ai1]), p2(points[ai2]),
i1(ai1), i2(ai2)
{
center = ::netgen::Center (p1, p2);
@ -182,11 +207,11 @@ public:
void BuildSearchTree();
void DeleteSearchTree();
int TestSeg(const Point<3> & p1, const Point<3> & p2, const Vec<3> & sn,
double sinchartangle, int divisions, NgArray<Point<3> >& points,
double eps);
int TestSegChartNV(const Point3d& p1, const Point3d& p2, const Vec3d& sn);
bool TestSeg(const Point<3> & p1, const Point<3> & p2, const Vec<3> & sn,
double sinchartangle, int divisions, Array<Point<3>,STLPointId>& points,
double eps);
bool TestSegChartNV(const Point3d& p1, const Point3d& p2, const Vec3d& sn);
};
@ -292,6 +317,6 @@ void STLSurfaceOptimization (STLGeometry & geom,
const MeshingParameters & mparam);
} // namespace netgen
#endif

View File

@ -515,9 +515,9 @@ void STLTopology :: InitSTLGeometry(const NgArray<STLReadTriangle> & readtrigs)
foundpos = AddPoint(p);
pointtree->Insert (p, foundpos);
}
if (Dist(p, points.Get(foundpos)) > 1e-10)
cout << "identify close points: " << p << " " << points.Get(foundpos)
<< ", dist = " << Dist(p, points.Get(foundpos))
if (Dist(p, points[foundpos]) > 1e-10)
cout << "identify close points: " << p << " " << points[foundpos]
<< ", dist = " << Dist(p, points[foundpos])
<< endl;
st[k] = foundpos;
}
@ -711,14 +711,14 @@ void STLTopology :: FindNeighbourTrigs()
for (STLTrigIndex ti = 0; ti < GetNT(); ti++)
for (STLTrigId ti = 0; ti < GetNT(); ti++)
{
STLTriangle & trig = trias[ti];
for (int k = 0; k < 3; k++)
{
STLPointIndex pi = trig[k] - STLBASE;
STLPointIndex pi2 = trig[(k+1)%3] - STLBASE;
STLPointIndex pi3 = trig[(k+2)%3] - STLBASE;
STLPointId pi = trig[k]; // - STLBASE;
STLPointId pi2 = trig[(k+1)%3]; // - STLBASE;
STLPointId pi3 = trig[(k+2)%3]; // - STLBASE;
// vector along edge
Vec<3> ve = points[pi2] - points[pi];
@ -736,24 +736,24 @@ void STLTopology :: FindNeighbourTrigs()
for (int j = 0; j < trigsperpoint[pi].Size(); j++)
{
STLTrigIndex ti2 = trigsperpoint[pi][j] - STLBASE;
STLTrigId ti2 = trigsperpoint[pi][j]; // - STLBASE;
const STLTriangle & trig2 = trias[ti2];
if (ti == ti2) continue;
bool hasboth = 0;
for (int l = 0; l < 3; l++)
if (trig2[l] - STLBASE == pi2)
if (trig2[l] /* - STLBASE */ == pi2)
{
hasboth = 1;
break;
}
if (!hasboth) continue;
STLPointIndex pi4(0);
STLPointId pi4(0);
for (int l = 0; l < 3; l++)
if (trig2[l] - STLBASE != pi && trig2[l] - STLBASE != pi2)
pi4 = trig2[l] - STLBASE;
if (trig2[l] /* - STLBASE */ != pi && trig2[l] /* - STLBASE */ != pi2)
pi4 = trig2[l] /* - STLBASE */;
Vec<3> vt2 = points[pi4] - points[pi];
@ -763,12 +763,12 @@ void STLTopology :: FindNeighbourTrigs()
if (phi < phimin)
{
phimin = phi;
trig.NBTrig (0, (k+2)%3) = ti2 + STLBASE;
trig.NBTrig (0, (k+2)%3) = ti2; // + STLBASE;
}
if (phi > phimax)
{
phimax = phi;
trig.NBTrig (1, (k+2)%3) = ti2 + STLBASE;
trig.NBTrig (1, (k+2)%3) = ti2; // + STLBASE;
}
}
}

View File

@ -14,41 +14,58 @@
*/
namespace netgen {
class STLGeometry;
#define STLBASE 1
// #define STLBASE 1
class STLPointIndex
class STLPointId
{
int i;
public:
STLPointIndex () { ; }
STLPointIndex (int ai) : i(ai) { ; }
STLPointIndex & operator= (const STLPointIndex & ai) { i = ai.i; return *this; }
STLPointIndex & operator= (int ai) { i = ai; return *this; }
STLPointId () { ; }
constexpr STLPointId (int ai) : i(ai) { ; }
STLPointId & operator= (const STLPointId & ai) { i = ai.i; return *this; }
STLPointId & operator= (int ai) { i = ai; return *this; }
operator int () const { return i; }
STLPointIndex operator++ (int) { return i++; }
STLPointIndex operator-- (int) { return i--; }
STLPointId operator++ (int) { return i++; }
STLPointId operator-- (int) { return i--; }
void DoArchive(Archive& ar) { ar & i; }
};
class STLTrigIndex
class STLTrigId
{
int i;
public:
STLTrigIndex () { ; }
STLTrigIndex (int ai) : i(ai) { ; }
STLTrigIndex & operator= (const STLTrigIndex & ai) { i = ai.i; return *this; }
STLTrigIndex & operator= (int ai) { i = ai; return *this; }
STLTrigId () { ; }
constexpr STLTrigId (int ai) : i(ai) { ; }
STLTrigId & operator= (const STLTrigId & ai) { i = ai.i; return *this; }
STLTrigId & operator= (int ai) { i = ai; return *this; }
operator int () const { return i; }
STLTrigIndex operator++ (int) { return i++; }
STLTrigIndex operator-- (int) { return i--; }
STLTrigId operator++ (int) { return i++; }
STLTrigId operator-- (int) { return i--; }
};
}
namespace ngcore
{
template<>
constexpr netgen::STLPointId IndexBASE<netgen::STLPointId> () { return netgen::STLPointId(1); }
template<>
constexpr netgen::STLTrigId IndexBASE<netgen::STLTrigId> () { return netgen::STLTrigId(1); }
}
namespace netgen {
// triangle structure for loading stl files
class STLReadTriangle
@ -73,7 +90,7 @@ class STLTriangle
// normalized stored normal vector ??
Vec<3> normal;
// point numbers of triangle
int pts[3];
STLPointId pts[3];
// front-side and back-side domains
int domains[2];
@ -93,22 +110,23 @@ public:
STLTriangle (const int * apts);
STLTriangle (const STLPointId * apts);
STLTriangle () {pts[0]=0;pts[1]=0;pts[2]=0;}
void DoArchive(Archive& ar)
{
ar.Do(&topedges[0],3);
ar.Do(&nbtrigs[0][0], 6);
ar.Do(&pts[0],3);
// ar.Do(&pts[0],3);
ar & pts[0] & pts[1] & pts[2];
ar.Do(&domains[0],2);
size_t i = flags.toperror;
ar & normal & box & center & rad & facenum & i;
flags.toperror = i;
}
int operator[] (int i) const { return pts[i]; }
int & operator[] (int i) { return pts[i]; }
STLPointId operator[] (int i) const { return pts[i]; }
STLPointId & operator[] (int i) { return pts[i]; }
int EdgeNum(int i) const { return topedges[(i-1)]; }
int & EdgeNum(int i) { return topedges[(i-1)]; }
@ -123,10 +141,10 @@ public:
// obsolete:
int PNum(int i) const { return pts[(i-1)]; }
int & PNum(int i) { return pts[(i-1)]; }
int PNumMod(int i) const { return pts[(i-1)%3]; }
int & PNumMod(int i) { return pts[(i-1)%3]; }
STLPointId PNum(int i) const { return pts[(i-1)]; }
STLPointId & PNum(int i) { return pts[(i-1)]; }
STLPointId PNumMod(int i) const { return pts[(i-1)%3]; }
STLPointId & PNumMod(int i) { return pts[(i-1)%3]; }
int EdgeNumMod(int i) const { return topedges[(i-1)%3]; }
int & EdgeNumMod(int i) { return topedges[(i-1)%3]; }
@ -149,7 +167,7 @@ public:
// NON-normalized geometry - normal vector
Vec<3> GeomNormal(const NgArray<Point<3> >& ap) const;
Vec<3> GeomNormal(const Array<Point<3>,STLPointId>& ap) const;
// Stored normal vector, normalized
void SetNormal (const Vec<3> & n);
@ -159,10 +177,10 @@ public:
void ChangeOrientation();
//project with a certain normal vector in plane
void ProjectInPlain(const NgArray<Point<3> >& ap,
void ProjectInPlain(const Array<Point<3>, STLPointId>& ap,
const Vec<3> & n, Point<3> & pp) const;
//project with the triangle's normal vector in plane
void ProjectInPlain(const NgArray<Point<3> > & ap, Point<3> & pp) const;
void ProjectInPlain(const Array<Point<3>, STLPointId> & ap, Point<3> & pp) const;
/*
@ -177,20 +195,20 @@ public:
pp(output) = P1 + lam1 v1 + lam2 v2
*/
int ProjectInPlain (const NgArray<Point<3> >& ap,
int ProjectInPlain (const Array<Point<3>,STLPointId>& ap,
const Vec<3> & nproj,
Point<3> & pp, Vec<3> & lam) const;
int PointInside(const NgArray<Point<3> >& ap, const Point<3> & pp) const;
int PointInside(const Array<Point<3>,STLPointId>& ap, const Point<3> & pp) const;
//get nearest point on triangle and distance to it
double GetNearestPoint(const NgArray<Point<3> >& ap,
double GetNearestPoint(const Array<Point<3>,STLPointId>& ap,
Point<3> & p3d) const;
double Area(const NgArray<Point<3> >& ap) const;
double Area(const Array<Point<3>,STLPointId>& ap) const;
double MinHeight(const NgArray<Point<3> >& ap) const;
double MaxLength(const NgArray<Point<3> >& ap) const;
double MinHeight(const Array<Point<3>,STLPointId>& ap) const;
double MaxLength(const Array<Point<3>,STLPointId>& ap) const;
//max length of a side of triangle
int GetFaceNum() {return facenum;}
@ -250,9 +268,9 @@ ostream& operator<<(ostream& os, const STLTriangle& t);
class STLTopology
{
protected:
NgArray<STLTriangle> trias;
Array<STLTriangle, STLTrigId> trias;
NgArray<STLTopEdge> topedges;
NgArray<Point<3> > points;
Array<Point<3>, STLPointId> points;
// mapping of sorted pair of points to topedge
INDEX_2_HASHTABLE<int> * ht_topedges;
@ -312,24 +330,24 @@ public:
int GetNP() const { return points.Size(); }
int AddPoint(const Point<3> & p) { points.Append(p); return points.Size(); }
const Point<3> & GetPoint(int nr) const { return points.Get(nr); }
const Point<3> & GetPoint(int nr) const { return points[nr]; } // .Get(nr); }
int GetPointNum (const Point<3> & p);
void SetPoint(int nr, const Point<3> & p) { points.Elem(nr) = p; }
const NgArray<Point<3> >& GetPoints() const { return points; }
void SetPoint(int nr, const Point<3> & p) { points[nr] = p; } // { points.Elem(nr) = p; }
auto & GetPoints() const { return points; }
const Point<3> & operator[] (STLPointIndex i) const { return points[i]; }
Point<3> & operator[] (STLPointIndex i) { return points[i]; }
const Point<3> & operator[] (STLPointId i) const { return points[i]; }
Point<3> & operator[] (STLPointId i) { return points[i]; }
int GetNT() const { return trias.Size(); }
void AddTriangle(const STLTriangle& t);
const STLTriangle & GetTriangle (int nr) const { return trias.Get(nr); }
STLTriangle & GetTriangle (int nr) { return trias.Elem(nr); }
const STLTriangle & GetTriangle (int nr) const { return trias[nr]; } // .Get(nr); }
STLTriangle & GetTriangle (int nr) { return trias[nr]; } // .Elem(nr); }
const STLTriangle & operator[] (STLTrigIndex i) const { return trias[i]; }
STLTriangle & operator[] (STLTrigIndex i) { return trias[i]; }
const STLTriangle & operator[] (STLTrigId i) const { return trias[i]; }
STLTriangle & operator[] (STLTrigId i) { return trias[i]; }
int GetNTE() const { return topedges.Size(); }
@ -376,5 +394,6 @@ public:
const Box<3> & GetBoundingBox () const { return boundingbox; }
};
} // namespace netgen
#endif