Sign error for dphi move of jet axis close to phi = 0 corrected.
authormorsch <morsch@f7af4fe6-9843-0410-8265-dc069ae4e863>
Thu, 27 Sep 2007 08:16:11 +0000 (08:16 +0000)
committermorsch <morsch@f7af4fe6-9843-0410-8265-dc069ae4e863>
Thu, 27 Sep 2007 08:16:11 +0000 (08:16 +0000)
JETAN/AliUA1JetFinderV1.cxx
JETAN/AliUA1JetFinderV2.cxx

index 80edb16..aacda75 100644 (file)
@@ -97,16 +97,17 @@ void AliUA1JetFinderV1::FindJets()
     etaT[i] = lv->Eta();
     phiT[i] = ((lv->Phi() < 0) ? (lv->Phi()) + 2 * TMath::Pi() : lv->Phi());
     if (fReader->GetCutFlag(i) != 1) continue;
-    fLego->Fill(etaT[i], phiT[i], ptT[i]);
+    fLego ->Fill(etaT[i], phiT[i], ptT[i]);
     hPtTotal->Fill(ptT[i]);
     etbgTotal+= ptT[i];
   }
+
   fJets->SetNinput(nIn);
 
   // calculate total energy and fluctuation in map
   Double_t meanpt = hPtTotal->GetMean();
-  Double_t ptRMS = hPtTotal->GetRMS();
-  Double_t npart = hPtTotal->GetEntries();
+  Double_t ptRMS  = hPtTotal->GetRMS();
+  Double_t npart  = hPtTotal->GetEntries();
   Double_t dEtTotal = (TMath::Sqrt(npart))*TMath::Sqrt(meanpt * meanpt + ptRMS*ptRMS);
 
   // arrays to hold jets
