[u/mrichter/AliRoot.git] / MUON / AliMUONContourMaker.cxx

/**************************************************************************
* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
*                                                                        *
* Author: The ALICE Off-line Project.                                    *
* Contributors are mentioned in the code where appropriate.              *
*                                                                        *
* Permission to use, copy, modify and distribute this software and its   *
* documentation strictly for non-commercial purposes is hereby granted   *
* without fee, provided that the above copyright notice appears in all   *
* copies and that both the copyright notice and this permission notice   *
* appear in the supporting documentation. The authors make no claims     *
* about the suitability of this software for any purpose. It is          *
* provided "as is" without express or implied warranty.                  *
**************************************************************************/

// $Id$

///
/// Maker/merger of contours. Can create contour from a set polygons or
/// merger a set of contours into a single one.
/// 
/// This is based on (one of the) algorithm found in 
/// Diane L. Souvaine and Iliana Bjorling-Sachs,
/// Proceedings of the IEEE, Vol. 80, No. 9, September 1992, p. 1449
///
/// Note that besides the AliMUON prefix, nothing is really MUON specific
/// in this class...
///
/// \author Laurent Aphecetche, Subatech
///

#include "AliMUONContourMaker.h"

#include "AliCodeTimer.h"
#include "AliLog.h"
#include "AliMUONContour.h"
#include "AliMUONPointWithRef.h"
#include "AliMUONPolygon.h"
#include "AliMUONSegment.h"
#include "AliMUONSegmentTree.h"
#include "Riostream.h"
#include "TArrayD.h"
#include "TMath.h"
#include <cassert>
#include "TArrayI.h"

/// \cond CLASSIMP
ClassImp(AliMUONContourMaker)
/// \endcond

namespace
{
  void PrintSegments(const TObjArray& array)
  {
    TIter next(&array);
    AliMUONSegment* s;
    Int_t i(0);
    while ( ( s = static_cast<AliMUONSegment*>(next()) ) )
    {
      cout << Form("i=%d %s",i,s->AsString()) << endl;
      ++i;
    }
  }
}

//_____________________________________________________________________________
AliMUONContourMaker::AliMUONContourMaker() 
{
/// Default constructor
}


//_____________________________________________________________________________
AliMUONContourMaker::~AliMUONContourMaker()
{
/// Destructor
}

//_____________________________________________________________________________
AliMUONContour* 
AliMUONContourMaker::CreateContour(const TObjArray& polygons, const char* name) const
{
  /// Create the contour of the polygon array
  /// and get back the intermediate verticals and horizontal segments
  /// both arrays are arrays of AliMUONSegment objects.
  
  AliCodeTimerAuto("",0);
  
  if ( polygons.IsEmpty() ) return 0x0; // protection against user error...
  
  // Sanity check : insure that all polygons are oriented counter-clockwise
  TIter next(&polygons);
  AliMUONPolygon* pol;
  while ( ( pol = static_cast<AliMUONPolygon*>(next()) ) )
  {
    if ( !pol->IsCounterClockwiseOriented() )
    {
      AliError(Form("Got a clockwise oriented polygon in CreateContour(%s). That's not OK !",name));
      StdoutToAliError(polygons.Print());
      return 0x0;
    }
  }
  
  AliMUONContour* contour(0x0);
  
  if ( polygons.GetLast() == 0 ) 
  {
    AliCodeTimerAuto("Trivial case",1);
    contour = new AliMUONContour(name);
    pol = static_cast<AliMUONPolygon*>(polygons.First());
    contour->Add(*pol);
    contour->AssertOrientation();
    return contour;
  }
  
  TObjArray polygonVerticalEdges; // arrray of AliMUONSegment objects  
  polygonVerticalEdges.SetOwner(kTRUE);
  // get vertical edges of input polygons
  GetVerticalEdges(polygons,polygonVerticalEdges);
  
  // sort them in ascending x order
  // if same x, insure that left edges are before right edges
  // within same x, order by increasing bottommost y (see AliMUONSegment::Compare method)
  polygonVerticalEdges.Sort();
  
  if ( polygonVerticalEdges.GetLast()+1 < 2 ) 
  {
    AliError(Form("Got too few edges here for createContour %s",name));
    polygons.Print();
    TObject* o(0x0);
    o->Print();
  }
  
  // Find the vertical edges of the merged contour. This is the meat of the algorithm...
  TObjArray contourVerticalEdges;
  contourVerticalEdges.SetOwner(kTRUE);
  Sweep(polygonVerticalEdges,contourVerticalEdges);
  
  TObjArray horizontals;
  horizontals.SetOwner(kTRUE);
  VerticalToHorizontal(contourVerticalEdges,horizontals);
  
  contour = FinalizeContour(contourVerticalEdges,horizontals);
  
  if ( contour && name ) contour->SetName(name);
  
  return contour;
}

