Sun, 09/04/2011 - 09:41
Forums:
Hello,
For those interested in doing Gauss curvature here is a Gaussian Curvature AIS object.
Attached is the Header,implementation.
Attachments:
Collaborative development portal
Sun, 09/04/2011 - 09:43
/// implementation file and an image of how it should look.
// AIS_Gauss.cpp: implementation of the AIS_Gauss class.
//
//////////////////////////////////////////////////////////////////////
//// following the code on this link:
//http://www.opencascade.org/org/forum/thread_1125/
#include "AIS_Gauss.hxx"
#ifdef _DEBUG
#undef THIS_FILE
static char THIS_FILE[]=__FILE__;
//#define new DEBUG_NEW
#endif
IMPLEMENT_STANDARD_HANDLE(AIS_Gauss,AIS_InteractiveObject)
IMPLEMENT_STANDARD_RTTI(AIS_Gauss,AIS_InteractiveObject)
//
// Foreach ancestors, we add a IMPLEMENT_STANDARD_SUPERTYPE and
// a IMPLEMENT_STANDARD_SUPERTYPE_ARRAY_ENTRY macro.
// We must respect the order: from the direct ancestor class
// to the base class.
//
IMPLEMENT_STANDARD_TYPE(AIS_Gauss)
IMPLEMENT_STANDARD_SUPERTYPE(AIS_InteractiveObject)
IMPLEMENT_STANDARD_SUPERTYPE(SelectMgr_SelectableObject)
IMPLEMENT_STANDARD_SUPERTYPE(PrsMgr_PresentableObject)
IMPLEMENT_STANDARD_SUPERTYPE(MMgt_TShared)
IMPLEMENT_STANDARD_SUPERTYPE(Standard_Transient)
IMPLEMENT_STANDARD_SUPERTYPE_ARRAY()
IMPLEMENT_STANDARD_SUPERTYPE_ARRAY_ENTRY(AIS_InteractiveObject)
IMPLEMENT_STANDARD_SUPERTYPE_ARRAY_ENTRY(SelectMgr_SelectableObject)
IMPLEMENT_STANDARD_SUPERTYPE_ARRAY_ENTRY(PrsMgr_PresentableObject)
IMPLEMENT_STANDARD_SUPERTYPE_ARRAY_ENTRY(MMgt_TShared)
IMPLEMENT_STANDARD_SUPERTYPE_ARRAY_ENTRY(Standard_Transient)
IMPLEMENT_STANDARD_SUPERTYPE_ARRAY_END()
IMPLEMENT_STANDARD_TYPE_END(AIS_Gauss)
#include
#include
#include
#include
#include
#include
#include "PrsMgr_PresentationManager2d.hxx"
#include "SelectMgr_Selection.hxx"
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
#include
//////////////////////////////////////////////////////////////////////
// Construction/Destruction
//////////////////////////////////////////////////////////////////////
//AIS_Gauss::AIS_Gauss(const Graphic3d_Array1OfVertex &points)
//:cloud_data(points)
//{
//
//}
AIS_Gauss::AIS_Gauss(TopoDS_Shape theshape, bool showedges)
{
myShape = theshape;
_showedges = showedges;
}
AIS_Gauss::~AIS_Gauss()
{
}
void AIS_Gauss::Compute(const Handle(PrsMgr_PresentationManager3d)& aPresentationManager,
const Handle(Prs3d_Presentation)& aPresentation,
const Standard_Integer aMode)
{
int gradationsteps = 20;
int numberColors = gradationsteps;
int redsteps = gradationsteps/3;
int bluesteps = gradationsteps/3;
int greensteps = gradationsteps -(redsteps + bluesteps);
float * curvatureColors = new float[(gradationsteps-1)*3]; //blue to red or whatever
int redcounter = 0;
int greencounter =0;
int bluecounter = 0;
int colorcount = 0;
for(int cs=0 ; cs < gradationsteps ; cs++)
{
double myR ;
double myG ;
double myB ;
if (cs <= (redsteps -1) )
{
myR = 1 - ((redcounter-1) / (double) redsteps);
myG = 0;
myB = 0;
if (myR = 0.1) myG = 1.0;
redcounter ++;
} else if (cs > (redsteps -1) && cs < (redsteps + greensteps -1) )
{
myR = 0;
myG = 1-((greencounter-1) / (double) greensteps);
myB = 0;
if (myG = 0.1) myB = 1.0;
greencounter++;
} else if (cs > (redsteps + greensteps -1) )
{
myR = 0;
myG = 0;
myB = 1-((bluecounter-1) / (double)bluesteps);
bluecounter++;
}
int rindex = colorcount ;
int gindex = colorcount+1;
int bindex = colorcount+2 ;
curvatureColors[rindex]= myR;
curvatureColors[gindex]= myG;
curvatureColors[bindex]= myB;
colorcount+=3;
}
Handle(Graphic3d_Structure) theStructure = Handle(Graphic3d_Structure)::DownCast(aPresentation);
Handle(Graphic3d_Group) mygroup= new Graphic3d_Group(theStructure);
bool visualizationDone;
bool flatShading = false;
//Explore the shape and get all the triangles of the corresponding mesh
TopExp_Explorer Ex;
for (Ex.Init(myShape,TopAbs_FACE); Ex.More(); Ex.Next())
{
TopoDS_Shape curshape = Ex.Current();
Handle(Geom_Surface) mysurf = BRep_Tool::Surface(TopoDS::Face(curshape));
TopoDS_Face myFace = TopoDS::Face(curshape);
TopLoc_Location L;
Handle(Poly_Triangulation) myT = BRep_Tool::Triangulation(myFace, L);
int nbnodes = myT->NbNodes();
int nbtriangles = myT->NbTriangles();
//Poly_Connect pc(myT);
const TColgp_Array1OfPnt& Nodes = myT->Nodes();
const TColgp_Array1OfPnt2d& UVNodes = myT->UVNodes();
const Poly_Array1OfTriangle& triangles = myT->Triangles();
Handle(Graphic3d_ArrayOfTriangles) trianglearr = new Graphic3d_ArrayOfTriangles(triangles.Length()* 3,0,false,true);
Graphic3d_Array1OfVertexNC Points(1,triangles.Length()* 3);
Handle(Graphic3d_ArrayOfPolylines) polylines = new Graphic3d_ArrayOfPolylines ( triangles.Length()* 6,triangles.Length()* 6 );
Standard_Real fUMin, fUMax, fVMin, fVMax;
BRepTools::UVBounds(TopoDS::Face(myFace), fUMin, fUMax, fVMin, fVMax);
Standard_Real fU = (fUMax-fUMin)/2.0;
Standard_Real fV = (fVMax-fVMin)/2.0;
GeomLProp_SLProps prop(mysurf, fU, fV, 1, Precision::Confusion());
double maxVisualizedCurvature = 0;
double minVisualizedCurvature = 0;
for (int i = 1; i < UVNodes.Length()+1; i++)
{
gp_Pnt2d uvval = UVNodes.Value(i);
prop.SetParameters(uvval.X(),uvval.Y());
double g = prop.MeanCurvature();//GaussianCurvature();
if (g < minVisualizedCurvature)minVisualizedCurvature = g;
if (g > maxVisualizedCurvature)maxVisualizedCurvature = g;
}
//you could use something like this:
double span = 1;
span = (maxVisualizedCurvature - minVisualizedCurvature) / numberColors;
int spanNumber = 0;
double curvatureStrength = 0;
//loop over all faces & triangles
double prevR;
double prevG;
double prevB;
double R;
double G;
double B;
//For each triangle get the three nodes
//Graphic3d_Array1OfVertexNC Points(1,3);
if (myT->HasUVNodes())
{
for (int i = 1; i < triangles.Length()+1; i++)
{
Poly_Triangle mytriangle = triangles.Value(i);
int n1,n2,n3;
mytriangle.Get(n1,n2,n3);
gp_Pnt p1 = Nodes.Value(n1);
gp_Pnt p2 = Nodes.Value(n2);
gp_Pnt p3 = Nodes.Value(n3);
if (_showedges)
{
// edge1
polylines->AddBound ( 2 );
polylines->AddVertex ( p1.X(), p1.Y(), p1.Z() );
polylines->AddVertex ( p2.X(), p2.Y(), p2.Z() );
// edge1
polylines->AddBound ( 2 );
polylines->AddVertex ( p2.X(), p2.Y(), p2.Z() );
polylines->AddVertex ( p3.X(), p3.Y(), p3.Z() );
// edge1
polylines->AddBound ( 2 );
polylines->AddVertex ( p3.X(), p3.Y(), p3.Z() );
polylines->AddVertex ( p1.X(), p1.Y(), p1.Z() );
}
//trianglearr->AddBound(3);
for(int k = 1; k <= 3 ; k++)
{
int vertindex = mytriangle.Value(k);
if(k == 1)
prop.SetParameters(UVNodes(n1).X(), UVNodes(n1).Y());
else if(k == 2)
prop.SetParameters(UVNodes(n2).X(), UVNodes(n2).Y());
else
prop.SetParameters(UVNodes(n3).X(), UVNodes(n3).Y());
if (prop.IsCurvatureDefined())
{
// switch curvature types later like the mean curvature or the max curvature.
double curvature = prop.MeanCurvature();//GaussianCurvature();
visualizationDone = true;
if(span != 0)
{
spanNumber = 1 + (int) ( (curvature - minVisualizedCurvature) /span);
curvatureStrength = (curvature - minVisualizedCurvature - (spanNumber - 1) * span) / span;
}
else
{
curvatureStrength = 1.0;
spanNumber = 0;
}
if( (spanNumber < 0) || (spanNumber > numberColors - 1) ||
(curvatureStrength > 1.0) || (curvatureStrength < 0) )
{
curvatureStrength = 1.0;
if(spanNumber > numberColors - 1)
spanNumber = numberColors - 1;
else if(spanNumber < 0)
spanNumber = 0;
}
if(spanNumber > 0 )
{
prevR = curvatureColors[(spanNumber - 1) + 0];
prevG = curvatureColors[(spanNumber - 1) + 1];
prevB = curvatureColors[(spanNumber - 1) + 2];
}
else
{
prevR = curvatureColors[0];
prevG = curvatureColors[1];
prevB = curvatureColors[2];
}
R = prevR + (float)curvatureStrength * ( curvatureColors[spanNumber + 0] - prevR );
G = prevG + (float)curvatureStrength * ( curvatureColors[spanNumber + 1] - prevG );
B = prevB + (float)curvatureStrength * ( curvatureColors[spanNumber + 2] - prevB );
if(flatShading == true)
{
R = ceil(R);
G = ceil(G);
B = ceil(B);
}
gp_Pnt mynodept = Nodes.Value(vertindex);
if (R > 1) R = 1;
if (G > 1) G = 1;
if (B > 1) B = 1;
Quantity_Color thecolor(R,G,B,Quantity_TOC_RGB);
//qDebug() << "color and gauss:" << R << "," << G << "," << B << "-" << curvature;
trianglearr->AddVertex(mynodept,thecolor);
//points(k).SetColor(Quantity_Color(R,G,B,Quantity_TOC_RGB));
// mygroup->Marker(Graphic3d_Vertex(mynodept.X(),mynodept.Y(),mynodept.Z()));
}
}
} // end of lopping thorugh points
} // end of has uv nodes
// draw mesh triangles
Handle(Graphic3d_AspectFillArea3d) aspect = new Graphic3d_AspectFillArea3d();
aspect->SetInteriorStyle(Aspect_IS_SOLID);
//mygroup->Clear();
mygroup->BeginPrimitives();
mygroup->SetPrimitivesAspect(aspect);
mygroup->AddPrimitiveArray(trianglearr);
mygroup->EndPrimitives();
if(_showedges)
{
// draw mesh edges inside polylines array
mygroup->SetPrimitivesAspect ( new Graphic3d_AspectLine3d( Quantity_NameOfColor::Quantity_NOC_BLACK , Aspect_TOL_SOLID, 0.1 ) );
mygroup->BeginPrimitives ();
mygroup->AddPrimitiveArray ( polylines );
mygroup->EndPrimitives ();
}
} // end of iterator on face
}
void AIS_Gauss::Compute(const Handle(Prs3d_Projector)& aProjector,
const Handle(Prs3d_Presentation)& aPresentation)
{
}
void AIS_Gauss::Compute(const Handle(PrsMgr_PresentationManager2d)& aPresentationManager,
const Handle(Graphic2d_GraphicObject)& aGrObj,
const Standard_Integer unMode)
{
}
void AIS_Gauss::ComputeSelection(const Handle(SelectMgr_Selection)& aSelection,
const Standard_Integer unMode)
{
}
Sun, 09/04/2011 - 09:45
// Usage:
// ic == interactive context.
// first parmeter is the TopoDS_Shape
// second Parameter is True or False for Displaying edges.
TopoDS_Shape myshape = makeshape;
Handle(AIS_Gauss) mygauss = new AIS_Gauss(myshape,true);
ic->Display(mygauss);
// best
// AlexP
Sun, 09/04/2011 - 09:50
by the way remove the following includes at the top of the implementation file (*.cpp) so its not Qt dependant.
I was using qt gradients to blend between one color and the next, not needed any more
my color scales are still a mess, if anyone knows a correct way of generating the color scales from red to blue, in a way you can specify how many gradation steps, please tag along.
#include
#include
#include
#include
Best,
AlexP
Sun, 09/04/2011 - 11:01
cool! thanks for sharing...
Sun, 09/04/2011 - 12:18
Thanks Jelle,
I actually used your version as guide during the process, im not a python user, so i was super confused with the tuple thing. how functions can return more than one value, but i got the overall principle pretty much from your python implementation.
i noticed you where using triangle fans made of 4 points, so i think you where creating the points on surface looping through the uv points.
i made it using the brep meshing.
in doing so, i noticed a weird phenomena, maybe its part of the design intent of the algorithm.
the size of the surface, changes the effect of the accuracy parameter fed to the brepmesh algorithm,
it seems the algorithm is creating the decimation based on euclidean distance, rather than using the ratio of the iso parm length.
i think on my brepmesh algorithm i will add a relation:
steps = 100;
accuracy = area of surface / steps;
note to everyone-else :
before i forget:
you need to mesh before you use the AIS_Gauss code so add the following just after creation of the TopoDS_Shape.
// create the shape ...
TopoDS_Shape myshape = makeshape;
// determine accuracy of meshing ... the more dense the mesh, the more accurate the gaussian map is.
accuracy = 10; // this value depends on what you need, ive been getting good results between 10 - 0.2
// if you want to make it constant use steps = 100; accuracy = measure area of myshape / steps;
BRepTools::Clean(myshape);
BRepMesh::Mesh(myshape, accuracy);
// then you visualize
Handle(AIS_Gauss) mygauss = new AIS_Gauss(myshape,true);
ic->Display(mygauss);
Best,
AlexP
Thu, 09/08/2011 - 23:27
The following code fails to display the AIS_Gauss shape:
gp_Pnt P1(0,0,0);
gp_Pnt P2(-100,0,0);
gp_Pnt P3(0,0,100);
gp_Pnt P4(-100,0,100);
TColgp_Array1OfPnt thePoles(1,4);
thePoles.SetValue(1, P1);
thePoles.SetValue(2, P2);
thePoles.SetValue(3, P3);
thePoles.SetValue(4, P4);
Handle(Geom_BezierCurve) theBez = new Geom_BezierCurve(thePoles);
Handle(Geom_SurfaceOfRevolution) theRev = new Geom_SurfaceOfRevolution(theBez, gp::OZ());
BRepBuilderAPI_MakeEdge MkEdge(theBez);
BRepBuilderAPI_MakeWire MkWire(MkEdge.Edge());
BRepBuilderAPI_MakeFace MkFace(MkWire.Wire());
BRepPrimAPI_MakeRevol MkRevol(MkWire.Wire(), gp::OZ());
myAISContext->Display(new AIS_Gauss(MkRevol.Shape(), false));
The failure is in the Compute method line no 204 where it tries to retrieve the nodes after triangulation:
Handle(Poly_Triangulation) myT = BRep_Tool::Triangulation(myFace, L);
int nbnodes = myT->NbNodes();
Am I missing something that needs to be done?
Fri, 09/09/2011 - 08:03
you need to mesh the shape prior to feeding it to the AIS_Gauss
TopoDS_Shape myrev = MkRevol.Shape();
BRepTools::Clean(myrev);
BRepMesh::Mesh(myshape, 5);
myAISContext->Display(new AIS_Gauss(MkRevol.Shape(), true));
Best,
AlexP
Fri, 09/09/2011 - 12:06
Correction to code above **************
TopoDS_Shape myrev = MkRevol.Shape();
BRepTools::Clean(myrev);
BRepMesh::Mesh(myrev, 5);
myAISContext->Display(new AIS_Gauss(myrev, true));
Best,
AlexP
Fri, 09/09/2011 - 12:12
By the way the color mapping of this file is a little goofy, i need to understand how to blend from one color to another.
the issue is Gaussian curvature ranges from negative to positive
the negative limit is red
the neutral or (flat) condition is orange
the positive limit is blue
you need to create a color scale that takes the gaussian values and divis-them-up in to as many color steps you need.
ill look into it later and will add it to the openshapefactory AIS_Gauss.cpp
Fri, 09/09/2011 - 21:00
fyi, a loooong time ago, I adapted the gradient class from this article ( http://www.codeproject.com/KB/miscctrl/gradient.aspx ) to be pure C++ without any Win32 related class to use in OpenGL.
(I cannot share the code but you can do the same, it's not that hard )
Stephane
Sat, 09/10/2011 - 08:23
Hi Stephane,
Great Link thanks for the tip.
Best,
AlexP