Subversion Repositories hepmc

//////////////////////////////////////////////////////////////////////////

// Matt.Dobbs@Cern.CH, September 1999

// GenVertex within an event

//////////////////////////////////////////////////////////////////////////

#include "HepMC/GenParticle.h"

#include "HepMC/GenVertex.h"

#include "HepMC/GenEvent.h"

#include <iomanip>       // needed for formatted output

namespace HepMC {

    GenVertex::GenVertex( const FourVector& position,

                          int id, const WeightContainer& weights )

        : m_position(position), m_id(id), m_weights(weights), m_event(0),

          m_barcode(0)

    {

        s_counter++;

    }

    GenVertex::GenVertex( const GenVertex& invertex ) : m_event(0),

                                                        m_barcode(0) {

        // Shallow copy: does not copy the FULL list of particle pointers.

        // Creates a copy of  - invertex

        //                    - outgoing particles of invertex, but sets the

        //                      decay vertex of these particles to NULL

        //                    - all incoming particles which do not have a

        //                      creation vertex.

        // (i.e. it creates copies of all particles which it owns)

        // (note - impossible to copy the FULL list of particle pointers 

        //         while having the vertex

        //         and particles in/out point-back to one another -- unless you

        //         copy the entire tree -- which we don't want to do)

        *this = invertex;

        s_counter++;

    }

    GenVertex::~GenVertex() {

        //

        // need to delete any particles previously owned by this vertex

        if ( parent_event() ) parent_event()->remove_barcode(this);

        delete_adopted_particles();

        s_counter--;

    }

    GenVertex& GenVertex::operator=( const GenVertex& invertex ) {

        // Shallow: does not copy the FULL list of particle pointers.

        // Creates a copy of  - invertex

        //                    - outgoing particles of invertex, but sets the

        //                      decay vertex of these particles to NULL

        //                    - all incoming particles which do not have a

        //                      creation vertex.

        //                    - it does not alter *this's m_event (!)

        // (i.e. it creates copies of all particles which it owns)

        // (note - impossible to copy the FULL list of particle pointers 

        //         while having the vertex

        //         and particles in/out point-back to one another -- unless you

        //         copy the entire tree -- which we don't want to do)

        //

        // first need to delete any particles previously owned by this vertex

        delete_adopted_particles();

        //

        for ( particles_in_const_iterator

                  part1 = invertex.particles_in_const_begin();

              part1 != invertex.particles_in_const_end(); part1++ ) {

            if ( !(*part1)->production_vertex() ) {

                GenParticle* pin = new GenParticle(**part1);

                add_particle_in(pin);

            }

        }

        for ( particles_out_const_iterator

                  part2 = invertex.particles_out_const_begin();

              part2 != invertex.particles_out_const_end(); part2++ ) {

            GenParticle* pin = new GenParticle(**part2);

            add_particle_out(pin);

        }

        set_position( invertex.position() );

        set_id( invertex.id() );

        weights() = invertex.weights();

        suggest_barcode( invertex.barcode() );

        return *this;

    }

    bool GenVertex::operator==( const GenVertex& a ) const {

        // Returns true if the positions and the particles in the lists of a 

        //  and this are identical. Does not compare barcodes.

        // Note that it is impossible for two vertices to point to the same 

        //  particle's address, so we need to do more than just compare the

        //  particle pointers

        //

        if ( a.position() !=  this->position() ) return 0;

        // if the size of the inlist differs, return false.

        if ( a.particles_in_size() !=  this->particles_in_size() ) return 0;

        // if the size of the outlist differs, return false.

        if ( a.particles_out_size() !=  this->particles_out_size() ) return 0;

        // loop over the inlist and ensure particles are identical

        //   (only do this if the lists aren't empty - we already know

        //   if one isn't, then other isn't either!)

        if ( a.particles_in_const_begin() != a.particles_in_const_end() ) {

            for ( GenVertex::particles_in_const_iterator

                      ia = a.particles_in_const_begin(),

                      ib = this->particles_in_const_begin();

                  ia != a.particles_in_const_end(); ia++, ib++ ){

                if ( **ia != **ib ) return 0;

            }

        }

        // loop over the outlist and ensure particles are identical

        //   (only do this if the lists aren't empty - we already know

        //   if one isn't, then other isn't either!)

        if ( a.particles_out_const_begin() != a.particles_out_const_end() ) {

            for ( GenVertex::particles_out_const_iterator

                      ia = a.particles_out_const_begin(),

                      ib = this->particles_out_const_begin();

                  ia != a.particles_out_const_end(); ia++, ib++ ){

                if ( **ia != **ib ) return 0;

            }

        }

        return 1;

    }