//_____________________________________________________________________________
AliMUONContour* 
AliMUONContourMaker::FinalizeContour(const TObjArray& verticals,
                                     const TObjArray& horizontals) const
{  
  /// For a list of vertical and horizontal edges, we build the final
  /// contour object.
  
  AliCodeTimerAuto("",0);
  
  TObjArray all; // array of AliMUONSegment
  TObjArray inorder; // array of AliMUONSegment

  all.SetOwner(kFALSE);
  inorder.SetOwner(kFALSE);
    
  for ( Int_t i = 0; i <= verticals.GetLast(); ++i ) 
  {
    all.Add(verticals.UncheckedAt(i));
    all.Add(horizontals.UncheckedAt(i));
  }

  TArrayI alreadyAdded(all.GetLast()+1);
  alreadyAdded.Reset();
  
  Int_t i(0);
  
  AliMUONContour* contour = new AliMUONContour;
  
  int total(0);
  
  while ( !all.IsEmpty() )
  {
    total++;
    
    if ( total > 1000 ) 
    {
      AliError("Total 1000 reached !!!!");
      return 0x0;
    }
    
    AliMUONSegment* si = static_cast<AliMUONSegment*>(all.UncheckedAt(i));
    inorder.Add(si);
    alreadyAdded[i] = 1;
    const AliMUONSegment* all0 = static_cast<const AliMUONSegment*>(all.First());
    if ( i != 0 && AliMUONSegment::AreEqual(si->EndX(),all0->StartX()) && AliMUONSegment::AreEqual(si->EndY(),all0->StartY()) )
    {
      Int_t n(-1);
      
      AliMUONPolygon polygon(inorder.GetLast()+2);

      // we got a cycle. Add it to the contour
      for ( Int_t j = 0; j <= inorder.GetLast(); ++j ) 
      {
        AliMUONSegment* s = static_cast<AliMUONSegment*>(inorder.UncheckedAt(j));
        polygon.SetVertex(++n,s->StartX(),s->StartY());
        all.Remove(s);
      }
      
      all.Compress();

      polygon.Close();
      
      contour->Add(polygon);
      
      if ( ! all.IsEmpty() )
      {
        i = 0;
        inorder.Clear();
        alreadyAdded.Set(all.GetLast()+1);
        alreadyAdded.Reset();
      }
      continue;
    }
    
    for ( Int_t j = 0; j <= all.GetLast(); ++j) 
    {
      if ( j != i && alreadyAdded[j] == 0 ) 
      {        
        const AliMUONSegment* sj = static_cast<const AliMUONSegment*>(all.UncheckedAt(j));
        if ( AliMUONSegment::AreEqual(si->EndX(),sj->StartX()) && AliMUONSegment::AreEqual(si->EndY(),sj->StartY()))
        {
          i = j;
          break;
        }
      }
    }    
  }
  
  contour->AssertOrientation(kTRUE);
  return contour;
}


//_____________________________________________________________________________
void 
AliMUONContourMaker::GetVerticalEdges(const TObjArray& polygons, TObjArray& polygonVerticalEdges) const
{
  /// From an array of polygons, extract the list of vertical edges.
  /// Output array polygonVerticalEdges should be empty before calling.
  
  AliCodeTimerAuto("",0);
  
  for ( Int_t i = 0; i <= polygons.GetLast(); ++i ) 
  {
    const AliMUONPolygon* g = static_cast<const AliMUONPolygon*>(polygons.UncheckedAt(i));
    for ( Int_t j = 0; j < g->NumberOfVertices()-1; ++j ) 
    {
      if ( AliMUONSegment::AreEqual(g->X(j),g->X(j+1)) ) // segment is vertical
      {
        polygonVerticalEdges.Add(new AliMUONSegment(g->X(j),g->Y(j),g->X(j+1),g->Y(j+1)));
      }
    }
  }
}


