Index: trunk/MagicSoft/Mars/mtemp/mifae/library/MIslandCalc.cc
===================================================================
--- trunk/MagicSoft/Mars/mtemp/mifae/library/MIslandCalc.cc	(revision 4222)
+++ trunk/MagicSoft/Mars/mtemp/mifae/library/MIslandCalc.cc	(revision 4286)
@@ -4,5 +4,5 @@
 !   Author(s): Ester Aliu, 2/2004 <aliu@ifae.es>
 !   
-!   Last Update: 5/2004
+!   Last Update: 6/2004
 !
 \* ======================================================================== */
@@ -124,5 +124,26 @@
 Int_t MIslandCalc::Process()
 {
-
+  
+  if (fIslandAlgorithm == 1)
+    Calc1();
+  
+  if (fIslandAlgorithm == 2)
+    Calc2();
+  return kTRUE;  
+   
+}
+
+
+Int_t MIslandCalc::Calc1(){
+  
+  
+  /////////////////////////////
+  //
+  //        ALGORITHM # 1
+  // counts the number of algorithms as you can see in 
+  // the event display after doing the std image cleaning
+  //
+  /////////////////////////////
+  
   Float_t  noise;
   Float_t  signal;
@@ -333,7 +354,252 @@
   fIsl->SetReadyToSave();
   
-  return 1;
-  
+  return kTRUE;  
 }
 
 
