ForkMaster/netgen
1
0
Fork 0
mirror of https://github.com/NGSolve/netgen.git synced 2024-12-25 05:20:34 +05:00
Code Issues Packages Projects Releases Wiki Activity

Merge branch 'visualization_data_geometries' into 'master'

Browse Source 
add functions to collect visualization data to python export of geometries

See merge request jschoeberl/netgen!92
This commit is contained in:
Joachim Schöberl 2018-07-12 17:21:02 +02:00
commit 6034f6ecc3
5 changed files with 270 additions and 0 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

62
libsrc/csg/python_csg.cpp
View File

@ -623,6 +623,68 @@ However, when r = 0, the top part becomes a point(tip) and meshing fails!
}) })
) )
.def_property_readonly ("ntlo", &CSGeometry::GetNTopLevelObjects) .def_property_readonly ("ntlo", &CSGeometry::GetNTopLevelObjects)
.def("_visualizationData", [](shared_ptr<CSGeometry> csg_geo)
{
std::vector<float> vertices;
std::vector<int> trigs;
std::vector<float> normals;
std::vector<float> min = {std::numeric_limits<float>::max(),
std::numeric_limits<float>::max(),
std::numeric_limits<float>::max()};
std::vector<float> max = {std::numeric_limits<float>::lowest(),
std::numeric_limits<float>::lowest(),
std::numeric_limits<float>::lowest()};
std::vector<string> surfnames;
for (int i = 0; i < csg_geo->GetNSurf(); i++)
{
auto surf = csg_geo->GetSurface(i);
surfnames.push_back(surf->GetBCName());
}
csg_geo->FindIdenticSurfaces(1e-6);
csg_geo->CalcTriangleApproximation(0.01,100);
auto nto = csg_geo->GetNTopLevelObjects();
size_t np = 0;
size_t ntrig = 0;
for (int i = 0; i < nto; i++){
np += csg_geo->GetTriApprox(i)->GetNP();
ntrig += csg_geo->GetTriApprox(i)->GetNT();
}
vertices.reserve(np*3);
trigs.reserve(ntrig*4);
normals.reserve(np*3);
int offset_points = 0;
for (int i = 0; i < nto; i++)
{
auto triapprox = csg_geo->GetTriApprox(i);
for (int j = 0; j < triapprox->GetNP(); j++)
for(int k = 0; k < 3; k++) {
float val = triapprox->GetPoint(j)[k];
vertices.push_back(val);
min[k] = min2(min[k], val);
max[k] = max2(max[k],val);
normals.push_back(triapprox->GetNormal(j)[k]);
}
for (int j = 0; j < triapprox->GetNT(); j++)
{
for(int k = 0; k < 3; k++)
trigs.push_back(triapprox->GetTriangle(j)[k]+offset_points);
trigs.push_back(triapprox->GetTriangle(j).SurfaceIndex());
}
offset_points += triapprox->GetNP();
}
py::gil_scoped_acquire ac;
py::dict res;
py::list snames;
for(auto name : surfnames)
snames.append(py::cast(name));
res["vertices"] = MoveToNumpy(vertices);
res["triangles"] = MoveToNumpy(trigs);
res["normals"] = MoveToNumpy(normals);
res["surfnames"] = snames;
res["min"] = MoveToNumpy(min);
res["max"] = MoveToNumpy(max);
return res;
}, py::call_guard<py::gil_scoped_release>())
; ;
m.def("GenerateMesh", FunctionPointer m.def("GenerateMesh", FunctionPointer

10
libsrc/general/ngpython.hpp
View File

@ -14,11 +14,21 @@ namespace pybind11 {
// END EVIL HACK // END EVIL HACK
#include <pybind11/pybind11.h> #include <pybind11/pybind11.h>
#include <pybind11/operators.h> #include <pybind11/operators.h>
#include <pybind11/numpy.h>
namespace py = pybind11; namespace py = pybind11;
#include <iostream> #include <iostream>
#include <sstream> #include <sstream>
template <typename T>
py::array MoveToNumpy(std::vector<T>& vec)
{
auto newvec = new std::vector<T>();
std::swap(*newvec, vec);
auto capsule = py::capsule(newvec, [](void *v) { delete reinterpret_cast<std::vector<T>*>(v); });
return py::array(newvec->size(), newvec->data(), capsule);
}
namespace PYBIND11_NAMESPACE { namespace PYBIND11_NAMESPACE {
template<typename T> template<typename T>
bool CheckCast( py::handle obj ) { bool CheckCast( py::handle obj ) {

73
libsrc/geom2d/python_geom2d.cpp
View File

@ -198,6 +198,79 @@ DLL_HEADER void ExportGeom2d(py::module &m)
return py::tuple(py::make_tuple(xlim, ylim, xpoints, ypoints)); return py::tuple(py::make_tuple(xlim, ylim, xpoints, ypoints));
})) }))
.def("_visualizationData", [](SplineGeometry2d &self)
{
Box<2> box(self.GetBoundingBox());
double xdist = box.PMax()(0) - box.PMin()(0);
double ydist = box.PMax()(1) - box.PMin()(1);
py::dict data;
py::dict segment_data;
auto min_val = py::make_tuple(box.PMin()(0), box.PMin()(1),0);
auto max_val = py::make_tuple(box.PMax()(1),box.PMax()(1),0);
py::list vertices;
py::list domains;
py::list segment_points;
py::list segment_normals;
py::list leftdom;
py::list rightdom;
int max_bcnr = 0;
for(int i = 0; i < self.splines.Size(); i++)
{
std::vector<netgen::GeomPoint<2>> lst;
if (self.splines[i]->GetType().compare("line") == 0)
lst = { self.splines[i]->StartPI(), self.splines[i]->EndPI() };
else if(self.splines[i]->GetType().compare("spline3") == 0)
{
double len = self.splines[i]->Length();
int n = floor(len/(0.05*min(xdist,ydist)));
lst.push_back(self.splines[i]->StartPI());
for (int j = 1; j < n; j++){
lst.push_back(self.splines[i]->GetPoint(j*1./n));
lst.push_back(self.splines[i]->GetPoint(j*1./n));
}
lst.push_back(self.splines[i]->EndPI());
}
else
{
throw NgException("Spline is neither line nor spline3");
}
for (auto point : lst)
{
for(auto val : {point(0), point(1), 0.})
vertices.append(val);
int bcnr = self.GetSpline(i).bc;
max_bcnr = max2(max_bcnr, bcnr);
domains.append(bcnr);
domains.append(self.GetSpline(i).leftdom);
domains.append(self.GetSpline(i).rightdom);
}
// segment data
auto pnt = self.splines[i]->GetPoint(0.5);
segment_points.append(py::make_tuple(pnt(0),pnt(1)));
auto normal = self.GetSpline(i).GetTangent(0.5);
std::swap(normal(0),normal(1));
normal(1) *= -1;
normal *= 1./sqrt(normal(0) * normal(0) + normal(1)*normal(1));
segment_normals.append(py::make_tuple(normal(0),normal(1)));
leftdom.append(self.GetSpline(i).leftdom);
rightdom.append(self.GetSpline(i).rightdom);
}
py::list bcnames;
for (int i = 1; i<max_bcnr + 1; i++)
bcnames.append(self.GetBCName(i));
segment_data["midpoints"] = segment_points;
segment_data["normals"] = segment_normals;
segment_data["leftdom"] = leftdom;
segment_data["rightdom"] = rightdom;
data["segment_data"] = segment_data;
data["vertices"] = vertices;
data["domains"] = domains;
data["min"] = min_val;
data["max"] = max_val;
data["bcnames"] = bcnames;
return data;
})
.def("PointData", FunctionPointer([](SplineGeometry2d &self) .def("PointData", FunctionPointer([](SplineGeometry2d &self)
{ {
py::list xpoints, ypoints, pointindex; py::list xpoints, ypoints, pointindex;

76
libsrc/occ/python_occ.cpp
View File

@ -30,6 +30,82 @@ DLL_HEADER void ExportNgOCC(py::module &m)
self.HealGeometry(); self.HealGeometry();
self.BuildFMap(); self.BuildFMap();
},py::arg("tolerance")=1e-3, py::arg("fixsmalledges")=true, py::arg("fixspotstripfaces")=true, py::arg("sewfaces")=true, py::arg("makesolids")=true, py::arg("splitpartitions")=false,R"raw_string(Heal the OCCGeometry.)raw_string",py::call_guard<py::gil_scoped_release>()) },py::arg("tolerance")=1e-3, py::arg("fixsmalledges")=true, py::arg("fixspotstripfaces")=true, py::arg("sewfaces")=true, py::arg("makesolids")=true, py::arg("splitpartitions")=false,R"raw_string(Heal the OCCGeometry.)raw_string",py::call_guard<py::gil_scoped_release>())
.def("_visualizationData", [] (shared_ptr<OCCGeometry> occ_geo)
{
std::vector<float> vertices;
std::vector<int> trigs;
std::vector<float> normals;
std::vector<float> min = {std::numeric_limits<float>::max(),
std::numeric_limits<float>::max(),
std::numeric_limits<float>::max()};
std::vector<float> max = {std::numeric_limits<float>::lowest(),
std::numeric_limits<float>::lowest(),
std::numeric_limits<float>::lowest()};
std::vector<string> surfnames;
auto box = occ_geo->GetBoundingBox();
for(int i = 0; i < 3; i++)
{
min[i] = box.PMin()[i];
max[i] = box.PMax()[i];
}
occ_geo->BuildVisualizationMesh(0.01);
gp_Pnt2d uv;
gp_Pnt pnt;
gp_Vec n;
gp_Pnt p[3];
int count = 0;
for (int i = 1; i <= occ_geo->fmap.Extent(); i++)
{
surfnames.push_back("occ_surface" + to_string(i));
auto face = TopoDS::Face(occ_geo->fmap(i));
auto surf = BRep_Tool::Surface(face);
TopLoc_Location loc;
BRepAdaptor_Surface sf(face, Standard_False);
BRepLProp_SLProps prop(sf, 1, 1e-5);
Handle(Poly_Triangulation) triangulation = BRep_Tool::Triangulation (face, loc);
if (triangulation.IsNull())
cout << "cannot visualize face " << i << endl;
trigs.reserve(trigs.size() + triangulation->NbTriangles()*4);
vertices.reserve(vertices.size() + triangulation->NbTriangles()*3*3);
normals.reserve(normals.size() + triangulation->NbTriangles()*3*3);
for (int j = 1; j < triangulation->NbTriangles()+1; j++)
{
auto triangle = (triangulation->Triangles())(j);
for (int k = 1; k < 4; k++)
p[k-1] = (triangulation->Nodes())(triangle(k)).Transformed(loc);
for (int k = 1; k < 4; k++)
{
vertices.insert(vertices.end(),{float(p[k-1].X()), float(p[k-1].Y()), float(p[k-1].Z())});
trigs.insert(trigs.end(),{count, count+1, count+2,i});
count += 3;
uv = (triangulation->UVNodes())(triangle(k));
prop.SetParameters(uv.X(), uv.Y());
if (prop.IsNormalDefined())
n = prop.Normal();
else
{
gp_Vec a(p[0], p[1]);
gp_Vec b(p[0], p[2]);
n = b^a;
}
if (face.Orientation() == TopAbs_REVERSED) n*= -1;
normals.insert(normals.end(),{float(n.X()), float(n.Y()), float(n.Z())});
}
}
}
py::gil_scoped_acquire ac;
py::dict res;
py::list snames;
for(auto name : surfnames)
snames.append(py::cast(name));
res["vertices"] = MoveToNumpy(vertices);
res["triangles"] = MoveToNumpy(trigs);
res["normals"] = MoveToNumpy(normals);
res["surfnames"] = snames;
res["min"] = MoveToNumpy(min);
res["max"] = MoveToNumpy(max);
return res;
}, py::call_guard<py::gil_scoped_release>())
; ;
m.def("LoadOCCGeometry",FunctionPointer([] (const string & filename) m.def("LoadOCCGeometry",FunctionPointer([] (const string & filename)
{ {

49
libsrc/stlgeom/python_stl.cpp
View File

@ -20,6 +20,55 @@ DLL_HEADER void ExportSTL(py::module & m)
{ {
py::class_<STLGeometry,shared_ptr<STLGeometry>, NetgenGeometry> (m,"STLGeometry") py::class_<STLGeometry,shared_ptr<STLGeometry>, NetgenGeometry> (m,"STLGeometry")
.def(py::init<>()) .def(py::init<>())
.def("_visualizationData", [](shared_ptr<STLGeometry> stl_geo)
{
std::vector<float> vertices;
std::vector<int> trigs;
std::vector<float> normals;
std::vector<float> min = {std::numeric_limits<float>::max(),
std::numeric_limits<float>::max(),
std::numeric_limits<float>::max()};
std::vector<float> max = {std::numeric_limits<float>::lowest(),
std::numeric_limits<float>::lowest(),
std::numeric_limits<float>::lowest()};
std::vector<string> surfnames;
surfnames.push_back("stl");
vertices.reserve(stl_geo->GetNT()*3*3);
trigs.reserve(stl_geo->GetNT()*4);
normals.reserve(stl_geo->GetNT()*3*3);
size_t ii = 0;
for(int i = 0; i < stl_geo->GetNT(); i++)
{
auto& trig = stl_geo->GetTriangle(i+1);
for(int k = 0; k < 3; k++)
{
trigs.push_back(ii++);
auto& pnt = stl_geo->GetPoint(trig[k]);
for (int l = 0; l < 3; l++)
{
float val = pnt[l];
vertices.push_back(val);
min[l] = min2(min[l], val);
max[l] = max2(max[l], val);
normals.push_back(trig.Normal()[l]);
}
}
trigs.push_back(0);
}
py::gil_scoped_acquire ac;
py::dict res;
py::list snames;
for(auto name : surfnames)
snames.append(py::cast(name));
res["vertices"] = MoveToNumpy(vertices);
res["triangles"] = MoveToNumpy(trigs);
res["normals"] = MoveToNumpy(normals);
res["surfnames"] = snames;
res["min"] = MoveToNumpy(min);
res["max"] = MoveToNumpy(max);
return res;
}, py::call_guard<py::gil_scoped_release>())
; ;
m.def("LoadSTLGeometry", FunctionPointer([] (const string & filename) m.def("LoadSTLGeometry", FunctionPointer([] (const string & filename)
{ {