//_____________________________________________________________________________
void
AliMUONContourMaker::GetYPositions(const TObjArray& polygonVerticalEdges,
                                   TArrayD& yPositions) const
{
  /// Fill the array yPositions with the different y positions found in 
  /// polygonVerticalEdges
  
  AliCodeTimerAuto("",0);
  
  Double_t* y = new Double_t[polygonVerticalEdges.GetSize()*2];
  Int_t n(0);
  
  for ( Int_t i = 0; i < polygonVerticalEdges.GetLast(); ++i ) 
  {
    AliMUONSegment* s = static_cast<AliMUONSegment*>(polygonVerticalEdges.UncheckedAt(i));
    y[n] = s->StartY();
    y[n+1] = s->EndY();
    n += 2;
  }
  Int_t* ix = new Int_t[n+1];
  
  TMath::Sort(n,y,ix,kFALSE);
  
  yPositions.Set(n+1);
  
  Int_t u(0);
  Double_t x(FLT_MAX);
  
  for ( Int_t i = 0; i < n; ++i ) 
  {
    if ( y[ix[i]] != x )
    {
      yPositions[u] = y[ix[i]];
      x = y[ix[i]];
      ++u;
    }
  }

  yPositions.Set(u);
  
  delete[] ix;
  delete[] y;
  
}

//_____________________________________________________________________________
AliMUONContour* 
AliMUONContourMaker::MergeContour(const TObjArray& contours, const char* name) const
{
  /// Merge all the polygons of all contours into a single contour
  
  AliCodeTimerAuto("",0);
  
  TObjArray polygons;
  polygons.SetOwner(kTRUE);
  
  TIter next(&contours);
  AliMUONContour* contour;
  while ( ( contour = static_cast<AliMUONContour*>(next()) ) )
  {
    const TObjArray* contourPolygons = contour->Polygons();
    TIter nextPol(contourPolygons);
    AliMUONPolygon* pol;
    while ( ( pol = static_cast<AliMUONPolygon*>(nextPol()) ) )
    {
      polygons.Add(new AliMUONPolygon(*pol));
    }
  }
  
  if ( polygons.IsEmpty() ) return 0x0;
  
  contour = CreateContour(polygons,name);
  
  return contour;
}

//_____________________________________________________________________________
void 
AliMUONContourMaker::SortPoints(const TObjArray& polygonVerticalEdges, 
                                TObjArray& sortedPoints) const
{
  /// Sort the point of the vertical edges in ascending order, first on ordinate, 
  /// then on abcissa, and put them in output vector sortedPoints.
  /// Output array sortedPoints should be empty before calling this method.
  
  AliCodeTimerAuto("",0);
  
  for ( Int_t i = 0; i <= polygonVerticalEdges.GetLast(); ++i )
  {
    const AliMUONSegment* e = static_cast<const AliMUONSegment*>(polygonVerticalEdges.UncheckedAt(i));
    sortedPoints.Add(new AliMUONPointWithRef(e->StartX(),e->StartY(),i));
    sortedPoints.Add(new AliMUONPointWithRef(e->EndX(),e->EndY(),i));
    // note that we keep track of the original edge, which is used
    // later on to deduce orientation of horizontal edges.
  }
  
  sortedPoints.Sort(); 
}