    bool GenVertex::operator!=( const GenVertex& a ) const {

        // Returns true if the positions and lists of a and this are not equal.

        return !( a == *this );

    }  

    void GenVertex::print( std::ostream& ostr ) const {

        if ( barcode()!=0 ) {

            if ( position() != FourVector(0,0,0,0) ) {

                ostr << "Vertex:";

                ostr.width(9);

                ostr << barcode();

                ostr << " ID:";

                ostr.width(5);

                ostr << id();

                ostr << " (X,cT)=";

                ostr.width(9);

                ostr.precision(2);

                ostr.setf(std::ios::scientific, std::ios::floatfield);

                ostr.setf(std::ios_base::showpos);

                ostr << position().x() << ",";

                ostr.width(9);

                ostr.precision(2);

                ostr << position().y() << ",";

                ostr.width(9);

                ostr.precision(2);

                ostr << position().z() << ",";

                ostr.width(9);

                ostr.precision(2);

                ostr << position().t();

                ostr.setf(std::ios::fmtflags(0), std::ios::floatfield);

                ostr.unsetf(std::ios_base::showpos);

                ostr << std::endl;

            } else {

                ostr << "GenVertex:";

                ostr.width(9);

                ostr << barcode();

                ostr << " ID:";

                ostr.width(5);

                ostr << id();

                ostr << " (X,cT):0";

                ostr << std::endl;

            }

        } else {

            // If the vertex doesn't have a unique barcode assigned, then

            //  we print its memory address instead... so that the

            //  print out gives us a unique tag for the particle.

            if ( position() != FourVector(0,0,0,0) ) {

                ostr << "Vertex:";

                ostr.width(9);

                ostr << (void*)this;

                ostr << " ID:";

                ostr.width(5);

                ostr << id();

                ostr << " (X,cT)=";

                ostr.width(9);

                ostr.precision(2);

                ostr.setf(std::ios::scientific, std::ios::floatfield);

                ostr.setf(std::ios_base::showpos);

                ostr << position().x();

                ostr.width(9);

                ostr.precision(2);

                ostr << position().y();

                ostr.width(9);

                ostr.precision(2);

                ostr << position().z();

                ostr.width(9);

                ostr.precision(2);

                ostr << position().t();

                ostr.setf(std::ios::fmtflags(0), std::ios::floatfield);

                ostr.unsetf(std::ios_base::showpos);

                ostr << std::endl;

            } else {

                ostr << "GenVertex:";

                ostr.width(9);

                ostr << (void*)this;

                ostr << " ID:";

                ostr.width(5);

                ostr << id();

                ostr << " (X,cT):0";

                ostr << std::endl;

            }

        }

        // print the weights if there are any

        if ( ! weights().empty() ) {

            ostr << " Wgts(" << weights().size() << ")=";

            for ( WeightContainer::const_iterator wgt = weights().begin();

                  wgt != weights().end(); wgt++ ) { ostr << *wgt << " "; }

            ostr << std::endl;

        }

        // print out all the incoming, then outgoing particles

        for ( particles_in_const_iterator part1 = particles_in_const_begin();

              part1 != particles_in_const_end(); part1++ ) {

            if ( part1 == particles_in_const_begin() ) {

                ostr << " I:";

                ostr.width(2);

                ostr << m_particles_in.size();

            } else { ostr << "     "; }

            //(*part1)->print( ostr );  //uncomment for long debugging printout

            ostr << **part1 << std::endl;

        }

        for ( particles_out_const_iterator part2 = particles_out_const_begin();

              part2 != particles_out_const_end(); part2++ ) {

            if ( part2 == particles_out_const_begin() ) {

                ostr << " O:";

                ostr.width(2);

                ostr << m_particles_out.size();

            } else { ostr << "     "; }

            //(*part2)->print( ostr ); // uncomment for long debugging printout

            ostr << **part2 << std::endl;

        }

    }