+Int_t MIslandCalc::Calc2(){
+  
+  
+  /////////////////////////////
+  //
+  //        ALGORITHM # 2
+  // counts the number of islands considering as the same 
+  // islands the ones separated for 2 or less pixels
+  //
+  /////////////////////////////
+  
+  Float_t  noise;
+  Float_t  signal;
+
+  Int_t    npix = 577;
+
+  Int_t    sflag;
+  Int_t    control;
+  
+  Int_t    nvect = 0;
+  Int_t    fIslNum = 0;
+ 
+  Int_t    vect[50][577];
+
+  Int_t    zeros[50];
+  
+  Int_t kk[577];
+
+  for(Int_t m = 0; m < 50 ; m++)
+    for(Int_t n = 0; n < npix ; n++)
+	vect[m][n] = 0;
+    
+  for(Int_t n = 0; n < 50 ; n++)
+    zeros[n] = 0;
+  
+  for(Int_t n = 0; n < 577 ; n++)
+    kk[n] = 0;
+  
+  MCerPhotPix *pix;
+
+  //1st loop over all pixels
+  MCerPhotEvtIter Next0(fEvt, kFALSE);
+ 
+  while ((pix=static_cast<MCerPhotPix*>(Next0())))
+    {
+      const Int_t idx = pix->GetPixId();
+      
+      const MGeomPix &gpix  = (*fCam)[idx]; 
+      const Int_t    nnmax  = gpix.GetNumNeighbors();
+
+      if( fEvt->IsPixelUsed(idx))
+	{
+	  kk[idx] = 1 ;
+	  for(Int_t j=0; j< nnmax; j++)
+	    {
+	      kk[gpix.GetNeighbor(j)] = 1;
+	    }
+	} 
+      
+    }
+      
+  //2nd loop over all pixels
+  MCerPhotEvtIter Next(fEvt, kFALSE);
+  
+  while ((pix=static_cast<MCerPhotPix*>(Next())))
+    {
+      const Int_t idx = pix->GetPixId();
+      
+      const MGeomPix &gpix  = (*fCam)[idx];
+      const Int_t    nnmax  = gpix.GetNumNeighbors();
+      
+      if ( kk[idx] > 0)
+     	{
+	  sflag = 0;
+	  
+	  for(Int_t j=0; j < nnmax ; j++)
+	    {
+	      const Int_t idx2 = gpix.GetNeighbor(j);
+	      
+	      if (idx2 < idx)
+		{
+		  for(Int_t k = 1; k <= nvect; k++)
+		    {
+		      if (vect[k][idx2] == 1)
+			{
+			  sflag = 1;
+			  vect[k][idx] = 1;
+			}
+		    }
+		}
+	    }
+	  
+	  if (sflag == 0)
+	    {
+	      nvect++;
+	      vect[nvect][idx] = 1;	     
+	    }
+	  
+	}
+    }
+  
+  fIslNum = nvect;
+  
+  
+  // Repeated Chain Corrections
+  
+  for(Int_t i = 1; i <= nvect; i++)
+    {
+      for(Int_t j = i+1; j <= nvect; j++)
+	{
+	  control = 0;
+	  for(Int_t k = 0; k < npix; k++)
+	    {
+	      
+	      if (vect[i][k] == 1 && vect[j][k] == 1)
+		{
+		  control = 1; 
+		  break;
+		}
+	    }
+	  if (control == 1)
+	    {
+	      for(Int_t k = 0; k < npix; k++)
+		{
+		  if(vect[j][k] == 1)
+		    vect[i][k] = 1;
+		  vect[j][k] = 0;
+		  zeros[j] = 1;
+		}	
+	      fIslNum = fIslNum-1;	    
+	    }
+	  
+	}
+    }
+  
+  
+  Int_t l = 1;
+  
+  for(Int_t i = 1;  i<= nvect ; i++)
+    {
+      if (zeros[i] == 0)
+	{
+	  for(Int_t k = 0; k<npix; k++)
+	    {
+	      vect[l][k] = vect[i][k];
+	    }
+	  l++;
+	}
+    }
+  
+  
+  //set the number of islands in one event
+  fIsl->SetIslNum(fIslNum);
+
+  //examine each island...
+  Int_t fPixNum[fIslNum];
+  Float_t fSigToNoise[fIslNum];
+  Float_t time[577];
+  Float_t timeVariance[fIslNum];
+  
+  //reset the "sets" functions
+  if (fIslNum <1)
+    fIsl->SetIslNum(0);
+
+  for(Int_t i = 0; i<20 ;i++)
+    {
+      fIsl->SetPixNum(i,-1);
+      fIsl->SetSigToNoise(i,-1);
+      fIsl->SetTimeSpread(i,-1);
+    
+      for(Int_t idx = 0; idx<npix; idx++)
+	{
+	  fIsl->SetIslId(idx, -1);
+	  fIsl->SetArrivalTime(i, idx, -1 );
+	}
+    }
+ 
+   
+   for(Int_t i = 1; i<=fIslNum ; i++)
+    {
+      Int_t n = 0;
+      Int_t ncore = 0;
+      
+      Float_t MIN = 10000;
+      Float_t MAX = 0;
+
+      signal = 0;
+      noise = 0;
+      fPixNum[i-1] = 0;
+      timeVariance[i-1] = 0;
+      
+      for(Int_t idx=0 ; idx<npix ; idx++)
+	{
+	  
+	  MCerPhotPix *pix = fEvt->GetPixById(idx);
+	  const MPedestalPix &ped  = (*fPed)[idx];
+	  const MArrivalTimePix &timepix = (*fTime)[idx];
+
+	  if (pix == NULL) break;
+	    
+	  if (vect[i][idx] ==1 && fEvt->IsPixelUsed(idx)){
+	    
+	    fPixNum[i-1]++;
+	    signal += pix->GetNumPhotons() * (fCam->GetPixRatio(idx));
+	    noise += pow(ped.GetPedestalRms(),2);
+	    
+	    time[n] = timepix.IsLoGainUsed() ? timepix.GetArrivalTimeLoGain() : timepix.GetArrivalTimeHiGain();
+	    
+	    
+	    if (fEvt->IsPixelCore(idx)){ 
+	      
+	      if (time[n] > MAX)
+		MAX = time[n];
+	      if (time[n] < MIN)
+		MIN = time[n];
+	      
+	      ncore++;
+	    }
+	    
+	    fIsl->SetIslId(idx, i-1);
+	    fIsl->SetArrivalTime(i-1, idx, time[n]);
+	   	
+	    n++;
+	  }
+	  
+	}
+
+      // Float_t mean = timeVariance[i-1]/ncore;
+           
+      timeVariance[i-1] = 0;
+     
+      timeVariance[i-1] = (MAX - MIN)/ncore; 
+      timeVariance[i-1] = (MAX - MIN)/ncore; 
+
+      fSigToNoise[i-1]= (Float_t)signal/(Float_t)sqrt(noise);
+
+      fIsl->SetPixNum(i-1,fPixNum[i-1]);
+      fIsl->SetSigToNoise(i-1,fSigToNoise[i-1]);
+      fIsl->SetTimeSpread(i-1, timeVariance[i-1]);
+      
+    } 
+  
+  fIsl->SetReadyToSave();
+  
+  return 1;  
+}
Index: trunk/MagicSoft/Mars/mtemp/mifae/library/MIslandCalc.h
===================================================================
--- trunk/MagicSoft/Mars/mtemp/mifae/library/MIslandCalc.h	(revision 4222)
+++ trunk/MagicSoft/Mars/mtemp/mifae/library/MIslandCalc.h	(revision 4286)
@@ -31,7 +31,11 @@
   
     TString           fIslName;     //    name of the 'MIslands' container
- 
+  
+    Int_t  fIslandAlgorithm;
+
     Int_t PreProcess(MParList *plist);
     Int_t Process();
+    Int_t Calc1();                  // algorithm of counting islands #1
+    Int_t Calc2();                  // algorithm of counting islands #2    
            
  public:
@@ -39,4 +43,6 @@
     void SetOutputName(TString outname)   { fIslName = outname; }
     
+    void SetAlgorithm(Int_t m)   {fIslandAlgorithm = m;}
+    
     ClassDef(MIslandCalc, 0)        // task doing the image cleaning
 }; 