//_____________________________________________________________________________
void 
AliMUONContourMaker::Sweep(const TObjArray& polygonVerticalEdges, 
                           TObjArray& contourVerticalEdges) const
{
  /// This is the meat of the algorithm of the contour merging...
  
  AliCodeTimerAuto("",0);
  
  TArrayD yPositions;
  GetYPositions(polygonVerticalEdges,yPositions);
  
  AliMUONSegmentTree segmentTree(yPositions);
  
  for ( Int_t i = 0; i <= polygonVerticalEdges.GetLast(); ++i )
  {
    const AliMUONSegment* edge = static_cast<const AliMUONSegment*>(polygonVerticalEdges.UncheckedAt(i));
    
    assert(edge!=0x0);
    
    if ( edge->IsLeftEdge() ) 
    {
      segmentTree.Contribution(edge->Bottom(),edge->Top());
      segmentTree.InsertInterval(edge->Bottom(),edge->Top());
    }
    else
    {
      segmentTree.DeleteInterval(edge->Bottom(),edge->Top());
      segmentTree.Contribution(edge->Bottom(),edge->Top());
    }
    
    AliMUONSegment e1(*edge);
    
    if ( i < polygonVerticalEdges.GetLast() ) 
    {
      const AliMUONSegment* next = static_cast<const AliMUONSegment*>(polygonVerticalEdges.UncheckedAt(i+1));
      e1 = *next;
    }

    if ( ( edge->IsLeftEdge() != e1.IsLeftEdge() ) ||
        ( !AliMUONSegment::AreEqual(edge->StartX(),e1.StartX() ) ) ||
        ( i == polygonVerticalEdges.GetLast() ) )
    {
      const TObjArray& stack = segmentTree.Stack();
      
      double x = edge->StartX();
      
      for ( Int_t j = 0; j <= stack.GetLast(); ++j )
      {
        AliMUONSegment* sj = static_cast<AliMUONSegment*>(stack.UncheckedAt(j));
        AliMUONSegment* s = new AliMUONSegment(x,sj->StartY(),x,sj->EndY());
        
        if  (s->IsAPoint()) 
        {
          delete s;
          continue;
        }
        
        if ( edge->IsLeftEdge() != s->IsLeftEdge() ) 
        {
          s->Set(x,sj->EndY(),x,sj->StartY());
        }
        contourVerticalEdges.Add(s);
      }
      segmentTree.ResetStack();
    }
  }
}

//_____________________________________________________________________________
void 
AliMUONContourMaker::VerticalToHorizontal(const TObjArray& polygonVerticalEdges,
                                          TObjArray& horizontalEdges) const
{
  /// Deduce the set of horizontal edges from the vertical edges
  /// Output array horizontalEdges should be empty before calling this method
  
  AliCodeTimerAuto("",0);
  
  TObjArray points; // array of AliMUONPointWithRef
  points.SetOwner(kTRUE);
  
  SortPoints(polygonVerticalEdges,points);
  
  for ( Int_t k = 0; k < (points.GetLast()+1)/2; ++k )
  {
    const AliMUONPointWithRef* p1 = static_cast<AliMUONPointWithRef*>(points.UncheckedAt(k*2));
    const AliMUONPointWithRef* p2 = static_cast<AliMUONPointWithRef*>(points.UncheckedAt(k*2+1));
    
    const AliMUONSegment* refEdge = static_cast<const AliMUONSegment*>(polygonVerticalEdges.UncheckedAt(p1->Ref()));
    
    // (p1,p2) is the horizontal edge.
    // refEdge is used to deduce the orientation of (p1,p2)
    
    if ( AliMUONSegment::AreEqual(p1->X(),refEdge->EndX()) && AliMUONSegment::AreEqual(p1->Y(),refEdge->EndY()) )
//    if ( AreEqual(p1,refEdge->End()) )
    {
      horizontalEdges.Add(new AliMUONSegment(p1->X(),p1->Y(),p2->X(),p2->Y()));
    }
    else
    {
      horizontalEdges.Add(new AliMUONSegment(p2->X(),p2->Y(),p1->X(),p1->Y()));
    }
  }
}
Commit	Line	Data
0b936dc0	1	/**************************************************************************
	2	* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
	3	* *
	4	* Author: The ALICE Off-line Project. *
	5	* Contributors are mentioned in the code where appropriate. *
	6	* *
	7	* Permission to use, copy, modify and distribute this software and its *
	8	* documentation strictly for non-commercial purposes is hereby granted *
	9	* without fee, provided that the above copyright notice appears in all *
	10	* copies and that both the copyright notice and this permission notice *
	11	* appear in the supporting documentation. The authors make no claims *
	12	* about the suitability of this software for any purpose. It is *
	13	* provided "as is" without express or implied warranty. *
	14	**************************************************************************/
	15
	16	// $Id$
	17
	18	///
	19	/// Maker/merger of contours. Can create contour from a set polygons or
	20	/// merger a set of contours into a single one.
	21	///
	22	/// This is based on (one of the) algorithm found in
	23	/// Diane L. Souvaine and Iliana Bjorling-Sachs,
	24	/// Proceedings of the IEEE, Vol. 80, No. 9, September 1992, p. 1449
	25	///
	26	/// Note that besides the AliMUON prefix, nothing is really MUON specific
	27	/// in this class...
	28	///
	29	/// \author Laurent Aphecetche, Subatech
	30	///
	31
	32	#include "AliMUONContourMaker.h"
	33
	34	#include "AliCodeTimer.h"
	35	#include "AliLog.h"
	36	#include "AliMUONContour.h"
	37	#include "AliMUONPointWithRef.h"
	38	#include "AliMUONPolygon.h"
	39	#include "AliMUONSegment.h"
	40	#include "AliMUONSegmentTree.h"
	41	#include "Riostream.h"
	42	#include "TArrayD.h"
	43	#include "TMath.h"
	44	#include <cassert>