    double GenVertex::check_momentum_conservation() const {

        // finds the difference between the total momentum out and the total

        // momentum in vectors, and returns the magnitude of this vector

        // i.e.         returns | \vec{p_in} - \vec{p_out} |

        double sumpx = 0, sumpy = 0, sumpz = 0;

        for ( particles_in_const_iterator part1 = particles_in_const_begin();

              part1 != particles_in_const_end(); part1++ ) {

            sumpx   += (*part1)->momentum().px();

            sumpy   += (*part1)->momentum().py();

            sumpz   += (*part1)->momentum().pz();

        }

        for ( particles_out_const_iterator part2 = particles_out_const_begin();

              part2 != particles_out_const_end(); part2++ ) {

            sumpx   -= (*part2)->momentum().px();

            sumpy   -= (*part2)->momentum().py();

            sumpz   -= (*part2)->momentum().pz();

        }

        return sqrt( sumpx*sumpx + sumpy*sumpy + sumpz*sumpz );

    }

    void GenVertex::add_particle_in( GenParticle* inparticle ) {

        if ( !inparticle ) return;

        // if inparticle previously had a decay vertex, remove it from that

        // vertex's list

        if ( inparticle->end_vertex() ) {

            inparticle->end_vertex()->m_particles_in.erase( inparticle );

        }

        m_particles_in.insert( inparticle );

        inparticle->set_end_vertex_( this );

    }

    void GenVertex::add_particle_out( GenParticle* outparticle ) {

        if ( !outparticle ) return;

        // if outparticle previously had a production vertex,

        // remove it from that vertex's list

        if ( outparticle->production_vertex() ) {

            outparticle->production_vertex()->m_particles_out.erase(

                outparticle );

        }

        m_particles_out.insert( outparticle );

        outparticle->set_production_vertex_( this );

    }

    GenParticle* GenVertex::remove_particle( GenParticle* particle ) {

        // this finds *particle in the in and/or out list and removes it from

        //  these lists ... it DOES NOT DELETE THE PARTICLE or its relations.

        // you could delete the particle too as follows:

        //      delete vtx->remove_particle( particle );

        // or if the particle has an end vertex, you could:

        //      delete vtx->remove_particle( particle )->end_vertex();

        // which would delete the particle's end vertex, and thus would

        // also delete the particle, since the particle would be 

        // owned by the end vertex.

        if ( !particle ) return 0;

        if ( particle->end_vertex() == this ) {

            particle->set_end_vertex_( 0 );

            m_particles_in.erase(particle);

        }

        if ( particle->production_vertex() == this ) {

            particle->set_production_vertex_(0);

            m_particles_out.erase(particle);

        }

        return particle;

    }

    void GenVertex::delete_adopted_particles() {

        // deletes all particles which this vertex owns

        // to be used by the vertex destructor and operator=

        //

        if ( m_particles_out.empty() && m_particles_in.empty() ) return;

        // 1. delete all outgoing particles which don't have decay vertices.

        //    those that do become the responsibility of the decay vertex

        //    and have their productionvertex pointer set to NULL

        for ( std::set<GenParticle*>::iterator part1 = m_particles_out.begin();

              part1 != m_particles_out.end(); ) {

            if ( !(*part1)->end_vertex() ) {

                delete *(part1++);

            } else {

                (*part1)->set_production_vertex_(0);

                ++part1;

            }

        }

        m_particles_out.clear();

        //

        // 2. delete all incoming particles which don't have production

        //    vertices. those that do become the responsibility of the 

        //    production vertex and have their decayvertex pointer set to NULL

        for ( std::set<GenParticle*>::iterator part2 = m_particles_in.begin();

              part2 != m_particles_in.end(); ) {

            if ( !(*part2)->production_vertex() ) {

                delete *(part2++);

            } else {

                (*part2)->set_end_vertex_(0);

                ++part2;

            }

        }

        m_particles_in.clear();

    }

