Hi J-S,
>This is the part I'm having trouble with. I actually don't mind using
>the scene bound I get with the ComputeBoundingBoxVisitor, but getting a
>usable camera frustum (both for main and light cameras) and then
>intersecting those is the hard part for me. Any tips there? How do you
>create the convex objects? How do you do the intersection itself?
This piece of code was actually written by Robert ;-). I use his utility
CustomPolytope class from osgSim::OverlayNode.cpp. It computes the
intersection of two polytopes. Its named CustomPolytope but it does work
well as convex hull representation. I attached this class repacked into
separate header. Its identical to CustomPolytope from OverlayNode with one
bug fixed in my version. I submitted this fix to OverlayNode.cpp but its
still pending. So I suggest to use my header before Robert returns.
Code for computing frustum polytope is obvious provided we have camera view
& projection matrices:
/////////////////////////////////
osg::Matrix mvp = camera->getViewMatrix() * camera->getProjectionMatrix();
CustomPolytope cpCameraFrustum;
cpCameraFrustum.setToUnitFrustum( );
// transform frustum from clip space into world coords
cpCameraFrustum.transform( osg::Matrix::inverse( mvp ), mvp );
/////////////////////////////////
Of course there is also a method to build CustomPolytope from bounding box:
/////////////////////////////////
CustomPolytope cpSceneBounds;
cp.setToBoundingBox( sceneBounds );
/////////////////////////////////
Computing intersection is straight forward
Suppose we have two polytopes for main camera frustum & coarse shadow camera
frustum and scene bounds polytope.
/////////////////////////////////
CustomPolytope cpMainCamera, cpCoarseShadowCamera, cpSceneBounds;
CustomPolytope cpIntersection( cpSceneBounds );
cpIntersection.cut( cpMainCamera );
cpIntersection.cut( cpCoarseShadowCamera );
/////////////////////////////////
CustomPolytope::getVertices returns corners of the convex hull so its easy
to find its bounding box. I compute precise shadow camera projection using
such bounding box around intersection polytope transformed to shadow space.
And thats it. Thank you Robert !
Wojtek Lewandowski
#ifndef CUSTOMPOLYTOPE_H
#define CUSTOMPOLYTOPE_H
#include<osg/Geometry>
#include<osg/Polytope>
#include<osg/CoordinateSystemNode>
// Class directly copied osgSim::OverlayNode (+prefixed few Vec3d occurences
with osg::)
class CustomPolytope
{
public:
CustomPolytope() {}
typedef std::vector<osg::Vec3d> Vertices;
struct Face
{
std::string name;
osg::Plane plane;
Vertices vertices;
};
Face& createFace() { _faces.push_back(Face()); return _faces.back(); }
/** Create a Polytope which is a cube, centered at 0,0,0, with sides of 2
units.*/
void setToUnitFrustum(bool withNear=true, bool withFar=true)
{
const osg::Vec3d v000(-1.0,-1.0,-1.0);
const osg::Vec3d v010(-1.0,1.0,-1.0);
const osg::Vec3d v001(-1.0,-1.0,1.0);
const osg::Vec3d v011(-1.0,1.0,1.0);
const osg::Vec3d v100(1.0,-1.0,-1.0);
const osg::Vec3d v110(1.0,1.0,-1.0);
const osg::Vec3d v101(1.0,-1.0,1.0);
const osg::Vec3d v111(1.0,1.0,1.0);
_faces.clear();
{ // left plane.
Face& face = createFace();
face.name = "left";
face.plane.set(1.0,0.0,0.0,1.0);
face.vertices.push_back(v000);
face.vertices.push_back(v001);
face.vertices.push_back(v011);
face.vertices.push_back(v010);
}
{ // right plane.
Face& face = createFace();
face.name = "right";
face.plane.set(-1.0,0.0,0.0,1.0);
face.vertices.push_back(v100);
face.vertices.push_back(v110);
face.vertices.push_back(v111);
face.vertices.push_back(v101);
}
{ // bottom plane.
Face& face = createFace();
face.name = "bottom";
face.plane.set(0.0,1.0,0.0,1.0);
face.vertices.push_back(v000);
face.vertices.push_back(v100);
face.vertices.push_back(v101);
face.vertices.push_back(v001);
}
{ // top plane.
Face& face = createFace();
face.name = "top";
face.plane.set(0.0,-1.0,0.0,1.0);
face.vertices.push_back(v111);
face.vertices.push_back(v011);
face.vertices.push_back(v010);
face.vertices.push_back(v110);
}
if (withNear)
{ // near plane
Face& face = createFace();
face.name = "near";
face.plane.set(0.0,0.0,1.0,1.0);
face.vertices.push_back(v000);
face.vertices.push_back(v010);
face.vertices.push_back(v110);
face.vertices.push_back(v100);
}
if (withFar)
{ // far plane
Face& face = createFace();
face.name = "far";
face.plane.set(0.0,0.0,-1.0,1.0);
face.vertices.push_back(v001);
face.vertices.push_back(v101);
face.vertices.push_back(v111);
face.vertices.push_back(v011);
}
}
void setToBoundingBox(const osg::BoundingBox& bb)
{
#if 0
osg::notify(osg::NOTICE)<<"setToBoundingBox xrange "<<bb.xMin()<<"
"<<bb.xMax()<<std::endl;
osg::notify(osg::NOTICE)<<" "<<bb.yMin()<<"
"<<bb.yMax()<<std::endl;
osg::notify(osg::NOTICE)<<" "<<bb.zMin()<<"
"<<bb.zMax()<<std::endl;
#endif
const osg::Vec3d v000(bb.xMin(),bb.yMin(), bb.zMin());
const osg::Vec3d v010(bb.xMin(),bb.yMax(), bb.zMin());
const osg::Vec3d v001(bb.xMin(),bb.yMin(), bb.zMax());
const osg::Vec3d v011(bb.xMin(),bb.yMax(), bb.zMax());
const osg::Vec3d v100(bb.xMax(),bb.yMin(), bb.zMin());
const osg::Vec3d v110(bb.xMax(),bb.yMax(), bb.zMin());
const osg::Vec3d v101(bb.xMax(),bb.yMin(), bb.zMax());
const osg::Vec3d v111(bb.xMax(),bb.yMax(), bb.zMax());
_faces.clear();
{ // x min plane
Face& face = createFace();
face.name = "xMin";
face.plane.set(1.0,0.0,0.0,-bb.xMin());
face.vertices.push_back(v000);
face.vertices.push_back(v001);
face.vertices.push_back(v011);
face.vertices.push_back(v010);
}
{ // x max plane.
Face& face = createFace();
face.name = "xMax";
face.plane.set(-1.0,0.0,0.0,bb.xMax());
face.vertices.push_back(v100);
face.vertices.push_back(v110);
face.vertices.push_back(v111);
face.vertices.push_back(v101);
}
{ // y min plane.
Face& face = createFace();
face.name = "yMin";
face.plane.set(0.0,1.0,0.0,-bb.yMin());
face.vertices.push_back(v000);
face.vertices.push_back(v100);
face.vertices.push_back(v101);
face.vertices.push_back(v001);
}
{ // y max plane.
Face& face = createFace();
face.name = "yMax";
face.plane.set(0.0,-1.0,0.0,bb.yMax());
face.vertices.push_back(v111);
face.vertices.push_back(v011);
face.vertices.push_back(v010);
face.vertices.push_back(v110);
}
{ // z min plane
Face& face = createFace();
face.name = "zMin";
face.plane.set(0.0,0.0,1.0,-bb.zMin());
face.vertices.push_back(v000);
face.vertices.push_back(v010);
face.vertices.push_back(v110);
face.vertices.push_back(v100);
}
{ // z max plane
Face& face = createFace();
face.name = "zMax";
face.plane.set(0.0,0.0,-1.0,bb.zMax());
face.vertices.push_back(v001);
face.vertices.push_back(v101);
face.vertices.push_back(v111);
face.vertices.push_back(v011);
}
}
void transform(const osg::Matrix& matrix, const osg::Matrix& inverse)
{
for(Faces::iterator itr = _faces.begin();
itr != _faces.end();
++itr)
{
Face& face = *itr;
face.plane.transformProvidingInverse(inverse);
for(Vertices::iterator vitr = face.vertices.begin();
vitr != face.vertices.end();
++vitr)
{
*vitr = *vitr * matrix;
}
}
}
void insertVertex(const osg::Vec3d& vertex, osg::EllipsoidModel* em=0,
double minHeight=0.0)
{
// osg::notify(osg::NOTICE)<<"Inserting vertex "<<vertex<<std::endl;
Faces removedFaces;
Faces::iterator itr;
for(itr = _faces.begin();
itr != _faces.end();
)
{
Face& face = *itr;
if (face.plane.distance(vertex)<0.0)
{
removedFaces.push_back(face);
itr = _faces.erase(itr);
}
else
{
++itr;
}
}
if (removedFaces.empty()) return;
typedef std::pair<osg::Vec3d, osg::Vec3d> Edge;
typedef std::map<Edge,int> Edges;
Edges edges;
double numVerticesAdded=0.0;
osg::Vec3d center;
for(itr = removedFaces.begin();
itr != removedFaces.end();
++itr)
{
Face& face = *itr;
Vertices& vertices = face.vertices;
for(unsigned int i=0; i<vertices.size(); ++i)
{
osg::Vec3d& a = vertices[i];
osg::Vec3d& b = vertices[ (i+1) % vertices.size()];
if (a<b) ++edges[Edge(a,b)];
else ++edges[Edge(b,a)];
center += a;
numVerticesAdded += 1.0;
}
}
if (numVerticesAdded==0.0) return;
center /= numVerticesAdded;
typedef std::set<osg::Vec3> VertexSet;
VertexSet uniqueVertices;
for(Edges::iterator eitr = edges.begin();
eitr != edges.end();
++eitr)
{
if (eitr->second==1)
{
// osg::notify(osg::NOTICE)<<" edge Ok"<<std::endl;
const Edge& edge = eitr->first;
Face face;
face.name = "baseSide";
face.plane.set(vertex, edge.first, edge.second);
face.vertices.push_back(vertex);
face.vertices.push_back(edge.first);
face.vertices.push_back(edge.second);
if (face.plane.distance(center)<0.0) face.plane.flip();
_faces.push_back(face);
uniqueVertices.insert(edge.first);
uniqueVertices.insert(edge.second);
}
}
// now trim the new polytope back the desired height
if (em)
{
// compute the base vertices at the new height
Vertices baseVertices;
for(VertexSet::iterator itr = uniqueVertices.begin();
itr != uniqueVertices.end();
++itr)
{
const osg::Vec3d& point = *itr;
double latitude, longitude, height;
em->convertXYZToLatLongHeight(point.x(), point.y(), point.z(),
latitude, longitude, height);
osg::Vec3d normal(point);
normal.normalize();
baseVertices.push_back(point - normal * (height - minHeight));
}
//compute centroid of the base vertices
osg::Vec3d center;
double totalArea = 0;
for(unsigned int i=0; i<baseVertices.size()-1; ++i)
{
const osg::Vec3d& a = baseVertices[i];
const osg::Vec3d& b = baseVertices[i+1];
const osg::Vec3d& c = baseVertices[(i+2)%baseVertices.size()];
double area = ((a-b)^(b-c)).length()*0.5;
osg::Vec3d localCenter = (a+b+c)/3.0;
center += localCenter*area;
totalArea += area;
}
center /= totalArea;
osg::Vec3d normal(center);
normal.normalize();
osg::Plane basePlane(normal, center);
cut(basePlane,"basePlane");
}
// osg::notify(osg::NOTICE)<<" Removed faces
"<<removedFaces.size()<<std::endl;
}
void projectDowntoBase(const osg::Vec3d& control, const osg::Vec3d& normal)
{
// osg::notify(osg::NOTICE)<<"CustomPolytope::projectDowntoBase not
implementated yet."<<std::endl;
Faces removedFaces;
Faces::iterator itr;
for(itr = _faces.begin();
itr != _faces.end();
)
{
Face& face = *itr;
if ((face.plane.getNormal()*normal)>=0.0)
{
removedFaces.push_back(face);
itr = _faces.erase(itr);
}
else
{
++itr;
}
}
if (removedFaces.empty()) return;
typedef std::pair<osg::Vec3d, osg::Vec3d> Edge;
typedef std::map<Edge,int> Edges;
Edges edges;
double numVerticesAdded=0.0;
osg::Vec3d center;
for(itr = removedFaces.begin();
itr != removedFaces.end();
++itr)
{
Face& face = *itr;
Vertices& vertices = face.vertices;
for(unsigned int i=0; i<vertices.size(); ++i)
{
osg::Vec3d& a = vertices[i];
osg::Vec3d& b = vertices[ (i+1) % vertices.size()];
if (a<b) ++edges[Edge(a,b)];
else ++edges[Edge(b,a)];
center += a;
numVerticesAdded += 1.0;
}
}
if (numVerticesAdded==0.0) return;
center /= numVerticesAdded;
typedef std::set<osg::Vec3> VertexSet;
VertexSet uniqueVertices;
for(Edges::iterator eitr = edges.begin();
eitr != edges.end();
++eitr)
{
if (eitr->second==1)
{
// osg::notify(osg::NOTICE)<<" edge Ok"<<std::endl;
const Edge& edge = eitr->first;
double h_first = (edge.first-control) * normal;
osg::Vec3d projected_first = edge.first - normal * h_first;
double h_second = (edge.second-control) * normal;
osg::Vec3d projected_second = edge.second - normal * h_second;
Face face;
face.name = "baseSide";
face.plane.set(projected_first, edge.first, edge.second);
face.vertices.push_back(projected_first);
face.vertices.push_back(projected_second);
face.vertices.push_back(edge.second);
face.vertices.push_back(edge.first);
if (face.plane.distance(center)<0.0) face.plane.flip();
_faces.push_back(face);
uniqueVertices.insert(projected_first);
uniqueVertices.insert(projected_second);
}
}
Face newFace;
newFace.name = "basePlane";
newFace.plane.set(normal,control);
osg::Vec3d side = ( fabs(normal.x()) < fabs(normal.y()) ) ?
osg::Vec3(1.0, 0.0, 0.0) :
osg::Vec3(0.0, 1.0, 0.0);
osg::Vec3 v = normal ^ side;
v.normalize();
osg::Vec3 u = v ^ normal;
u.normalize();
typedef std::map<double, osg::Vec3d> AnglePositions;
AnglePositions anglePositions;
for(VertexSet::iterator vitr = uniqueVertices.begin();
vitr != uniqueVertices.end();
++vitr)
{
osg::Vec3d delta = *vitr - center;
double angle = atan2(delta * u, delta * v);
if (angle<0.0) angle += 2.0*osg::PI;
anglePositions[angle] = *vitr;
}
for(AnglePositions::iterator aitr = anglePositions.begin();
aitr != anglePositions.end();
++aitr)
{
newFace.vertices.push_back(aitr->second);
}
_faces.push_back(newFace);
}
void computeSilhoette(const osg::Vec3d& normal, Vertices& vertices)
{
typedef std::pair<osg::Vec3d, osg::Vec3d> EdgePair;
typedef std::vector<Face*> EdgeFaces;
typedef std::map<EdgePair, EdgeFaces> EdgeMap;
EdgeMap edgeMap;
for(Faces::iterator itr = _faces.begin();
itr != _faces.end();
++itr)
{
Face& face = *itr;
for(unsigned int i=0; i<face.vertices.size(); ++i)
{
osg::Vec3d& va = face.vertices[i];
osg::Vec3d& vb = face.vertices[(i+1)%face.vertices.size()];
if (va < vb) edgeMap[EdgePair(va,vb)].push_back(&face);
else edgeMap[EdgePair(vb,va)].push_back(&face);
}
}
typedef std::set<osg::Vec3> VertexSet;
VertexSet uniqueVertices;
for(EdgeMap::iterator eitr = edgeMap.begin();
eitr != edgeMap.end();
++eitr)
{
const EdgePair& edge = eitr->first;
const EdgeFaces& edgeFaces = eitr->second;
if (edgeFaces.size()==1)
{
// osg::notify(osg::NOTICE)<<"Open edge found
"<<edgeFaces.front()->name<<std::endl;
}
else if (edgeFaces.size()==2)
{
double dotProduct0 = edgeFaces[0]->plane.getNormal() * normal;
double dotProduct1 = edgeFaces[1]->plane.getNormal() * normal;
if (dotProduct0 * dotProduct1 <0.0)
{
// osg::notify(osg::NOTICE)<<" Silhoette edge found
"<<edgeFaces[0]->name<<" "<<edgeFaces[1]->name<<std::endl;
uniqueVertices.insert(edge.first);
uniqueVertices.insert(edge.second);
}
else
{
// osg::notify(osg::NOTICE)<<" Non silhoette edge found
"<<edgeFaces[0]->name<<" "<<edgeFaces[1]->name<<std::endl;
}
}
else
{
// osg::notify(osg::NOTICE)<<"Confused edge found
"<<edgeFaces.size()<<std::endl;
}
}
// compute center
osg::Vec3d center;
VertexSet::iterator vitr;
for(vitr = uniqueVertices.begin();
vitr != uniqueVertices.end();
++vitr)
{
center += *vitr;
}
center /= (double)(uniqueVertices.size());
// compute the ordered points around silhoette
osg::Vec3d side = ( fabs(normal.x()) < fabs(normal.y()) ) ?
osg::Vec3(1.0, 0.0, 0.0) :
osg::Vec3(0.0, 1.0, 0.0);
osg::Vec3 v = normal ^ side;
v.normalize();
osg::Vec3 u = v ^ normal;
u.normalize();
typedef std::map<double, osg::Vec3d> AnglePositions;
AnglePositions anglePositions;
for(vitr = uniqueVertices.begin();
vitr != uniqueVertices.end();
++vitr)
{
osg::Vec3d delta = *vitr - center;
double angle = atan2(delta * u, delta * v);
if (angle<0.0) angle += 2.0*osg::PI;
anglePositions[angle] = *vitr;
}
for(AnglePositions::iterator aitr = anglePositions.begin();
aitr != anglePositions.end();
++aitr)
{
vertices.push_back(aitr->second);
}
}
void cut(const osg::Polytope& polytope)
{
char * apc[6] = { "left", "right", "bottom", "top", "near", "far" };
char ac[16];
int i = 0;
// osg::notify(osg::NOTICE)<<"Cutting with polytope "<<std::endl;
for(osg::Polytope::PlaneList::const_iterator itr =
polytope.getPlaneList().begin();
itr != polytope.getPlaneList().end();
++itr)
{
cut(*itr, std::string( i < 6 ? apc[i] : itoa( i, ac, 10 ) ) );
i++;
}
}
void cut(const CustomPolytope& polytope)
{
// osg::notify(osg::NOTICE)<<"Cutting with polytope "<<std::endl;
for(Faces::const_iterator itr = polytope._faces.begin();
itr != polytope._faces.end();
++itr)
{
cut(itr->plane, itr->name);
}
}
void cut(const osg::Plane& plane, const std::string& name=std::string())
{
// osg::notify(osg::NOTICE)<<" Cutting plane "<<plane<<std::endl;
Face newFace;
for(Faces::iterator itr = _faces.begin();
itr != _faces.end();
)
{
Face& face = *itr;
int intersect = plane.intersect(face.vertices);
if (intersect==1)
{
// osg::notify(osg::NOTICE)<<" Face inside"<<std::endl;
++itr;
}
else if (intersect==0)
{
// osg::notify(osg::NOTICE)<<" Face intersecting - before
"<<face.vertices.size()<<std::endl;
Vertices& vertices = face.vertices;
Vertices newVertices;
typedef std::vector<double> Distances;
Distances distances;
distances.reserve(face.vertices.size());
for(Vertices::iterator vitr = vertices.begin();
vitr != vertices.end();
++vitr)
{
distances.push_back(plane.distance(*vitr));
}
for(unsigned int i=0; i<vertices.size(); ++i)
{
osg::Vec3d& va = vertices[i];
osg::Vec3d& vb = vertices[(i+1)%vertices.size()];
double distance_a = distances[i];
double distance_b = distances[(i+1)%vertices.size()];
// is first edge point inside plane?
if (distance_a>=0.0) newVertices.push_back(vertices[i]);
if (distance_a==0.0)
newFace.vertices.push_back(vertices[i]);
// does edge intersect plane
if (distance_a * distance_b<0.0)
{
// inserting vertex
osg::Vec3d intersection = (vb*distance_a -
va*distance_b)/(distance_a-distance_b);
newVertices.push_back(intersection);
newFace.vertices.push_back(intersection);
// osg::notify(osg::NOTICE)<<" intersection distance
"<<plane.distance(intersection)<<std::endl;
}
}
vertices.swap(newVertices);
// osg::notify(osg::NOTICE)<<" intersecting - after
"<<face.vertices.size()<<std::endl;
++itr;
}
else if (intersect==-1)
{
// osg::notify(osg::NOTICE)<<" Face
outside"<<_faces.size()<<std::endl;
itr = _faces.erase(itr);
}
}
if (!newFace.vertices.empty())
{
// osg::notify(osg::NOTICE)<<" inserting newFace intersecting
"<<newFace.vertices.size()<<std::endl;
newFace.name = name;
newFace.plane = plane;
Vertices& vertices = newFace.vertices;
osg::Vec3d side = ( fabs(plane.getNormal().x()) <
fabs(plane.getNormal().y()) ) ?
osg::Vec3(1.0, 0.0, 0.0) :
osg::Vec3(0.0, 1.0, 0.0);
osg::Vec3 v = plane.getNormal() ^ side;
v.normalize();
osg::Vec3 u = v ^ plane.getNormal();
u.normalize();
osg::Vec3d center;
Vertices::iterator vitr;
for(vitr = vertices.begin();
vitr != vertices.end();
++vitr)
{
center += *vitr;
}
center /= vertices.size();
typedef std::map<double, osg::Vec3d> AnglePositions;
AnglePositions anglePositions;
for(vitr = vertices.begin();
vitr != vertices.end();
++vitr)
{
osg::Vec3d delta = *vitr - center;
double angle = atan2(delta * u, delta * v);
if (angle<0.0) angle += 2.0*osg::PI;
anglePositions[angle] = *vitr;
}
Vertices newVertices;
for(AnglePositions::iterator aitr = anglePositions.begin();
aitr != anglePositions.end();
++aitr)
{
newVertices.push_back(aitr->second);
}
newVertices.swap(vertices);
// osg::notify(osg::NOTICE)<<" after angle sort
"<<newFace.vertices.size()<<std::endl;
_faces.push_back(newFace);
}
}
void getPolytope(osg::Polytope& polytope)
{
for(Faces::const_iterator itr = _faces.begin();
itr != _faces.end();
++itr)
{
polytope.add(itr->plane);
}
}
void getPoints(Vertices& vertices)
{
typedef std::set<osg::Vec3d> VerticesSet;
VerticesSet verticesSet;
for(Faces::iterator itr = _faces.begin();
itr != _faces.end();
++itr)
{
Face& face = *itr;
for(Vertices::iterator vitr = face.vertices.begin();
vitr != face.vertices.end();
++vitr)
{
verticesSet.insert(*vitr);
}
}
for(VerticesSet::iterator sitr = verticesSet.begin();
sitr != verticesSet.end();
++sitr)
{
vertices.push_back(*sitr);
}
}
osg::Drawable* createDrawable(const osg::Vec4d& colour)
{
osg::Geometry* geometry = new osg::Geometry;
osg::Vec3Array* vertices = new osg::Vec3Array;
geometry->setVertexArray(vertices);
for(Faces::iterator itr = _faces.begin();
itr != _faces.end();
++itr)
{
Face& face = *itr;
geometry->addPrimitiveSet( new osg::DrawArrays(GL_LINE_LOOP,
vertices->size(), face.vertices.size()) );
for(Vertices::iterator vitr = face.vertices.begin();
vitr != face.vertices.end();
++vitr)
{
vertices->push_back(*vitr);
}
}
osg::Vec4Array* colours = new osg::Vec4Array;
colours->push_back(colour);
geometry->setColorArray(colours);
geometry->setColorBinding(osg::Geometry::BIND_OVERALL);
osg::StateSet* stateset = geometry->getOrCreateStateSet();
stateset->setMode(GL_LIGHTING, osg::StateAttribute::OFF);
stateset->setTextureMode(0, GL_TEXTURE_2D, osg::StateAttribute::OFF);
stateset->setTextureMode(1, GL_TEXTURE_2D, osg::StateAttribute::OFF);
return geometry;
}
typedef std::list<Face> Faces;
Faces _faces;
protected:
};
#endif_______________________________________________
osg-users mailing list
osg-users@lists.openscenegraph.org
http://lists.openscenegraph.org/listinfo.cgi/osg-users-openscenegraph.org