fef32488	45	#include "TArrayI.h"
0b936dc0	46
	47	/// \cond CLASSIMP
	48	ClassImp(AliMUONContourMaker)
	49	/// \endcond
	50
	51	namespace
	52	{
	53	void PrintSegments(const TObjArray& array)
	54	{
	55	TIter next(&array);
	56	AliMUONSegment* s;
	57	Int_t i(0);
	58	while ( ( s = static_cast<AliMUONSegment*>(next()) ) )
	59	{
	60	cout << Form("i=%d %s",i,s->AsString()) << endl;
	61	++i;
	62	}
	63	}
	64	}
	65
	66	//_____________________________________________________________________________
	67	AliMUONContourMaker::AliMUONContourMaker()
	68	{
cddcc1f3	69	/// Default constructor
0b936dc0	70	}
0b936dc0	71
cddcc1f3	72
0b936dc0	73	//_____________________________________________________________________________
	74	AliMUONContourMaker::~AliMUONContourMaker()
	75	{
cddcc1f3	76	/// Destructor
0b936dc0	77	}
	78
	79	//_____________________________________________________________________________
	80	AliMUONContour*
	81	AliMUONContourMaker::CreateContour(const TObjArray& polygons, const char* name) const
	82	{
	83	/// Create the contour of the polygon array
	84	/// and get back the intermediate verticals and horizontal segments
	85	/// both arrays are arrays of AliMUONSegment objects.
	86
99c136e1	87	AliCodeTimerAuto("",0);
0b936dc0	88
	89	if ( polygons.IsEmpty() ) return 0x0; // protection against user error...
	90
	91	// Sanity check : insure that all polygons are oriented counter-clockwise
	92	TIter next(&polygons);
	93	AliMUONPolygon* pol;
	94	while ( ( pol = static_cast<AliMUONPolygon*>(next()) ) )
	95	{
	96	if ( !pol->IsCounterClockwiseOriented() )
	97	{
	98	AliError(Form("Got a clockwise oriented polygon in CreateContour(%s). That's not OK !",name));
	99	StdoutToAliError(polygons.Print());
	100	return 0x0;
	101	}
	102	}
	103
	104	AliMUONContour* contour(0x0);
	105
	106	if ( polygons.GetLast() == 0 )
	107	{
dd0be8a7	108	AliCodeTimerAuto("Trivial case",1);
0b936dc0	109	contour = new AliMUONContour(name);
7da191b6	110	pol = static_cast<AliMUONPolygon*>(polygons.First());
0b936dc0	111	contour->Add(*pol);
	112	contour->AssertOrientation();
	113	return contour;
	114	}
	115
	116	TObjArray polygonVerticalEdges; // arrray of AliMUONSegment objects
	117	polygonVerticalEdges.SetOwner(kTRUE);
	118	// get vertical edges of input polygons
	119	GetVerticalEdges(polygons,polygonVerticalEdges);
	120
	121	// sort them in ascending x order
	122	// if same x, insure that left edges are before right edges
	123	// within same x, order by increasing bottommost y (see AliMUONSegment::Compare method)
	124	polygonVerticalEdges.Sort();
	125
	126	if ( polygonVerticalEdges.GetLast()+1 < 2 )
	127	{
	128	AliError(Form("Got too few edges here for createContour %s",name));
	129	polygons.Print();
	130	TObject* o(0x0);
	131	o->Print();
	132	}
	133
	134	// Find the vertical edges of the merged contour. This is the meat of the algorithm...
	135	TObjArray contourVerticalEdges;
	136	contourVerticalEdges.SetOwner(kTRUE);
	137	Sweep(polygonVerticalEdges,contourVerticalEdges);
	138
	139	TObjArray horizontals;
	140	horizontals.SetOwner(kTRUE);
	141	VerticalToHorizontal(contourVerticalEdges,horizontals);
	142
	143	contour = FinalizeContour(contourVerticalEdges,horizontals);
	144
	145	if ( contour && name ) contour->SetName(name);
	146
	147	return contour;
	148	}
	149
	150	//_____________________________________________________________________________
	151	AliMUONContour*
	152	AliMUONContourMaker::FinalizeContour(const TObjArray& verticals,
	153	const TObjArray& horizontals) const
	154	{
	155	/// For a list of vertical and horizontal edges, we build the final
	156	/// contour object.
	157
99c136e1	158	AliCodeTimerAuto("",0);
0b936dc0	159
	160	TObjArray all; // array of AliMUONSegment
	161	TObjArray inorder; // array of AliMUONSegment
	162
	163	all.SetOwner(kFALSE);
	164	inorder.SetOwner(kFALSE);
	165
	166	for ( Int_t i = 0; i <= verticals.GetLast(); ++i )
	167	{
	168	all.Add(verticals.UncheckedAt(i));
	169	all.Add(horizontals.UncheckedAt(i));
	170	}
fef32488	171
	172	TArrayI alreadyAdded(all.GetLast()+1);
	173	alreadyAdded.Reset();
0b936dc0	174
	175	Int_t i(0);
	176
	177	AliMUONContour* contour = new AliMUONContour;
	178
	179	int total(0);
	180
	181	while ( !all.IsEmpty() )
	182	{
	183	total++;
	184
	185	if ( total > 1000 )
	186	{
	187	AliError("Total 1000 reached !!!!");
	188	return 0x0;
	189	}
	190
	191	AliMUONSegment* si = static_cast<AliMUONSegment*>(all.UncheckedAt(i));
	192	inorder.Add(si);
fef32488	193	alreadyAdded[i] = 1;
0b936dc0	194	const AliMUONSegment* all0 = static_cast<const AliMUONSegment*>(all.First());
	195	if ( i != 0 && AliMUONSegment::AreEqual(si->EndX(),all0->StartX()) && AliMUONSegment::AreEqual(si->EndY(),all0->StartY()) )
	196	{
	197	Int_t n(-1);
	198
	199	AliMUONPolygon polygon(inorder.GetLast()+2);
	200
	201	// we got a cycle. Add it to the contour
	202	for ( Int_t j = 0; j <= inorder.GetLast(); ++j )
	203	{
	204	AliMUONSegment* s = static_cast<AliMUONSegment*>(inorder.UncheckedAt(j));
	205	polygon.SetVertex(++n,s->StartX(),s->StartY());
	206	all.Remove(s);
	207	}
	208
	209	all.Compress();
	210
	211	polygon.Close();
	212
	213	contour->Add(polygon);
	214
	215	if ( ! all.IsEmpty() )
	216	{
	217	i = 0;
	218	inorder.Clear();
fef32488	219	alreadyAdded.Set(all.GetLast()+1);
fef32488	220	alreadyAdded.Reset();
0b936dc0	221	}
	222	continue;
	223	}
	224
	225	for ( Int_t j = 0; j <= all.GetLast(); ++j)
	226	{
fef32488	227	if ( j != i && alreadyAdded[j] == 0 )
0b936dc0	228	{
	229	const AliMUONSegment* sj = static_cast<const AliMUONSegment*>(all.UncheckedAt(j));
	230	if ( AliMUONSegment::AreEqual(si->EndX(),sj->StartX()) && AliMUONSegment::AreEqual(si->EndY(),sj->StartY()))
	231	{
	232	i = j;
	233	break;
	234	}
	235	}
	236	}
	237	}
	238
	239	contour->AssertOrientation(kTRUE);
	240	return contour;
	241	}
	242
	243
	244	//_____________________________________________________________________________
	245	void
	246	AliMUONContourMaker::GetVerticalEdges(const TObjArray& polygons, TObjArray& polygonVerticalEdges) const
	247	{
	248	/// From an array of polygons, extract the list of vertical edges.
	249	/// Output array polygonVerticalEdges should be empty before calling.
	250
99c136e1	251	AliCodeTimerAuto("",0);
0b936dc0	252
	253	for ( Int_t i = 0; i <= polygons.GetLast(); ++i )
	254	{
	255	const AliMUONPolygon* g = static_cast<const AliMUONPolygon*>(polygons.UncheckedAt(i));
	256	for ( Int_t j = 0; j < g->NumberOfVertices()-1; ++j )
	257	{
	258	if ( AliMUONSegment::AreEqual(g->X(j),g->X(j+1)) ) // segment is vertical
	259	{
	260	polygonVerticalEdges.Add(new AliMUONSegment(g->X(j),g->Y(j),g->X(j+1),g->Y(j+1)));
	261	}
	262	}
	263	}
	264	}
	265
	266
	267	//_____________________________________________________________________________
	268	void
	269	AliMUONContourMaker::GetYPositions(const TObjArray& polygonVerticalEdges,
	270	TArrayD& yPositions) const
	271	{
	272	/// Fill the array yPositions with the different y positions found in
	273	/// polygonVerticalEdges
	274
99c136e1	275	AliCodeTimerAuto("",0);
0b936dc0	276
	277	Double_t* y = new Double_t[polygonVerticalEdges.GetSize()*2];
	278	Int_t n(0);
	279
	280	for ( Int_t i = 0; i < polygonVerticalEdges.GetLast(); ++i )
	281	{
	282	AliMUONSegment* s = static_cast<AliMUONSegment*>(polygonVerticalEdges.UncheckedAt(i));
	283	y[n] = s->StartY();
	284	y[n+1] = s->EndY();
	285	n += 2;
	286	}
	287	Int_t* ix = new Int_t[n+1];
	288
	289	TMath::Sort(n,y,ix,kFALSE);
	290
	291	yPositions.Set(n+1);
	292
	293	Int_t u(0);
	294	Double_t x(FLT_MAX);
	295
	296	for ( Int_t i = 0; i < n; ++i )
	297	{
	298	if ( y[ix[i]] != x )
	299	{
	300	yPositions[u] = y[ix[i]];
	301	x = y[ix[i]];
	302	++u;
	303	}
	304	}
	305
	306	yPositions.Set(u);
	307
	308	delete[] ix;
	309	delete[] y;
	310
	311	}
	312
	313	//_____________________________________________________________________________
	314	AliMUONContour*
	315	AliMUONContourMaker::MergeContour(const TObjArray& contours, const char* name) const
	316	{
	317	/// Merge all the polygons of all contours into a single contour
	318
99c136e1	319	AliCodeTimerAuto("",0);
0b936dc0	320
	321	TObjArray polygons;
	322	polygons.SetOwner(kTRUE);
	323
	324	TIter next(&contours);
	325	AliMUONContour* contour;
	326	while ( ( contour = static_cast<AliMUONContour*>(next()) ) )
	327	{
	328	const TObjArray* contourPolygons = contour->Polygons();
	329	TIter nextPol(contourPolygons);
	330	AliMUONPolygon* pol;
	331	while ( ( pol = static_cast<AliMUONPolygon*>(nextPol()) ) )
	332	{
	333	polygons.Add(new AliMUONPolygon(*pol));
	334	}
	335	}
	336
	337	if ( polygons.IsEmpty() ) return 0x0;
	338
	339	contour = CreateContour(polygons,name);
	340
	341	return contour;
	342	}
	343
	344	//_____________________________________________________________________________
	345	void
	346	AliMUONContourMaker::SortPoints(const TObjArray& polygonVerticalEdges,
	347	TObjArray& sortedPoints) const
	348	{
	349	/// Sort the point of the vertical edges in ascending order, first on ordinate,
	350	/// then on abcissa, and put them in output vector sortedPoints.
	351	/// Output array sortedPoints should be empty before calling this method.
	352
99c136e1	353	AliCodeTimerAuto("",0);
0b936dc0	354
	355	for ( Int_t i = 0; i <= polygonVerticalEdges.GetLast(); ++i )
	356	{
	357	const AliMUONSegment* e = static_cast<const AliMUONSegment*>(polygonVerticalEdges.UncheckedAt(i));
	358	sortedPoints.Add(new AliMUONPointWithRef(e->StartX(),e->StartY(),i));
	359	sortedPoints.Add(new AliMUONPointWithRef(e->EndX(),e->EndY(),i));
	360	// note that we keep track of the original edge, which is used
	361	// later on to deduce orientation of horizontal edges.
	362	}
	363
	364	sortedPoints.Sort();
	365	}
	366
	367	//_____________________________________________________________________________
	368	void
	369	AliMUONContourMaker::Sweep(const TObjArray& polygonVerticalEdges,
	370	TObjArray& contourVerticalEdges) const
	371	{
	372	/// This is the meat of the algorithm of the contour merging...
	373
99c136e1	374	AliCodeTimerAuto("",0);
0b936dc0	375
	376	TArrayD yPositions;
	377	GetYPositions(polygonVerticalEdges,yPositions);
	378
	379	AliMUONSegmentTree segmentTree(yPositions);
	380
	381	for ( Int_t i = 0; i <= polygonVerticalEdges.GetLast(); ++i )
	382	{
	383	const AliMUONSegment* edge = static_cast<const AliMUONSegment*>(polygonVerticalEdges.UncheckedAt(i));
	384
	385	assert(edge!=0x0);
	386
	387	if ( edge->IsLeftEdge() )
	388	{
	389	segmentTree.Contribution(edge->Bottom(),edge->Top());
	390	segmentTree.InsertInterval(edge->Bottom(),edge->Top());
	391	}
	392	else
	393	{
	394	segmentTree.DeleteInterval(edge->Bottom(),edge->Top());
	395	segmentTree.Contribution(edge->Bottom(),edge->Top());
	396	}
	397
	398	AliMUONSegment e1(*edge);
	399
	400	if ( i < polygonVerticalEdges.GetLast() )
	401	{
	402	const AliMUONSegment* next = static_cast<const AliMUONSegment*>(polygonVerticalEdges.UncheckedAt(i+1));
	403	e1 = *next;
	404	}
	405
	406	if ( ( edge->IsLeftEdge() != e1.IsLeftEdge() ) \|\|
	407	( !AliMUONSegment::AreEqual(edge->StartX(),e1.StartX() ) ) \|\|
	408	( i == polygonVerticalEdges.GetLast() ) )
	409	{
	410	const TObjArray& stack = segmentTree.Stack();
	411
	412	double x = edge->StartX();
	413
	414	for ( Int_t j = 0; j <= stack.GetLast(); ++j )
	415	{
	416	AliMUONSegment* sj = static_cast<AliMUONSegment*>(stack.UncheckedAt(j));
	417	AliMUONSegment* s = new AliMUONSegment(x,sj->StartY(),x,sj->EndY());
	418
	419	if (s->IsAPoint())
	420	{
	421	delete s;
	422	continue;
	423	}
	424
	425	if ( edge->IsLeftEdge() != s->IsLeftEdge() )
	426	{
	427	s->Set(x,sj->EndY(),x,sj->StartY());
	428	}
	429	contourVerticalEdges.Add(s);
	430	}
	431	segmentTree.ResetStack();
	432	}
	433	}
	434	}
	435
	436	//_____________________________________________________________________________
	437	void
	438	AliMUONContourMaker::VerticalToHorizontal(const TObjArray& polygonVerticalEdges,
439	TObjArray& horizontalEdges) const
440	{
441	/// Deduce the set of horizontal edges from the vertical edges
442	/// Output array horizontalEdges should be empty before calling this method
443
99c136e1	444	AliCodeTimerAuto("",0);
0b936dc0	445
	446	TObjArray points; // array of AliMUONPointWithRef
	447	points.SetOwner(kTRUE);
	448
	449	SortPoints(polygonVerticalEdges,points);
	450
	451	for ( Int_t k = 0; k < (points.GetLast()+1)/2; ++k )
	452	{
	453	const AliMUONPointWithRef* p1 = static_cast<AliMUONPointWithRef>(points.UncheckedAt(k2));
	454	const AliMUONPointWithRef* p2 = static_cast<AliMUONPointWithRef>(points.UncheckedAt(k2+1));
	455
	456	const AliMUONSegment* refEdge = static_cast<const AliMUONSegment*>(polygonVerticalEdges.UncheckedAt(p1->Ref()));
	457
	458	// (p1,p2) is the horizontal edge.
	459	// refEdge is used to deduce the orientation of (p1,p2)
	460
	461	if ( AliMUONSegment::AreEqual(p1->X(),refEdge->EndX()) && AliMUONSegment::AreEqual(p1->Y(),refEdge->EndY()) )
	462	// if ( AreEqual(p1,refEdge->End()) )
	463	{
	464	horizontalEdges.Add(new AliMUONSegment(p1->X(),p1->Y(),p2->X(),p2->Y()));
	465	}
	466	else
	467	{
	468	horizontalEdges.Add(new AliMUONSegment(p2->X(),p2->Y(),p1->X(),p1->Y()));
	469	}
	470	}
	471	}
	472