    bool GenVertex::suggest_barcode( int the_bar_code )

    {

        // allows a barcode to be suggested for this vertex.

        // In general it is better to let the event pick the barcode for

        // you, which is automatic.

        // Returns TRUE if the suggested barcode has been accepted (i.e. the

        //  suggested barcode has not already been used in the event, 

        //  and so it was used).

        // Returns FALSE if the suggested barcode was rejected, or if the

        //  vertex is not yet part of an event, such that it is not yet

        //  possible to know if the suggested barcode will be accepted).

        if ( the_bar_code >0 ) {

            std::cerr << "GenVertex::suggest_barcode WARNING, vertex bar codes"

                      << "\n MUST be negative integers. Positive integers "

                      << "\n are reserved for particles only. Your suggestion "

                      << "\n has been rejected." << std::endl;

            return false;

        }

        bool success = false;

        if ( parent_event() ) {

            success = parent_event()->set_barcode( this, the_bar_code );

        } else { set_barcode_( the_bar_code ); }

        return success;

    }

    void GenVertex::set_parent_event_( GenEvent* new_evt )

    {

        GenEvent* orig_evt = m_event;

        m_event = new_evt;

        //

        // every time a vertex's parent event changes, the map of barcodes

        //   in the new and old parent event needs to be modified to 

        //   reflect this

        if ( orig_evt != new_evt ) {

            if (new_evt) new_evt->set_barcode( this, barcode() );

            if (orig_evt) orig_evt->remove_barcode( this );

            // we also need to loop over all the particles which are owned by 

            //  this vertex, and remove their barcodes from the old event.

            for ( particles_in_const_iterator part1=particles_in_const_begin();

                  part1 != particles_in_const_end(); part1++ ) {

                if ( !(*part1)->production_vertex() ) {

                    if ( orig_evt ) orig_evt->remove_barcode( *part1 );

                    if ( new_evt ) new_evt->set_barcode( *part1,

                                                         (*part1)->barcode() );

                }

            }

            for ( particles_out_const_iterator

                      part2 = particles_out_const_begin();

                  part2 != particles_out_const_end(); part2++ ) {

                    if ( orig_evt ) orig_evt->remove_barcode( *part2 );

                    if ( new_evt ) new_evt->set_barcode( *part2,

                                                         (*part2)->barcode() );

            }

        }

    }

    /////////////

    // Static  //

    /////////////

    unsigned int GenVertex::counter() { return s_counter; }

    unsigned int GenVertex::s_counter = 0;

    /////////////

    // Friends //

    /////////////

    std::ostream& operator<<( std::ostream& ostr, const GenVertex& vtx ) {

        if ( vtx.barcode()!=0 ) ostr << "BarCode " << vtx.barcode();

        else ostr << "Address " << &vtx;

        ostr << " (X,cT)=";

        if ( vtx.position() != FourVector(0,0,0,0)) {

            ostr << vtx.position().x() << ","

                 << vtx.position().y() << ","

                 << vtx.position().z() << ","

                 << vtx.position().t();

        } else { ostr << 0; }

        ostr << " #in:" << vtx.particles_in_size()

             << " #out:" << vtx.particles_out_size();

        return ostr;

    }

    /////////////////////////////

    // edge_iterator           // (protected - for internal use only)

    /////////////////////////////

    // If the user wants the functionality of the edge_iterator, he should

    // use particle_iterator with IteratorRange = family, parents, or children

    //

