source: trunk/MagicSoft/Mars/manalysis/MCerPhotEvt.cc@ 646

#include "MCerPhotEvt.h"

#include <math.h>
#include <TClonesArray.h>
#include <TCanvas.h>

#include "MCamGeom.h"
#include "MCamNeighbor.h"
#include "MCamDisplay.h"
#include "MHexagon.h"

ClassImp(MCerPhotPix)
ClassImp(MCerPhotEvt)

MCerPhotPix::MCerPhotPix(Int_t pix, Float_t phot  , Float_t errphot ) 
{ 
  //  default constructor
  fPixId    = pix ; 
  fIsUsed   = kTRUE ; 
  fIsCore   = kFALSE ; 
  fPhot     = phot ; 
  fErrPhot  = errphot ; 
} 

void MCerPhotPix::SetPixelContent(Int_t pix, Float_t phot  , Float_t errphot)
{
  fPixId    = pix ; 
  fIsUsed   = kTRUE ; 
  fIsUsed   = kFALSE ; 
  fPhot     = phot ; 
  fErrPhot  = errphot ; 
}

void MCerPhotPix::Print() 
{ 
  //   information about a pixel
  cout << "MCerPhotPix: Pixel: "<< fPixId ; 
  
  if ( fIsUsed == kTRUE )
    cout << "   Used " ; 
  else 
    if ( fIsUsed == kFALSE )
      cout << " UnUsed " ;

  if ( fIsCore == kTRUE )
    cout << "   Core " ; 
  else 
    if ( fIsCore == kFALSE )
      cout << "        " ;

  cout << "  Nphot= " << fPhot
       << "  Error(Nphot) = " << fErrPhot
       << endl ; 
}

// ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
// ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
// ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++


MCerPhotEvt::MCerPhotEvt(const char *name, const char *title ) 
{
  //   the default constructor 

  
  *fName  = name  ? name  : "MCerPhotEvt";
  *fTitle = name  ? name  : "(Number of Photon)-Event Information";
  
  fType = 0 ; 
  fNbPixels = 0 ; 
  
  fPixels = new TClonesArray ("MCerPhotPix", 577) ;
  
  fNN = new MCamNeighbor() ; 

  fPixels->Clear() ; 
}

void MCerPhotEvt::Draw(Option_t* option) 
{
  //   FIXME!!!
  // 
  
  MCamDisplay *disp = new MCamDisplay(fType)  ; 
  
  disp->Draw( this ) ; 
  
  //  disp->Draw() ; 

}




Int_t MCerPhotEvt::GetNbPixels() 
{
  return fNbPixels ; 
} 

void MCerPhotEvt::AddPixel(Int_t id, Float_t nph, Float_t err)
{
  TClonesArray &caP = *fPixels ;
  new ( caP[fNbPixels++] ) MCerPhotPix( id, nph, err ) ;
}

void MCerPhotEvt::Clear() 
{
  fNbPixels = 0 ; 
  fPixels->Clear() ; 
}

void MCerPhotEvt::Print() 
{
  cout << "MCerPhotEvt::Print()"  
       << "Number of Pixels: " << fNbPixels
       << "(" << fPixels->GetEntries() << ")" 
       << endl ; 

  for (Int_t il=0; il<fPixels->GetEntries(); il++ ) 
    {
      ((MCerPhotPix *) fPixels->At(il))->Print() ; 
    }
}

void MCerPhotEvt::CleanLevel1()
{ 
  //  This method looks for all pixels with an entry (photons) 
  //  that is three times bigger than the noise of the pixel
  
  Float_t entry, noise ; 
  
  // first look for pixels above some noise level 
  
  for (Int_t il=0; il<fPixels->GetEntries(); il++ ) 
    {
      entry = ((MCerPhotPix *) fPixels->At(il))->GetPhotons() ; 
      noise = ((MCerPhotPix *) fPixels->At(il))->GetErrorPhot() ; 
      
      if (entry < 3 * noise ) 
        ((MCerPhotPix *) fPixels->At(il))->SetPixelUnused() ; 
        
    }
} 

void MCerPhotEvt::CleanLevel2()
{
  //  check if the  survived pixel have a neighbor, that also 
  //  survived
  
  Int_t id, id2 ; 
  Int_t itest ;
 
  for (Int_t il=0; il<fPixels->GetEntries(); il++ ) 
    {
      if ( ((MCerPhotPix *) fPixels->At(il))->IsPixelUsed() == kTRUE )
        { 
          id = ((MCerPhotPix *) fPixels->At(il))->GetPixId() ; 
          
          itest = 0 ; 
          for (Int_t in=0 ; in < 6 ; in++ ) { 
             
            id2 = fNN->GetNN(id, in ) ; 
            
            if (id2 >=0 ) { 
              if ( PixelIsUsed(id2) == kTRUE ) 
                itest++ ; 
            }   
          } 
          
          // 
          //   check if no neighbor, then set unused 
          //
          if ( itest == 0 ) 
            { 
              ((MCerPhotPix *) fPixels->At(il))->SetPixelUnused() ; 
            }

        } 
    } 
  
  // now we declare all pixels that survive as CorePixels
  
  for (Int_t il=0; il<fPixels->GetEntries(); il++ ) 
    {
      if ( ((MCerPhotPix *) fPixels->At(il))->IsPixelUsed() == kTRUE )
        { 
          ((MCerPhotPix *) fPixels->At(il))->SetCorePixel() ; 
        }
    }

} 