@@ -246,10 +247,7 @@ void AliUA1JetFinderV1::RunAlgoritm(Float_t etbgTotal, Double_t dEtTotal, Int_t&
   // check enough space! *to be done*
   AliUA1JetHeaderV1* header = (AliUA1JetHeaderV1*) fHeader;
   const Int_t nBinsMax = 70000;
-  if(header->GetLegoNbinEta()*header->GetLegoNbinPhi()>nBinsMax){
-    printf("%s:%d \n >> Too Large Eta-Phi Grid Canceling Jet search \n",(char*)__FILE__,__LINE__);
-    return;
-  }
+  
   Float_t etCell[nBinsMax];   //! Cell Energy
   Float_t etaCell[nBinsMax];  //! Cell eta
   Float_t phiCell[nBinsMax];  //! Cell phi
@@ -294,22 +292,23 @@ void AliUA1JetFinderV1::RunAlgoritm(Float_t etbgTotal, Double_t dEtTotal, Int_t&
   Int_t * index  = new Int_t[nCell];
   TMath::Sort(nCell, etCell, index);
   // variable used in centroide loop
-  Float_t eta = 0.0;
-  Float_t phi = 0.0;
-  Float_t eta0 = 0.0;
-  Float_t phi0 = 0.0;
-  Float_t etab = 0.0;
-  Float_t phib = 0.0;
-  Float_t etas = 0.0;
-  Float_t phis = 0.0;
-  Float_t ets = 0.0;
-  Float_t deta = 0.0;
-  Float_t dphi = 0.0;
-  Float_t dr = 0.0;
-  Float_t etsb = 0.0;
+  Float_t eta   = 0.0;
+  Float_t phi   = 0.0;
+  Float_t eta0  = 0.0;
+  Float_t phi0  = 0.0;
+  Float_t etab  = 0.0;
+  Float_t phib  = 0.0;
+  Float_t etas  = 0.0;
+  Float_t phis  = 0.0;
+  Float_t ets   = 0.0;
+  Float_t deta  = 0.0;
+  Float_t dphi  = 0.0;
+  Float_t dr    = 0.0;
+  Float_t etsb  = 0.0;
   Float_t etasb = 0.0;
   Float_t phisb = 0.0;
-
+  Float_t dphib = 0.0;
+  
 
   for(Int_t icell = 0; icell < nCell; icell++){
         Int_t jcell = index[icell];
@@ -329,21 +328,21 @@ void AliUA1JetFinderV1::RunAlgoritm(Float_t etbgTotal, Double_t dEtTotal, Int_t&
         phisb = 0.0;
         for(Int_t kcell =0; kcell < nCell; kcell++){
             Int_t lcell = index[kcell];
-            if(lcell == jcell) continue; // cell itself
-            if(flagCell[lcell] != 0) continue; // cell used before
-            if(etCell[lcell] > etCell[jcell]) continue;
+            if(lcell == jcell)                 continue; // cell itself
+            if(flagCell[lcell] != 0)           continue; // cell used before
+            if(etCell[lcell] > etCell[jcell])  continue; // can this happen
             //calculate dr
             deta = etaCell[lcell] - eta;
-                dphi = phiCell[lcell] - phi;
-                if (dphi < -TMath::Pi()) dphi= -dphi - 2.0 * TMath::Pi();
-                if (dphi > TMath::Pi()) dphi = 2.0 * TMath::Pi() - dphi;
-                dr = TMath::Sqrt(deta * deta + dphi * dphi);
+           dphi = TMath::Abs(phiCell[lcell] - phi);
+           if (dphi > TMath::Pi()) dphi = 2. * TMath::Pi() - dphi;
+           dr = TMath::Sqrt(deta * deta + dphi * dphi);
             if(dr <= rc){
                // calculate offset from initiate cell
                deta = etaCell[lcell] - eta0;
                dphi = phiCell[lcell] - phi0;
-               if (dphi < -TMath::Pi()) dphi= -dphi - 2.0 * TMath::Pi();
-              if (dphi > TMath::Pi()) dphi = 2.0 * TMath::Pi() - dphi;
+               if (dphi < - TMath::Pi()) dphi=  dphi + 2.0 * TMath::Pi();
+              if (dphi >   TMath::Pi()) dphi = dphi - 2.0 * TMath::Pi();
+              
                etas = etas + etCell[lcell]*deta;
                phis = phis + etCell[lcell]*dphi;
                ets = ets + etCell[lcell];
@@ -351,24 +350,26 @@ void AliUA1JetFinderV1::RunAlgoritm(Float_t etbgTotal, Double_t dEtTotal, Int_t&
                eta = eta0 + etas/ets;
                phi = phi0 + phis/ets;
                // if cone does not move much, just go to next step
-               dr = TMath::Sqrt((eta-etab)*(eta-etab) + (phi-phib)*(phi-phib));
+              dphib = TMath::Abs(phi - phib);
+              if (dphib > TMath::Pi()) dphib = 2. * TMath::Pi() - dphib;
+              dr = TMath::Sqrt((eta-etab)*(eta-etab) + dphib * dphib);
                if(dr <= minmove) break;
                // cone should not move more than max_mov
                dr = TMath::Sqrt((etas/ets)*(etas/ets) + (phis/ets)*(phis/ets));
                if(dr > maxmove){
-                  eta = etab;
-                  phi = phib;
-                  ets = etsb;
-                  etas = etasb;
-                  phis = phisb;
-               }else{ // store this loop information
-                 etab=eta;
-                 phib=phi;
-                 etsb = ets;
-                 etasb = etas;
-                 phisb = phis;
+                  eta = etab;
+                  phi = phib;
+                  ets = etsb;
+                  etas = etasb;
+                  phis = phisb;
+              } else { // store this loop information
+                  etab  = eta;
+                  phib  = phi;
+                  etsb  = ets;
+                  etasb = etas;
+                  phisb = phis;
                }
-            }
+            } // inside cone
         }//end of cells loop looking centroide
 
         //avoid cones overloap (to be implemented in the future)
@@ -433,8 +434,8 @@ void AliUA1JetFinderV1::RunAlgoritm(Float_t etbgTotal, Double_t dEtTotal, Int_t&
 
 
   //delete
-  delete [] index;
-  delete [] idx;
+  delete[] index;
+  delete[] idx;
 
 }
 ////////////////////////////////////////////////////////////////////////
@@ -832,4 +833,5 @@ void AliUA1JetFinderV1::Init()
         header->GetLegoEtaMax(),  header->GetLegoNbinPhi(),
         header->GetLegoPhiMin(),  header->GetLegoPhiMax());
 
+
 }
index 81da6b2..0182689 100644 (file)
@@ -331,64 +331,67 @@ void AliUA1JetFinderV2::RunAlgoritm(Int_t nIn, Float_t* etCell, Float_t* etaCell
   TMath::Sort(nCell, etCell, index);
 
   // variable used in centroide loop
-  Float_t eta = 0.0;
-  Float_t phi = 0.0;
-  Float_t eta0 = 0.0;
-  Float_t phi0 = 0.0;
-  Float_t etab = 0.0;
-  Float_t phib = 0.0;
-  Float_t etas = 0.0;
-  Float_t phis = 0.0;
-  Float_t ets = 0.0;
-  Float_t deta = 0.0;
-  Float_t dphi = 0.0;
-  Float_t dr = 0.0;
-  Float_t etsb = 0.0;
+  Float_t eta   = 0.0;
+  Float_t phi   = 0.0;
+  Float_t eta0  = 0.0;
+  Float_t phi0  = 0.0;
+  Float_t etab  = 0.0;
+  Float_t phib  = 0.0;
+  Float_t etas  = 0.0;
+  Float_t phis  = 0.0;
+  Float_t ets   = 0.0;
+  Float_t deta  = 0.0;
+  Float_t dphi  = 0.0;
+  Float_t dr    = 0.0;
+  Float_t etsb  = 0.0;
   Float_t etasb = 0.0;
   Float_t phisb = 0.0;
+  Float_t dphib = 0.0;
 
 
   for(Int_t icell = 0; icell < nCell; icell++){
-        Int_t jcell = index[icell];
-        if(etCell[jcell] <= etseed) continue; // if cell energy is low et seed
-        if(flagCell[jcell] != 0) continue; // if cell was used before
-        eta  = etaCell[jcell];
-        phi  = phiCell[jcell];
-        eta0 = eta;
-        phi0 = phi;
-        etab = eta;
-        phib = phi;
-        ets  = etCell[jcell];
-        etas = 0.0;
-        phis = 0.0;
-        etsb = ets;
-        etasb = 0.0;
-        phisb = 0.0;
-        for(Int_t kcell =0; kcell < nCell; kcell++){
-            Int_t lcell = index[kcell];
-            if(lcell == jcell) continue; // cell itself
-            if(flagCell[lcell] != 0) continue; // cell used before
-            if(etCell[lcell] > etCell[jcell]) continue;
-            //calculate dr
-            deta = etaCell[lcell] - eta;
-                dphi = phiCell[lcell] - phi;
-                if (dphi < -TMath::Pi()) dphi= -dphi - 2.0 * TMath::Pi();
-                if (dphi > TMath::Pi()) dphi = 2.0 * TMath::Pi() - dphi;
-                dr = TMath::Sqrt(deta * deta + dphi * dphi);
-            if(dr <= rc){
-               // calculate offset from initiate cell
-               deta = etaCell[lcell] - eta0;
-               dphi = phiCell[lcell] - phi0;
-               if (dphi < -TMath::Pi()) dphi= -dphi - 2.0 * TMath::Pi();
-                   if (dphi > TMath::Pi()) dphi = 2.0 * TMath::Pi() - dphi;
-               etas = etas + etCell[lcell]*deta;
-               phis = phis + etCell[lcell]*dphi;
-               ets = ets + etCell[lcell];
+      Int_t jcell = index[icell];
+      if(etCell[jcell] <= etseed) continue; // if cell energy is low et seed
+      if(flagCell[jcell] != 0) continue; // if cell was used before
+      eta  = etaCell[jcell];
+      phi  = phiCell[jcell];
+      eta0 = eta;
+      phi0 = phi;
+      etab = eta;
+      phib = phi;
+      ets  = etCell[jcell];
+      etas = 0.0;
+      phis = 0.0;
+      etsb = ets;
+      etasb = 0.0;
+      phisb = 0.0;
+      for(Int_t kcell =0; kcell < nCell; kcell++){
+         Int_t lcell = index[kcell];
+         if(lcell == jcell) continue; // cell itself
+         if(flagCell[lcell] != 0) continue; // cell used before
+         if(etCell[lcell] > etCell[jcell]) continue;
+         //calculate dr
+         deta = etaCell[lcell] - eta;
+         dphi = phiCell[lcell] - phi;
+         if (dphi < -TMath::Pi()) dphi= -dphi - 2.0 * TMath::Pi();
+         if (dphi > TMath::Pi()) dphi = 2.0 * TMath::Pi() - dphi;
+         dr = TMath::Sqrt(deta * deta + dphi * dphi);
+         if(dr <= rc){
+             // calculate offset from initiate cell
+             deta = etaCell[lcell] - eta0;
+             dphi = phiCell[lcell] - phi0;
+             if (dphi < -TMath::Pi()) dphi = dphi + 2.0 * TMath::Pi();
+             if (dphi >  TMath::Pi()) dphi = dphi - 2.0 * TMath::Pi();
+             etas = etas + etCell[lcell]*deta;
+             phis = phis + etCell[lcell]*dphi;
+             ets = ets + etCell[lcell];
                //new weighted eta and phi including this cell
                eta = eta0 + etas/ets;
                phi = phi0 + phis/ets;
                // if cone does not move much, just go to next step
-               dr = TMath::Sqrt((eta-etab)*(eta-etab) + (phi-phib)*(phi-phib));
+              dphib = TMath::Abs(phi - phib);
+              if (dphib > TMath::Pi()) dphib = 2. * TMath::Pi() - dphib;
+              dr = TMath::Sqrt((eta-etab)*(eta-etab) + dphib * dphib);
                if(dr <= minmove) break;
                // cone should not move more than max_mov
                dr = TMath::Sqrt((etas/ets)*(etas/ets) + (phis/ets)*(phis/ets));