    GenVertex::edge_iterator::edge_iterator() : m_vertex(0), m_range(family),

        m_is_inparticle_iter(0), m_is_past_end(1)

    {}

    GenVertex::edge_iterator::edge_iterator( const GenVertex& vtx,

                                          IteratorRange range ) :

        m_vertex(&vtx), m_range(family)

    {

        // Note: (26.1.2000) the original version of edge_iterator inheritted

        //       from set<GenParticle*>::const_iterator() rather than using

        //       composition as it does now.

        //       The inheritted version suffered from a strange bug, which

        //       I have not fully understood --- it only occurred after many

        //       events were processed and only when I called the delete 

        //       function on past events. I believe it had something to do with

        //       the past the end values, which are now robustly coded in this

        //       version as boolean members.

        //

        // default range is family, only other choices are children/parents

        //    descendants/ancestors not allowed & recasted ot children/parents

        if ( range == descendants || range == children ) m_range = children;

        if ( range == ancestors   || range == parents  ) m_range = parents;

        //

        if ( m_vertex->m_particles_in.empty() &&

             m_vertex->m_particles_out.empty() ) {

            // Case: particles_in and particles_out is empty.

            m_is_inparticle_iter = 0;

            m_is_past_end = 1;

        } else if ( m_range == parents && m_vertex->m_particles_in.empty() ){

            // Case: particles in is empty and parents is requested.

            m_is_inparticle_iter = 1;

            m_is_past_end = 1;

        } else if ( m_range == children && m_vertex->m_particles_out.empty() ){

            // Case: particles out is empty and children is requested.

            m_is_inparticle_iter = 0;

            m_is_past_end = 1;

        } else if ( m_range == children ) {

            // Case: particles out is NOT empty, and children is requested

            m_set_iter = m_vertex->m_particles_out.begin();

            m_is_inparticle_iter = 0;

            m_is_past_end = 0;

        } else if ( m_range == family && m_vertex->m_particles_in.empty() ) {

            // Case: particles in is empty, particles out is NOT empty,

            //       and family is requested. Then skip ahead to partilces out.

            m_set_iter = m_vertex->m_particles_out.begin();

            m_is_inparticle_iter = 0;

            m_is_past_end = 0;

        } else {

            // Normal scenario: start with the first incoming particle

            m_set_iter = m_vertex->m_particles_in.begin();

            m_is_inparticle_iter = 1;

            m_is_past_end = 0;

        }

    }

    GenVertex::edge_iterator::edge_iterator( const edge_iterator& p ) {

        *this = p;

    }

    GenVertex::edge_iterator::~edge_iterator() {}

    GenVertex::edge_iterator& GenVertex::edge_iterator::operator=(

        const edge_iterator& p ) {

        m_vertex = p.m_vertex;

        m_range = p.m_range;

        m_set_iter = p.m_set_iter;

        m_is_inparticle_iter = p.m_is_inparticle_iter;

        m_is_past_end = p.m_is_past_end;

        return *this;

    }

    GenParticle* GenVertex::edge_iterator::operator*(void) const {

        if ( !m_vertex || m_is_past_end ) return 0;

        return *m_set_iter;

    }

    GenVertex::edge_iterator& GenVertex::edge_iterator::operator++(void){

        // Pre-fix increment

        //

        // increment the set iterator (unless we're past the end value)

        if ( m_is_past_end ) return *this;

        ++m_set_iter;

        // handle cases where m_set_iter points past the end

        if ( m_range == family && m_is_inparticle_iter &&

             m_set_iter == m_vertex->m_particles_in.end() ) {

            // at the end on in particle set, and range is family, so move to

            // out particle set

            m_set_iter = m_vertex->m_particles_out.begin();

            m_is_inparticle_iter = 0;

        } else if ( m_range == parents &&

                    m_set_iter == m_vertex->m_particles_in.end() ) {

            // at the end on in particle set, and range is parents only, so

            // move into past the end state

            m_is_past_end = 1;

        }

        // the following is not else if because we might have range=family

        // with an empty particles_out set.     

        if ( m_set_iter == m_vertex->m_particles_out.end() ) {

            //whenever out particles end is reached, go into past the end state

            m_is_past_end = 1;

        }

        return *this;

    }

