C++ Triangles类代码示例

void Mesh::Load(
	lpxmlnode						pNode)
{
	lpxmlnode pCurrNode = pNode->first_node();

	while(pCurrNode != NULL) {
		std::string Name = pCurrNode->name();
		if(Name == "source") {
			Source src;
			src.Load(pCurrNode);
			m_source.push_back(src);
		} else if(Name == "polylist") {
			Polylist pl;
			pl.Load(pCurrNode);
			m_polylist.push_back(pl);
		} else if(Name == "triangles") {
			Triangles triangles;
			triangles.Load(pCurrNode);
			m_triangles.push_back(triangles);
		} else if(Name == "vertices") {
			Vertices vertices;
			m_vertices.Load(pCurrNode);
		}
		pCurrNode = pCurrNode->next_sibling();
	};
}

double
NBHeightMapper::getZ(const Position& geo) const {
    if (!ready()) {
        WRITE_WARNING("Cannot supply height since no height data was loaded");
        return 0;
    }
    if (myRaster != 0) {
        double result = -1e6;
        if (myBoundary.around(geo)) {
            const int xSize = int((myBoundary.xmax() - myBoundary.xmin()) / mySizeOfPixel.x() + .5);
            const double normX = (geo.x() - myBoundary.xmin()) / mySizeOfPixel.x();
            const double normY = (geo.y() - myBoundary.ymax()) / mySizeOfPixel.y();
            PositionVector corners;
            corners.push_back(Position(floor(normX) + 0.5, floor(normY) + 0.5, myRaster[(int)normY * xSize + (int)normX]));
            if (normX - floor(normX) > 0.5) {
                corners.push_back(Position(floor(normX) + 1.5, floor(normY) + 0.5, myRaster[(int)normY * xSize + (int)normX + 1]));
            } else {
                corners.push_back(Position(floor(normX) - 0.5, floor(normY) + 0.5, myRaster[(int)normY * xSize + (int)normX - 1]));
            }
            if (normY - floor(normY) > 0.5) {
                corners.push_back(Position(floor(normX) + 0.5, floor(normY) + 1.5, myRaster[((int)normY + 1) * xSize + (int)normX]));
            } else {
                corners.push_back(Position(floor(normX) + 0.5, floor(normY) - 0.5, myRaster[((int)normY - 1) * xSize + (int)normX]));
            }
            result = Triangle(corners).getZ(Position(normX, normY));
        }
        if (result > -1e5 && result < 1e5) {
            return result;
        }
    }
    // coordinates in degrees hence a small search window
    float minB[2];
    float maxB[2];
    minB[0] = (float)geo.x() - 0.00001f;
    minB[1] = (float)geo.y() - 0.00001f;
    maxB[0] = (float)geo.x() + 0.00001f;
    maxB[1] = (float)geo.y() + 0.00001f;
    QueryResult queryResult;
    int hits = myRTree.Search(minB, maxB, queryResult);
    Triangles result = queryResult.triangles;
    assert(hits == (int)result.size());
    UNUSED_PARAMETER(hits); // only used for assertion

    for (Triangles::iterator it = result.begin(); it != result.end(); it++) {
        const Triangle* triangle = *it;
        if (triangle->contains(geo)) {
            return triangle->getZ(geo);
        }
    }
    WRITE_WARNING("Could not get height data for coordinate " + toString(geo));
    return 0;
}

int main(int argc, char*argv[])
{
  std::ifstream in((argc>1)?argv[1]:"data/triangles.xyz");
  Triangles triangles;
  Triangle_3 t;
  while(in >> t){
    triangles.push_back(t);
  }
  // Create the corresponding vector of bounding boxes
  std::vector<Box> boxes;
  for ( Iterator i = triangles.begin(); i != triangles.end(); ++i)
    boxes.push_back( Box( i->bbox(), i));
  
  // Create the corresponding vector of pointers to bounding boxes
  std::vector<Box *> ptr;
  for ( std::vector<Box>::iterator i = boxes.begin(); i != boxes.end(); ++i)
    ptr.push_back( &*i);
  
  // Run the self intersection algorithm with all defaults on the
  // indirect pointers to bounding boxes. Avoids copying the boxes.
  CGAL::box_self_intersection_d( ptr.begin(), ptr.end(), Report(triangles));
  return 0;
}

