CGAL users discussion list

[<>]

Hello, 
I am a newbie in CGAL and I have a question concerning Alphashapes.

I have set of about 1 million points, containing several roof planes. i would
like to get the shape of each roof plane. is it possible to get the shape by
using alphashapes.

I wrote the following function, but I don know if it is the right way to get
the shape of each roof plane.

I don`t have much experience using alphashapes and CGAL, but maybe somebody of
you has.

void alphashapes(struct points input[], const char *outputfile, int
pointCount){

//open file
                ofstream outfile(outputfile, ios::app);
                
                if(!outfile) { // file couldn't be opened
                cerr << "Error: file could not be opened" << endl;
                exit(1);
                }


//variables

        Delaunay_hierarchy dt;
        Point p;
        std::list<Edge> edges;
        std::list<Edge>::iterator it;
        
        for(int i = 0; i<pointCount; i++){
                int out = i;
                p=Point((double)input[out].X, (double)input[out].Y,
(double)input[out].Z); 
                
                dt.insert(p);
        }

        std::cout << "Delaunay computed." << std::endl;

        // compute alpha shape  
        Alpha_shape_3 as(dt,);
        
        std::cout << "Alpha shape computed in REGULARIZED mode by defaut." 
                << std::endl;

// find optimal alpha values
        Alpha_shape_3::NT alpha_solid = as.find_alpha_solid();
        Alpha_iterator opt = as.find_optimal_alpha(1);
        std::cout << "Smallest alpha value to get a solid through data points
is " 
                << alpha_solid << std::endl;
        std::cout << "Optimal alpha value to get one connected component is " 
                <<  *opt    << std::endl;
        as.set_alpha(*opt+5.5);
        as.get_alpha_shape_edges(std::back_inserter(edges),
Alpha_shape_3::REGULAR);

        cout << "edges Size: " << edges.size() << endl;  
        

        for(it = edges.begin(); it != edges.end(); ++it) {
                Point p1=(*it).first->vertex((*it).second)->point();
                Point q=(*it).first->vertex((*it).third)->point();
                cout << fixed;
                cout << setw(2) << setprecision(2);
                cout<< p1 << endl;
                cout << q << endl;
                
                
                
                }

        
        as.clear();
        edges.clear();
        dt.clear();
        outfile.close();


}       



Thanks.

Andi