    GenVertex::edge_iterator GenVertex::edge_iterator::operator++(int){

        // Post-fix increment

        edge_iterator returnvalue = *this;

        ++*this;

        return returnvalue;

    }

    bool GenVertex::edge_iterator::is_parent() const {

        if ( **this && (**this)->end_vertex() == m_vertex ) return 1;

        return 0;

    }

    bool GenVertex::edge_iterator::is_child() const {

        if ( **this && (**this)->production_vertex() == m_vertex ) return 1;

        return 0;

    }

    int GenVertex::edges_size( IteratorRange range ) const {

        if ( range == children ) return m_particles_out.size();

        if ( range == parents ) return m_particles_in.size();

        if ( range == family ) return m_particles_out.size()

                                      + m_particles_in.size();

        return 0;

    }

    /////////////////////

    // vertex_iterator //

    /////////////////////

    GenVertex::vertex_iterator::vertex_iterator()

        : m_vertex(0), m_visited_vertices(0), m_it_owns_set(0),

          m_recursive_iterator(0)

    {}

    GenVertex::vertex_iterator::vertex_iterator( GenVertex& vtx_root,

                                              IteratorRange range )

        : m_vertex(&vtx_root), m_range(range)

    {

        // standard public constructor

        //

        m_visited_vertices = new std::set<const GenVertex*>;

        m_it_owns_set = 1;

        m_visited_vertices->insert( m_vertex );

        m_recursive_iterator = 0;

        m_edge = m_vertex->edges_begin( m_range );

        // advance to the first good return value

        if ( !follow_edge_() &&

             m_edge != m_vertex->edges_end( m_range )) ++*this;

    }

    GenVertex::vertex_iterator::vertex_iterator( GenVertex& vtx_root,

        IteratorRange range, std::set<const GenVertex*>& visited_vertices ) :

        m_vertex(&vtx_root), m_range(range),

        m_visited_vertices(&visited_vertices), m_it_owns_set(0),

        m_recursive_iterator(0)

    {

        // This constuctor is only to be called internally by this class

        //   for use with its recursive member pointer (m_recursive_iterator).

        // Note: we do not need to insert m_vertex_root in the vertex - that is

        //  the responsibility of this iterator's mother, which is normally

        //  done just before calling this protected constructor.

        m_edge = m_vertex->edges_begin( m_range );

        // advance to the first good return value

        if ( !follow_edge_() &&

             m_edge != m_vertex->edges_end( m_range )) ++*this;

     }

    GenVertex::vertex_iterator::vertex_iterator( const vertex_iterator& v_iter)

        : m_vertex(0), m_visited_vertices(0), m_it_owns_set(0),

          m_recursive_iterator(0)

    {

        *this = v_iter;

    }

    GenVertex::vertex_iterator::~vertex_iterator() {

        if ( m_recursive_iterator ) delete m_recursive_iterator;

        if ( m_it_owns_set ) delete m_visited_vertices;

    }

    GenVertex::vertex_iterator& GenVertex::vertex_iterator::operator=(

        const vertex_iterator& v_iter )