Mesh Mesh::fromOBJ(const std::string &filename, bool centralizeLoadedMesh)
{
    std::ifstream in(filename);
    if (!in.is_open())
        throw FACELIB_EXCEPTION("Can't open file " + filename);

    VectorOfPoints points;
    Triangles triangles;
    std::string line;
    while (std::getline(in, line))
    {
        if (line.empty()) continue;

        if (line[0] == 'v')
        {
            Poco::StringTokenizer tokens(line, " ");
            double x = Poco::NumberParser::parseFloat(tokens[1]);
            double y = Poco::NumberParser::parseFloat(tokens[2]);
            double z = Poco::NumberParser::parseFloat(tokens[3]);

            points.push_back(cv::Point3d(x,y,z));
        }
        else if (line[0] == 'f')
        {
            Poco::StringTokenizer tokens(line, " ");
            int t1 = Poco::NumberParser::parse(tokens[1]) - 1;
            int t2 = Poco::NumberParser::parse(tokens[2]) - 1;
            int t3 = Poco::NumberParser::parse(tokens[3]) - 1;

            triangles.push_back(cv::Vec3i(t1, t2, t3));
        }
    }

    Mesh result = Mesh::fromPointcloud(points, centralizeLoadedMesh, false);
    result.triangles = triangles;
    return result;
}

      static IGL_INLINE bool intersect_other_helper(
        const Eigen::PlainObjectBase<DerivedVA> & VA,
        const Eigen::PlainObjectBase<DerivedFA> & FA,
        const Eigen::PlainObjectBase<DerivedVB> & VB,
        const Eigen::PlainObjectBase<DerivedFB> & FB,
        const RemeshSelfIntersectionsParam & params,
        Eigen::PlainObjectBase<DerivedIF> & IF,
        Eigen::PlainObjectBase<DerivedVVAB> & VVAB,
        Eigen::PlainObjectBase<DerivedFFAB> & FFAB,
        Eigen::PlainObjectBase<DerivedJAB>  & JAB,
        Eigen::PlainObjectBase<DerivedIMAB> & IMAB)
      {

        using namespace std;
        using namespace Eigen;

        typedef typename DerivedFA::Index Index;
        // 3D Primitives
        typedef CGAL::Point_3<Kernel>    Point_3;
        typedef CGAL::Segment_3<Kernel>  Segment_3; 
        typedef CGAL::Triangle_3<Kernel> Triangle_3; 
        typedef CGAL::Plane_3<Kernel>    Plane_3;
        typedef CGAL::Tetrahedron_3<Kernel> Tetrahedron_3; 
        // 2D Primitives
        typedef CGAL::Point_2<Kernel>    Point_2;
        typedef CGAL::Segment_2<Kernel>  Segment_2; 
        typedef CGAL::Triangle_2<Kernel> Triangle_2; 
        // 2D Constrained Delaunay Triangulation types
        typedef CGAL::Triangulation_vertex_base_2<Kernel>  TVB_2;
        typedef CGAL::Constrained_triangulation_face_base_2<Kernel> CTAB_2;
        typedef CGAL::Triangulation_data_structure_2<TVB_2,CTAB_2> TDS_2;
        typedef CGAL::Exact_intersections_tag Itag;
        // Axis-align boxes for all-pairs self-intersection detection
        typedef std::vector<Triangle_3> Triangles;
        typedef typename Triangles::iterator TrianglesIterator;
        typedef typename Triangles::const_iterator TrianglesConstIterator;
        typedef 
          CGAL::Box_intersection_d::Box_with_handle_d<double,3,TrianglesIterator> 
          Box;
        typedef 
          std::map<Index,std::vector<std::pair<Index,CGAL::Object> > > 
          OffendingMap;
        typedef std::map<std::pair<Index,Index>,std::vector<Index> >  EdgeMap;
        typedef std::pair<Index,Index> EMK;

        Triangles TA,TB;
        // Compute and process self intersections
        mesh_to_cgal_triangle_list(VA,FA,TA);
        mesh_to_cgal_triangle_list(VB,FB,TB);
        // http://www.cgal.org/Manual/latest/doc_html/cgal_manual/Box_intersection_d/Chapter_main.html#Section_63.5 
        // Create the corresponding vector of bounding boxes
        std::vector<Box> A_boxes,B_boxes;
        const auto box_up = [](Triangles & T, std::vector<Box> & boxes) -> void
        {
          boxes.reserve(T.size());
          for ( 
            TrianglesIterator tit = T.begin(); 
            tit != T.end(); 
            ++tit)
          {
            boxes.push_back(Box(tit->bbox(), tit));
          }
        };
        box_up(TA,A_boxes);
        box_up(TB,B_boxes);
        OffendingMap offendingA,offendingB;
        //EdgeMap edge2facesA,edge2facesB;

        std::list<int> lIF;
        const auto cb = [&](const Box &a, const Box &b) -> void
        {
          using namespace std;
          // index in F and T
          int fa = a.handle()-TA.begin();
          int fb = b.handle()-TB.begin();
          const Triangle_3 & A = *a.handle();
          const Triangle_3 & B = *b.handle();
          if(CGAL::do_intersect(A,B))
          {
            // There was an intersection
            lIF.push_back(fa);
            lIF.push_back(fb);
            if(params.first_only)
            {
              throw IGL_FIRST_HIT_EXCEPTION;
            }
            if(!params.detect_only)
            {
              CGAL::Object result = CGAL::intersection(A,B);

              push_result(FA,fa,fb,result,offendingA);
              push_result(FB,fb,fa,result,offendingB);
            }
          }
        };
        try{
          CGAL::box_intersection_d(
            A_boxes.begin(), A_boxes.end(),
            B_boxes.begin(), B_boxes.end(),
            cb);
//.........这里部分代码省略.........

//----------------------------------------------------------------------------
void Picker::ExecuteRecursive (Movable* object, bool &hasMeshPicked)
{
	if (object)
	{
		if (object->Culling == Movable::CULL_ALWAYS)
			return;
	}

	Triangles* mesh = DynamicCast<Triangles>(object);
	if (mesh)
	{
		if (mesh->WorldBound.TestIntersection(mOrigin, mDirection, mTMin, mTMax))
		{
			// 将射线从世界坐标系转换到模型坐标系。
			APoint ptmp = mesh->WorldTransform.Inverse()*mOrigin;
			Vector3f modelOrigin(ptmp[0], ptmp[1], ptmp[2]);

			AVector vtmp = mesh->WorldTransform.Inverse()*mDirection;
			Vector3f modelDirection(vtmp[0], vtmp[1], vtmp[2]);

			Line3f line(modelOrigin, modelDirection);

			// 访问方位数据
			VertexBufferAccessor vba(mesh);

			int numTriangles = mesh->GetNumTriangles();
			for (int i = 0; i < numTriangles; ++i)
			{
				int v0, v1, v2;
				if (!mesh->GetTriangle(i, v0, v1, v2))
				{
					continue;
				}

				Vector3f vertex0 = vba.Position<Vector3f>(v0);
				Vector3f vertex1 = vba.Position<Vector3f>(v1);
				Vector3f vertex2 = vba.Position<Vector3f>(v2);
				Triangle3f triangle(vertex0, vertex1, vertex2);

				IntrLine3Triangle3f calc(line, triangle);
				if (calc.Find() && mTMin <= calc.GetLineParameter()
					&&  calc.GetLineParameter() <= mTMax)
				{
					PickRecord record;
					record.Intersected = mesh;
					record.T = calc.GetLineParameter();
					record.Triangle = i;
					record.Bary[0] = calc.GetTriBary0();
					record.Bary[1] = calc.GetTriBary1();
					record.Bary[2] = calc.GetTriBary2();
					Records.push_back(record);

					if (mIsDoMovPickCall)
					{
						hasMeshPicked = true;
						mesh->OnPicked(mPickInfo);
					}
				}
			}
		}
		else
		{
			if (mIsDoMovPickCall)
				mesh->OnNotPicked(mPickInfo);
		}
		return;
	}

	SwitchNode* switchNode = DynamicCast<SwitchNode>(object);
	if (switchNode)
	{
		int activeChild = switchNode->GetActiveChild();
		if (activeChild != SwitchNode::SN_INVALID_CHILD)
		{
			if (switchNode->WorldBound.TestIntersection(mOrigin, mDirection, mTMin, mTMax))
			{
				Movable* child = switchNode->GetChild(activeChild);
				if (child)
				{
					ExecuteRecursive(child, hasMeshPicked);
				}

				if (mIsDoMovPickCall)
				{
					if (hasMeshPicked)
					{
						switchNode->OnPicked(mPickInfo);
					}
					else
					{
						switchNode->OnNotPicked(mPickInfo);
					}
				}
			}
			else
			{
				if (mIsDoMovPickCall)
					switchNode->OnNotPicked(mPickInfo);
			}
//.........这里部分代码省略.........

void KeyEdge::triangulate_(double width, Time time, Triangles & out) const
{
    out.clear();
    if (exists(time))
        geometry()->triangulate(width, out);
}

//----------------------------------------------------------------------------
void Picker::ExecuteRecursive (Spatial* object)
{
    Triangles* mesh = DynamicCast<Triangles>(object);
    if (mesh)
    {
        if (mesh->WorldBound.TestIntersection(mOrigin, mDirection, mTMin,
            mTMax))
        {
            // Convert the linear component to model-space coordinates.
            APoint ptmp = mesh->WorldTransform.Inverse()*mOrigin;
            Vector3f modelOrigin(ptmp[0], ptmp[1], ptmp[2]);

            AVector vtmp = mesh->WorldTransform.Inverse()*mDirection;
            Vector3f modelDirection(vtmp[0], vtmp[1], vtmp[2]);

            Line3f line(modelOrigin, modelDirection);

            // Get the position data.
            VertexBufferAccessor vba(mesh);

            // Compute intersections with the model-space triangles.
            int numTriangles = mesh->GetNumTriangles();
            for (int i = 0; i < numTriangles; ++i)
            {
                int v0, v1, v2;
                if (!mesh->GetTriangle(i, v0, v1, v2))
                {
                    continue;
                }

                Vector3f vertex0 = vba.Position<Vector3f>(v0);
                Vector3f vertex1 = vba.Position<Vector3f>(v1);
                Vector3f vertex2 = vba.Position<Vector3f>(v2);
                Triangle3f triangle(vertex0, vertex1, vertex2);

                IntrLine3Triangle3f calc(line, triangle);
                if (calc.Find() && mTMin <= calc.GetLineParameter()
                &&  calc.GetLineParameter() <= mTMax)
                {
                    PickRecord record;
                    record.Intersected = mesh;
                    record.T = calc.GetLineParameter();
                    record.Triangle = i;
                    record.Bary[0] = calc.GetTriBary0();
                    record.Bary[1] = calc.GetTriBary1();
                    record.Bary[2] = calc.GetTriBary2();
                    Records.push_back(record);
                }
            }
        }
        return;
    }

    SwitchNode* switchNode = DynamicCast<SwitchNode>(object);
    if (switchNode)
    {
        int activeChild = switchNode->GetActiveChild();
        if (activeChild != SwitchNode::SN_INVALID_CHILD)
        {
            if (switchNode->WorldBound.TestIntersection(mOrigin,
                mDirection, mTMin, mTMax))
            {
                Spatial* child = switchNode->GetChild(activeChild);
                if (child)
                {
                    ExecuteRecursive(child);
                }
            }
        }
        return;
    }

    Node* node = DynamicCast<Node>(object);
    if (node)
    {
        if (node->WorldBound.TestIntersection(mOrigin, mDirection, mTMin,
            mTMax))
        {
            for (int i = 0; i < node->GetNumChildren(); ++i)
            {
                Spatial* child = node->GetChild(i);
                if (child)
                {
                    ExecuteRecursive(child);
                }
            }
        }
    }
}

NiSkinPartition::NiSkinPartition(Ref<NiTriBasedGeom> shape) {
   NiSkinInstanceRef skinInst = shape->GetSkinInstance();
   if ( skinInst == NULL ) {
      throw runtime_error( "You must bind a skin before setting generating skin partitions.  No NiSkinInstance found." );
   }
   NiSkinDataRef skinData = skinInst->GetSkinData();
   if ( skinData == NULL ) {
      throw runtime_error( "You must bind a skin before setting generating skin partitions.  No NiSkinData found." );
   }
   NiTriBasedGeomDataRef geomData = DynamicCast<NiTriBasedGeomData>( shape->GetData() );
   if ( geomData == NULL ) {
      throw runtime_error( "Attempted to generate a skin partition on a mesh with no geometry data." );
   }

   int nWeightsPerVertex = 4;
   vector<WeightList> vertexWeights;
   BoneList boneMap;
   vector<unsigned short> vertexMap;
   Strips strips;
   vector<BoneList> boneIndexList;
   Triangles triangles;

   int totalBones = skinInst->GetBoneCount();
   boneMap.resize(totalBones);

   int nv = geomData->GetVertexCount();
   vertexMap.resize(nv);
   vertexWeights.resize(nv);
   boneIndexList.resize(nv);

   for (int i=0; i<totalBones; ++i) {
      boneMap[i] = i;

      vector<SkinWeight> skinWeights = skinData->GetBoneWeights(i);
      for (vector<SkinWeight>::const_iterator skinWeight = skinWeights.begin(); skinWeight != skinWeights.end(); ++skinWeight) {
         WeightList& vertexWeight = vertexWeights[skinWeight->index];
         BoneList& boneIndex = boneIndexList[skinWeight->index];

         vertexWeight.push_back(skinWeight->weight);
         boneIndex.push_back(i);

         // Adjust upper limit on number of weights per vertex if necessary.
         int nWeights = vertexWeight.size();
         if (nWeights > nWeightsPerVertex)
            nWeightsPerVertex = nWeights;
      }
   }
   if (nWeightsPerVertex == 0) {
      throw runtime_error( "Attempted to generate a skin partition on a mesh with no weights specified." );
   }
   for (int i=0; i<nv; ++i) {
      vertexMap[i] = i;

      WeightList& vertexWeight = vertexWeights[i];
      BoneList& boneIndex = boneIndexList[i];
      vertexWeight.reserve(nWeightsPerVertex);
      boneIndex.reserve(nWeightsPerVertex);
      for (int j = nWeightsPerVertex - vertexWeight.size(); j>0; --j) {
         vertexWeight.push_back(0.0f);
         boneIndex.push_back(0);
      }
   }

   SetNumPartitions(1);
   SetWeightsPerVertex(0, nWeightsPerVertex);
   SetBoneMap(0, boneMap);
   SetNumVertices(0, (unsigned short)(vertexMap.size()) );
   SetVertexMap(0, vertexMap);
   EnableVertexWeights(0, true);
   EnableVertexBoneIndices(0, true);
   for (int i=0; i<nv; ++i) {
      SetVertexWeights(0, i, vertexWeights[i]);
      SetVertexBoneIndices(0, i, boneIndexList[i]);
   }

   // Special case for pre-stripped data
   if (NiTriStripsDataRef stripData = DynamicCast<NiTriStripsData>(geomData)) {
      unsigned short nstrips = stripData->GetStripCount();
      SetStripCount(0, nstrips);
      for (int i=0; i<int(nstrips); ++i) {
         SetStrip(0, i, stripData->GetStrip(i));
      }
   } else {

      Triangles triangles = geomData->GetTriangles();
      SetTriangles(0, triangles);

      unsigned short *data = new unsigned short[triangles.size() * 3 * 2];
      for (size_t i=0; i< triangles.size(); i++) {
         data[i * 3 + 0] = triangles[i][0];
         data[i * 3 + 1] = triangles[i][1];
         data[i * 3 + 2] = triangles[i][2];
      }
      PrimitiveGroup * groups = 0;
      unsigned short numGroups = 0;

      // GF 3+
      SetCacheSize(CACHESIZE_GEFORCE3);
      // don't generate hundreds of strips
      SetStitchStrips(true);
//.........这里部分代码省略.........

void EdgeCell::triangulate(double /*width*/, Time /*time*/, Triangles & out)
{
    out.clear();
}

int main()
{
   Triangles application; // create Triangles object 
   application.drawTriangles(); // call its drawTriangles function
} // end main

//----------------------------------------------------------------------------
void Picker::ExecuteRecursive (Movable* object, bool &hasMeshPicked)
{
	if (object)
	{
		if (!object->IsDoPick())
			return;

		if (!object->IsShow() && !object->IsPickIngoreCullingMode())
			return;
	}

	Triangles* mesh = DynamicCast<Triangles>(object);
	if (mesh)
	{
		if (!mesh->GetVertexBuffer())
			return;

		if (mesh->WorldBound.TestIntersection(mOrigin, mDirection, mTMin, mTMax))
		{
			if (mesh->IsUseBoundPick())
			{
				AVector dir = mesh->GetWorldTransform().GetTranslate() - mOrigin;
				float length = dir.Length();

				PickRecord record;
				record.Intersected = mesh;
				record.T = length;
				Records.push_back(record);
			}
			else
			{
				// 将射线从世界坐标系转换到模型坐标系。
				APoint ptmp;
				if (!mesh->IsSkinCtrlSetWroldTrans)
					ptmp = mesh->WorldTransform.Inverse()*mOrigin;
				else
					ptmp = mesh->BoundWorldTransform.Inverse()*mOrigin;
				Vector3f modelOrigin(ptmp[0], ptmp[1], ptmp[2]);

				AVector vtmp;
				if (!mesh->IsSkinCtrlSetWroldTrans)
					vtmp = mesh->WorldTransform.Inverse()*mDirection;
				else
					vtmp = mesh->BoundWorldTransform.Inverse()*mDirection;
				Vector3f modelDirection(vtmp[0], vtmp[1], vtmp[2]);

				Line3f line(modelOrigin, modelDirection);

				// 访问方位数据
				VertexBufferAccessor vba(mesh);

				int numTriangles = mesh->GetNumTriangles();
				for (int i = 0; i < numTriangles; ++i)
				{
					int v0, v1, v2;
					if (!mesh->GetTriangle(i, v0, v1, v2))
					{
						continue;
					}

					Vector3f vertex0 = vba.Position<Vector3f>(v0);
					Vector3f vertex1 = vba.Position<Vector3f>(v1);
					Vector3f vertex2 = vba.Position<Vector3f>(v2);
					Triangle3f triangle(vertex0, vertex1, vertex2);

					IntrLine3Triangle3f calc(line, triangle);
					if (calc.Find())
					{
						float lineParameter = calc.GetLineParameter();
						if (mTMin<=lineParameter && lineParameter<=mTMax)
						{
							PickRecord record;
							record.Intersected = mesh;
							record.T = calc.GetLineParameter();
							record.Triangle = i;
							record.Bary[0] = calc.GetTriBary0();
							record.Bary[1] = calc.GetTriBary1();
							record.Bary[2] = calc.GetTriBary2();

							Vector3f edg1 = vertex1 - vertex0;
							Vector3f edg2 = vertex2 - vertex0;
							Vector3f normal = edg1.UnitCross(edg2);
							float dotVal = line.Direction.Dot(normal);
							if (dotVal > Mathf::ZERO_TOLERANCE)
							{
								normal *= -1.0f;
							}
							record.LocalNormal = normal;
							record.WorldPos = mOrigin + mDirection * record.T;

							Records.push_back(record);

							if (mIsDoMovPickCall)
							{
								hasMeshPicked = true;
								mesh->OnPicked(mPickInfo);
							}
						}
					}
//.........这里部分代码省略.........

// callback function that reports all truly intersecting triangles
void report_inters( const Box* a, const Box* b) {
    std::cout << "Box " << (a->handle() - triangles.begin()) << " and "
              << (b->handle() - triangles.begin()) << " intersect";
    if ( ! a->handle()->is_degenerate() && ! b->handle()->is_degenerate()
         && CGAL::do_intersect( *(a->handle()), *(b->handle()))) {
        std::cout << ", and the triangles intersect also";
    }
    std::cout << '.' << std::endl;
}

bool detect(vector<Triangle_3> &a, vector<Triangle_3> &b)
{
    std::vector<Box> boxes;
    triangles.clear();
    for ( Iterator i = a.begin(); i != a.end(); ++i)
        triangles.push_back(*i);
    for ( Iterator i = b.begin(); i != b.end(); ++i)
        triangles.push_back(*i);

    for(Iterator i = triangles.begin(); i!= triangles.end(); ++i)
        boxes.push_back( Box( i->bbox(), i));

    // Run the self intersection algorithm with all defaults
    CGAL::box_self_intersection_d( boxes.begin(), boxes.end(), report_inters);
   return true;
}

int main() {
    // Create 10 random triangles
    typedef CGAL::Random_points_in_cube_3<Point_3>           Pts;
    typedef CGAL::Creator_uniform_3< Point_3, Triangle_3>    Creator;
    typedef CGAL::Join_input_iterator_3<Pts,Pts,Pts,Creator> Triangle_gen;
    Pts    points( 1); // in centered cube [-1,1)^3
    Triangle_gen triangle_gen( points, points, points);
    CGAL::cpp11::copy_n( triangle_gen, 10, std::back_inserter(triangles));

    // Create the corresponding vector of bounding boxes
    std::vector<Box> boxes;
    for ( Iterator i = triangles.begin(); i != triangles.end(); ++i)
        boxes.push_back( Box( i->bbox(), i));

    // Create the corresponding vector of pointers to bounding boxes
    std::vector<Box *> ptr;
    for ( std::vector<Box>::iterator i = boxes.begin(); i != boxes.end(); ++i)
        ptr.push_back( &*i);

    // Run the self intersection algorithm with all defaults on the
    // indirect pointers to bounding boxes. Avoids copying the boxes.
    CGAL::box_self_intersection_d( ptr.begin(), ptr.end(), report_inters);
    return 0;
}