void MCerPhotEvt::CleanLevel3()
{
  //   Look for the boundary pixels around the core pixels
  //   if a pixel has more than 2.5 sigma, and a core neigbor
  //   it is declared as used. 
  
  Float_t entry, noise ; 
  Int_t   id, id2 ; 
  for (Int_t il=0; il<fPixels->GetEntries(); il++ ) 
    {
      if ( ((MCerPhotPix *) fPixels->At(il))->IsCorePixel() == kFALSE )
        {  
          entry = ((MCerPhotPix *) fPixels->At(il))->GetPhotons() ; 
          noise = ((MCerPhotPix *) fPixels->At(il))->GetErrorPhot() ; 
      
          if (entry > 2.5 * noise ) { 
            id = ((MCerPhotPix *) fPixels->At(il))->GetPixId()  ; 
            for (Int_t in=0 ; in < 6 ; in++ )
              { 
                id2 = fNN->GetNN(id, in ) ; 
                if (id2 >=0 ) 
                  { 
                    if ( PixelIsCore(id2) == kTRUE ) { 
                      ((MCerPhotPix *) fPixels->At(il))->SetPixelUsed() ; 
                      break ; 
                    } 
                  } 
              } 
          } 
        } 
    } 
} 




Bool_t MCerPhotEvt::PixelExist(Int_t id ) 
{ 
  //   Checks if in the pixel list is an entry with pixel id
  
  for (Int_t il=0; il<fPixels->GetEntries(); il++ ) 
    {
      if ( id == ((MCerPhotPix *) fPixels->At(il))->GetPixId() ) { 
        
        // cout << " PixelExist " << il ; 
        return kTRUE ; 
      } 
    } 

  return kFALSE ; 

} 

Bool_t MCerPhotEvt::PixelIsUsed(Int_t id ) 
{ 
  //   Checks if in the pixel list is an entry with pixel id
  
  for (Int_t il=0; il<fPixels->GetEntries(); il++ ) 
    {
      if ( id == ((MCerPhotPix *) fPixels->At(il))->GetPixId() && 
           ((MCerPhotPix *) fPixels->At(il))->IsPixelUsed() == kTRUE ) { 
        
        // cout << " PixelIsUsed  " << il ; 
        return kTRUE ; 
      } 
    } 

  return kFALSE ; 

} 

Bool_t MCerPhotEvt::PixelIsCore(Int_t id ) 
{ 
  //   Checks if in the pixel list is an entry with pixel id
  
  for (Int_t il=0; il<fPixels->GetEntries(); il++ ) 
    {
      if ( id == ((MCerPhotPix *) fPixels->At(il))->GetPixId() && 
           ((MCerPhotPix *) fPixels->At(il))->IsCorePixel() == kTRUE ) { 
        
        return kTRUE ; 
      } 
    } 

  return kFALSE ; 

} 

Int_t MCerPhotEvt::GetPixelId(Int_t i ) 
{ 
  return ( ( (MCerPhotPix *) fPixels->At(i))->GetPixId() ) ; 
} 

Bool_t  MCerPhotEvt::IsPixelUsed(Int_t i ) 
{
  return ( ( (MCerPhotPix *) fPixels->At(i))->IsPixelUsed() ) ;
}
 
Float_t MCerPhotEvt::GetPhotons(Int_t i ) 
{ 
  return ( ( (MCerPhotPix *) fPixels->At(i))->GetPhotons() ) ; 
}
 
Float_t MCerPhotEvt::GetErrorPhot(Int_t i ) 
{ 
  return ( ( (MCerPhotPix *) fPixels->At(i))->GetErrorPhot() ) ; 
} 


Float_t MCerPhotEvt::GetMinimumPhoton() 
{ 
  if ( fNbPixels <= 0 ) 
    return -5. ; 
  
  Float_t minWert ; 
  minWert = ((MCerPhotPix *) fPixels->At(0))->GetPhotons() ; 

  Float_t testWert ; 

  for ( Int_t i =0 ; i<fNbPixels ; i++ ) { 
    testWert = ((MCerPhotPix *) fPixels->At(i))->GetPhotons() ; 
    
    if ( minWert >= testWert ) 
      minWert = testWert ; 
  }

  return  minWert ;
} 

Float_t MCerPhotEvt::GetMaximumPhoton() 
{ 
  if ( fNbPixels <= 0 ) 
    return 50. ; 
  
  Float_t maxWert ; 
  maxWert = ((MCerPhotPix *) fPixels->At(0))->GetPhotons() ; 

  Float_t testWert ; 

  for ( Int_t i =0 ; i<fNbPixels ; i++ ) { 
    testWert = ((MCerPhotPix *) fPixels->At(i))->GetPhotons() ; 
    
    if ( maxWert <= testWert ) 
      maxWert = testWert ; 
  }

  return  maxWert ;
}