    {

        // Note: when copying a vertex_iterator that is NOT the owner

        // of its set container, the pointer to the set is copied. Beware!

        // (see copy_with_own_set() if you want a different set pointed to)

        // In practise the user never needs to worry

        // since such iterators are only intended to be used internally.

        //

        // destruct data member pointers

        if ( m_recursive_iterator ) delete m_recursive_iterator;

        m_recursive_iterator = 0;

        if ( m_it_owns_set ) delete m_visited_vertices;

        m_visited_vertices = 0;

        m_it_owns_set = 0;

        // copy the target vertex_iterator to this iterator 

        m_vertex = v_iter.m_vertex;

        m_range = v_iter.m_range;

        if ( v_iter.m_it_owns_set ) {

            // i.e. this vertex will own its set if v_iter points to any

            // vertex set regardless of whether v_iter owns the set or not!

            m_visited_vertices =

                new std::set<const GenVertex*>(*v_iter.m_visited_vertices);

            m_it_owns_set = 1;

        } else {

            m_visited_vertices = v_iter.m_visited_vertices;

            m_it_owns_set = 0;

        }

        //

        // Note: m_vertex_root is already included in the set of 

        //  tv_iter.m_visited_vertices, we do not need to insert it.

        //

        m_edge = v_iter.m_edge;

        copy_recursive_iterator_( v_iter.m_recursive_iterator );

        return *this;

    }

    GenVertex* GenVertex::vertex_iterator::operator*(void) const {

        // de-reference operator

        //

        // if this iterator has an iterator_node, then we return the iterator

        // node.

        if ( m_recursive_iterator ) return  **m_recursive_iterator;

        //

        // an iterator can only return its m_vertex -- any other vertex

        //  is returned by means of a recursive iterator_node 

        //  (so this is the only place in the iterator code that a vertex 

        //   is returned!) 

        if ( m_vertex ) return m_vertex;

        return 0;

    }

    GenVertex::vertex_iterator& GenVertex::vertex_iterator::operator++(void) {

        // Pre-fix incremental operator

        //

        // check for "past the end condition" denoted by m_vertex=0

        if ( !m_vertex ) return *this;

        // if at the last edge, move into the "past the end condition"

        if ( m_edge == m_vertex->edges_end( m_range ) ) {

            m_vertex = 0;

            return *this;

        }

        // check to see if we need to create a new recursive iterator by

        // following the current edge only if a recursive iterator doesn't

        // already exist. If a new recursive_iterator is created, return it.

        if ( follow_edge_() ) {

              return *this;

        }

        //

        // if a recursive iterator already exists, increment it, and return its

        // value (unless the recursive iterator has null root_vertex [its

        // root vertex is set to null if it has already returned its root] 

        // - in which case we delete it) 

        // and return the vertex pointed to by the edge.

        if ( m_recursive_iterator ) {

            ++(*m_recursive_iterator);

            if ( **m_recursive_iterator ) {

                return *this;

            } else {

                delete m_recursive_iterator;

                m_recursive_iterator = 0;

            }

        }

        //

        // increment to the next particle edge

        ++m_edge;

        // if m_edge is at the end, then we have incremented through all

        // edges, and it is time to return m_vertex, which we accomplish

        // by returning *this

        if ( m_edge == m_vertex->edges_end( m_range ) ) return *this;

        // otherwise we follow the current edge by recursively ++ing.

        return ++(*this);

    }

    GenVertex::vertex_iterator GenVertex::vertex_iterator::operator++(int) {

        // Post-fix increment

        vertex_iterator returnvalue(*this);

        ++(*this);

        return returnvalue;

    }

    void GenVertex::vertex_iterator::copy_with_own_set(

        const vertex_iterator& v_iter,

        std::set<const GenVertex*>& visited_vertices ) {

        // intended for internal use only. (use with care!)

        // this is the same as the operator= method, but it allows the

        // user to specify which set container m_visited_vertices points to.

        // in all cases, this vertex will NOT own its set.

        //

        // destruct data member pointers

        if ( m_recursive_iterator ) delete m_recursive_iterator;

        m_recursive_iterator = 0;

        if ( m_it_owns_set ) delete m_visited_vertices;

        m_visited_vertices = 0;

        m_it_owns_set = 0;

        // copy the target vertex_iterator to this iterator 

        m_vertex = v_iter.m_vertex;

        m_range = v_iter.m_range;

        m_visited_vertices = &visited_vertices;

        m_it_owns_set = 0;

        m_edge = v_iter.m_edge;

        copy_recursive_iterator_( v_iter.m_recursive_iterator );

    }

    GenVertex* GenVertex::vertex_iterator::follow_edge_() {