RALICE/icepack/IceCleanHits.cxx

   1 /*******************************************************************************
   2  * Copyright(c) 2003, IceCube Experiment at the South Pole. All rights reserved.
   3  *
   4  * Author: The IceCube RALICE-based Offline 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.
  12  * The authors make no claims about the suitability of this software for
  13  * any purpose. It is provided "as is" without express or implied warranty.
  14  *******************************************************************************/
  15
  16 // $Id$
  17
  18 ///////////////////////////////////////////////////////////////////////////
  19 // Class IceCleanHits
  20 // TTask derived class to perform hit cleaning.
  21 //
  22 // In case an event has been rejected by an AliEventSelector (based) processor,
  23 // this task (and its sub-tasks) is not executed.
  24 //
  25 // The code in this processor is based on the algorithms as developed
  26 // by Oladipo Fadiran and George Japaridze (Clark Atlanta University, USA).
  27 //
  28 // Criteria applied for Amanda MuDaq data :
  29 // ----------------------------------------
  30 // 1) ADC within [min,max]  Default : [0.3,999999] PE
  31 // 2) TOT within [min,max]  Default : electrical [125,2000] optical [20,2000] ns
  32 // 3) abs(LE-Ttrig)<=win    Default : win=2250 TDC counts
  33 //    where : LE=uncalibrated hit LE (i.e. TDC counts)   Ttrig=trigger pulse LE in TDC counts
  34 // 4) At least one hit in another OM within radius R and time difference dt
  35 //    to remove isolated hits. Defaults : R=70 m  dt=500 ns
  36 //
  37 // Criteria applied for Amanda TWRDaq data :
  38 // ----------------------------------------
  39 // 1) ADC within [min,max]  Default : [0.3,999999] PE
  40 // 2) TOT within [min,max]  Default : electrical [125,2000] optical [20,2000] ns
  41 // 3) abs(LE-Ttrig)<=win    Default : win=3000 ns
  42 //    where : LE=uncalibrated hit LE    Ttrig=trigger pulse LE
  43 // 4) At least one hit in another OM within radius R and time difference dt
  44 //    to remove isolated hits. Defaults : R=70 m  dt=500 ns
  45 //
  46 // The actual DAQ system is obtained automatically from the IceEvent structure
  47 // via the device called "Daq".
  48 //
  49 // The defaults of the various parameters can be changed by the corresponding
  50 // Set memberfunctions.
  51 //
  52 // Information about the actual parameter settings can be found in the event
  53 // structure itself via the device named "IceCleanHits".
  54 //
  55 // Concerning the trigger time :
  56 // -----------------------------
  57 // By default the trigger time is obtained automatically from the IceEvent structure
  58 // via the device called "Trigger".
  59 // The uncalibrated LE of a specified trigger pulse is used.
  60 // The user can impose a specific trigger name or time to be used
  61 // by invokation of the memberfunctions SetTnameA or SetTtimeA, respectively.
  62 // Specification of a negative trigger time will result in the automatic
  63 // trigger time setting corresponding to the "main" trigger.
  64 // By default the trigger time of the "main" trigger will be used.
  65 //
  66 // The hits which do not fullfill the criteria are flagged "dead" for the
  67 // corresponding signal slot. This means they are still present in the
  68 // IceEvent structure and are as such still accessible.
  69 // It is left to the user to decide (based on the various "dead" flag settings)
  70 // whether or not to use these hits in his/her reconstruction or analysis.
  71 //
  72 // Note : This processor only works properly on Time and ADC calibrated data.
  73 //
  74 //--- Author: Nick van Eijndhoven 13-oct-2005 Utrecht University
  75 //- Modified: NvE $Date$ Utrecht University
  76 ///////////////////////////////////////////////////////////////////////////
  77
  78 #include "IceCleanHits.h"
  79 #include "Riostream.h"
  80
  81 ClassImp(IceCleanHits) // Class implementation to enable ROOT I/O
  82
  83 IceCleanHits::IceCleanHits(const char* name,const char* title) : TTask(name,title)
  84 {
  85 // Default constructor.
  86  fEvt=0;
  87  fAdcminAM=0.3;
  88  fAdcmaxAM=999999;
  89  fAdcminAT=0.3;
  90  fAdcmaxAT=999999;
  91  fTotminAEM=125;
  92  fTotmaxAEM=2000;
  93  fTotminAOM=20;
  94  fTotmaxAOM=2000;
  95  fTotminAET=125;
  96  fTotmaxAET=2000;
  97  fTotminAOT=20;
  98  fTotmaxAOT=2000;
  99  fRmaxA=70;
 100  fDtmaxA=500;
 101  fTwinAM=2250;
 102  fTwinAT=3000;
 103  fTtimAM=-1;
 104  fTnamAM="main";
 105  fTtimAT=-1;
 106  fTnamAT="main";
 107 }
 108 ///////////////////////////////////////////////////////////////////////////
 109 IceCleanHits::~IceCleanHits()
 110 {
 111 // Default destructor.
 112 }
 113 ///////////////////////////////////////////////////////////////////////////
 114 void IceCleanHits::SetAdcRangeA(Float_t min,Float_t max,TString s)
 115 {
 116 // Set Amanda ADC range in PE.
 117 // The default for the maximum is 999999.
 118 // The argument "s" specifies the DAQ system (i.e. "MuDaq" or "TWRDaq").
 119 // For backward compatibility the default is s="MuDaq".
 120
 121  if (s=="MuDaq")
 122  {
 123   fAdcminAM=min;
 124   fAdcmaxAM=max;
 125  }
 126  if (s=="TWRDaq")
 127  {
 128   fAdcminAT=min;
 129   fAdcmaxAT=max;
 130  }
 131 }
 132 ///////////////////////////////////////////////////////////////////////////
 133 void IceCleanHits::SetTotRangeAE(Float_t min,Float_t max,TString s)
 134 {
 135 // Set Amanda electrical TOT range in ns.
 136 // The default for the maximum is 2000.
 137 // The argument "s" specifies the DAQ system (i.e. "MuDaq" or "TWRDaq").
 138 // For backward compatibility the default is s="MuDaq".
 139
 140  if (s=="MuDaq")
 141  {
 142   fTotminAEM=min;
 143   fTotmaxAEM=max;
 144  }
 145  if (s=="TWRDaq")
 146  {
 147   fTotminAET=min;
 148   fTotmaxAET=max;
 149  }
 150 }
 151 ///////////////////////////////////////////////////////////////////////////
 152 void IceCleanHits::SetTotRangeAO(Float_t min,Float_t max,TString s)
 153 {
 154 // Set Amanda optical TOT range in ns.
 155 // The default for the maximum is 2000.
 156 // The argument "s" specifies the DAQ system (i.e. "MuDaq" or "TWRDaq").
 157 // For backward compatibility the default is s="MuDaq".
 158
 159  if (s=="MuDaq")
 160  {
 161   fTotminAOM=min;
 162   fTotmaxAOM=max;
 163  }
 164  if (s=="TWRDaq")
 165  {
 166   fTotminAOT=min;
 167   fTotmaxAOT=max;
 168  }
 169 }
 170 ///////////////////////////////////////////////////////////////////////////
 171 void IceCleanHits::SetIsolationA(Float_t rmax,Float_t dtmax)
 172 {
 173 // Set Amanda isolation radius (in m) and time difference (in ns).
 174  fRmaxA=rmax;
 175  fDtmaxA=dtmax;
 176 }
 177 ///////////////////////////////////////////////////////////////////////////
 178 void IceCleanHits::SetTwindowA(Float_t dtmax,TString s)
 179 {
 180 // Set Amanda maximal trigger window.
 181 // Only hits which occur in [T-dtmax,T+dtmax] will be kept,
 182 // where T indicates the trigger time.
 183 // For the Amanda MuDaq hardware, the times are all in TDC counts,
 184 // where 1 TDC corresponds to about 1.04 ns.
 185 // For the TWRDaq, the times are all in nanoseconds.
 186 // The argument "s" specifies the DAQ system (i.e. "MuDaq" or "TWRDaq").
 187 // For backward compatibility the default is s="MuDaq".
 188
 189  if (s=="MuDaq") fTwinAM=dtmax;
 190  if (s=="TWRDaq") fTwinAT=dtmax;
 191 }
 192 ///////////////////////////////////////////////////////////////////////////
 193 void IceCleanHits::SetTtimeA(Float_t t,TString s)
 194 {
 195 // Set Amanda trigger time.
 196 // For the Amanda MuDaq hardware, the times are all in TDC counts,
 197 // where 1 TDC corresponds to about 1.04 ns.
 198 // For the TWRDaq, the times are all in nanoseconds.
 199 // The argument "s" specifies the DAQ system (i.e. "MuDaq" or "TWRDaq").
 200 // For backward compatibility the default is s="MuDaq".
 201 // A negative value will induce automatic trigger time setting based
 202 // on "main" trigger as recorded in the IceEvent structure.
 203
 204  if (s=="MuDaq")
 205  {
 206   fTtimAM=t;
 207   fTnamAM="user";
 208   if (t<0) fTnamAM="main";
 209  }
 210  if (s=="TWRDaq")
 211  {
 212   fTtimAT=t;
 213   fTnamAT="user";
 214   if (t<0) fTnamAT="main";
 215  }
 216 }
 217 ///////////////////////////////////////////////////////////////////////////
 218 void IceCleanHits::SetTnameA(TString name,TString s)
 219 {
 220 // Set Amanda trigger name.
 221 // The argument "s" specifies the DAQ system (i.e. "MuDaq" or "TWRDaq").
 222 // For backward compatibility the default is s="MuDaq".
 223 // Specification of a non-existing trigger name will result in a trigger time
 224 // value of 0.
 225
 226  if (s=="MuDaq")
 227  {
 228   fTtimAM=0;
 229   fTnamAM=name;
 230  }
 231  if (s=="TWRDaq")
 232  {
 233   fTtimAT=0;
 234   fTnamAT=name;
 235  }
 236 }
 237 ///////////////////////////////////////////////////////////////////////////
 238 void IceCleanHits::Exec(Option_t* opt)
 239 {
 240 // Implementation of the hit cleaning procedures.
 241
 242  TString name=opt;
 243  AliJob* parent=(AliJob*)(gROOT->GetListOfTasks()->FindObject(name.Data()));
 244
 245  if (!parent) return;
 246
 247  fEvt=(IceEvent*)parent->GetObject("IceEvent");
 248  if (!fEvt) return;
 249
 250  // Only process accepted events
 251  AliDevice* seldev=(AliDevice*)fEvt->GetDevice("AliEventSelector");
 252  if (seldev)
 253  {
 254   if (seldev->GetSignal("Select") < 0.1) return;
 255  }
 256
 257  // Storage of the used parameters in the IceCleanHits device
 258  AliSignal params;
 259  params.SetNameTitle("IceCleanHits","IceCleanHits processor parameters");
 260  params.SetSlotName("AdcminAM",1);
 261  params.SetSlotName("AdcmaxAM",2);
 262  params.SetSlotName("TotminAEM",3);
 263  params.SetSlotName("TotmaxAEM",4);
 264  params.SetSlotName("TotminAOM",5);
 265  params.SetSlotName("TotmaxAOM",6);
 266  params.SetSlotName("RmaxA",7);
 267  params.SetSlotName("DtmaxA",8);
 268  params.SetSlotName("TwinAM",9);
 269  params.SetSlotName("TtimAM",10);
 270  params.SetSlotName("AdcminAT",11);
 271  params.SetSlotName("AdcmaxAT",12);
 272  params.SetSlotName("TotminAET",13);
 273  params.SetSlotName("TotmaxAET",14);
 274  params.SetSlotName("TotminAOT",15);
 275  params.SetSlotName("TotmaxAOT",16);
 276  params.SetSlotName("TwinAT",17);
 277  params.SetSlotName("TtimAT",18);
 278
 279  params.SetSignal(fAdcminAM,1);
 280  params.SetSignal(fAdcmaxAM,2);
 281  params.SetSignal(fTotminAEM,3);
 282  params.SetSignal(fTotmaxAEM,4);
 283  params.SetSignal(fTotminAOM,5);
 284  params.SetSignal(fTotmaxAOM,6);
 285  params.SetSignal(fRmaxA,7);
 286  params.SetSignal(fDtmaxA,8);
 287  params.SetSignal(fTwinAM,9);
 288  params.SetSignal(fTtimAM,10);
 289  params.SetSignal(fAdcminAT,11);
 290  params.SetSignal(fAdcmaxAT,12);
 291  params.SetSignal(fTotminAET,13);
 292  params.SetSignal(fTotmaxAET,14);
 293  params.SetSignal(fTotminAOT,15);
 294  params.SetSignal(fTotmaxAOT,16);
 295  params.SetSignal(fTwinAT,17);
 296  params.SetSignal(fTtimAT,18);
 297
 298  fEvt->AddDevice(params);
 299
 300  Amanda();
 301  InIce();
 302  IceTop();
 303 }
 304 ///////////////////////////////////////////////////////////////////////////
 305 void IceCleanHits::Amanda()
 306 {
 307 // Hit cleaning for Amanda modules.
 308
 309  AliDevice* daq=(AliDevice*)fEvt->GetDevice("Daq");
 310  if (!daq) return;
 311
 312  if (daq->GetSignal("Muon")) MuDaq();
 313  if (daq->GetSignal("TWR")) TWRDaq();
 314 }
 315 ///////////////////////////////////////////////////////////////////////////
 316 void IceCleanHits::MuDaq()
 317 {
 318 // Hit cleaning for Amanda MuDaq data.
 319
 320  // Trigger time setting according to the user selection.
 321  // By default the "main" trigger of the event will be used.
 322  if (fTnamAM != "user")
 323  {
 324   AliDevice* tdev=(AliDevice*)fEvt->GetDevice("Trigger");
 325   if (tdev)
 326   {
 327    AliSignal* trig=tdev->GetHit(fTnamAM);
 328    if (trig) fTtimAM=trig->GetSignal("trig_pulse_le");
 329   }
 330  }
 331
 332  // All Amanda OMs with a signal
 333  TObjArray* aoms=fEvt->GetDevices("IceAOM");
 334  if (!aoms) return;
 335
 336  // Local OM array with bad/dead OMs (as indicated via IceCalibrate) discarded
 337  TObjArray oms;
 338  IceAOM* omx=0;
 339  for (Int_t i=0; i<aoms->GetEntries(); i++)
 340  {
 341   omx=(IceAOM*)aoms->At(i);
 342   if (!omx) continue;
 343   if (omx->GetDeadValue("ADC") || omx->GetDeadValue("LE") || omx->GetDeadValue("TOT")) continue;
 344   oms.Add(omx);
 345  }
 346
 347  // Local array with clean hits
 348  TObjArray hits;
 349  Int_t omid=0;
 350  Int_t clean=1;
 351  AliSignal* sx=0;
 352  Float_t adc,le,tot;
 353  Int_t readout;
 354  for (Int_t iom=0; iom<oms.GetEntries(); iom++)
 355  {
 356   omx=(IceAOM*)oms.At(iom);
 357   if (!omx) continue;
 358   omid=omx->GetUniqueID();
 359   readout=int(omx->GetSignal("READOUT"));
 360   // General readout setting in case info was missing
 361   if (!readout)
 362   {
 363    readout=1;
 364    if (omid>=303) readout=2; // Optical OMs
 365   }
 366   for (Int_t ih=1; ih<=omx->GetNhits(); ih++)
 367   {
 368    sx=omx->GetHit(ih);
 369    if (!sx) continue;
 370    adc=sx->GetSignal("ADC",7);
 371    le=sx->GetSignal("LE",7);
 372    tot=sx->GetSignal("TOT",7);
 373
 374    // Remove hits with an ADC value outside the range
 375    if (adc<fAdcminAM || adc>fAdcmaxAM)
 376    {
 377     sx->SetDead("ADC");
 378     clean=0;
 379    }
 380    // Remove hits with a TOT value outside the range
 381    // Note : Different ranges for electrical and optical modules
 382    if (readout==1) // Electrical OMs
 383    {
 384     if (tot<fTotminAEM || tot>fTotmaxAEM)
 385     {
 386      sx->SetDead("TOT");
 387      clean=0;
 388     }
 389    }
 390    else // Optical OMs
 391    {
 392     if (tot<fTotminAOM || tot>fTotmaxAOM)
 393     {
 394      sx->SetDead("TOT");
 395      clean=0;
 396     }
 397    }
 398    // Remove hits that are outside the trigger time window.
 399    // Since the trigger time was determined from uncalibrated LE's
 400    // (to include cable length effects) the uncalibrated LE of each
 401    // hit should be used here as well.
 402    le=sx->GetSignal("LE",-7);
 403    if (fabs(le-fTtimAM)>fTwinAM)
 404    {
 405      sx->SetDead("LE");
 406      clean=0;
 407    }
 408    // Store only the current clean hits in our local hit array
 409    // This will save CPU time for the isolation criterion
 410    if (clean) hits.Add(sx);
 411   }
 412  }
 413
 414  // Isolation cut
 415  // Only retain hits that have at least one hit of another OM within a certain
 416  // radius and within a certain time window
 417  Int_t nhits=hits.GetEntries();
 418  AliSignal* sx1=0;
 419  AliSignal* sx2=0;
 420  Float_t t1,t2;
 421  IceAOM* omx1=0;
 422  IceAOM* omx2=0;
 423  AliPosition r1;
 424  AliPosition r2;
 425  Float_t dt,dr;
 426  Int_t iso;
 427  for (Int_t jh1=0; jh1<nhits; jh1++)
 428  {
 429   sx1=(AliSignal*)hits.At(jh1);
 430   if (!sx1) continue;
 431
 432   iso=1;
 433   for (Int_t jh2=0; jh2<nhits; jh2++)
 434   {
 435    sx2=(AliSignal*)hits.At(jh2);
 436    if (!sx2) continue;
 437
 438    omx1=(IceAOM*)sx1->GetDevice();
 439    omx2=(IceAOM*)sx2->GetDevice();
 440    if (omx1==omx2) continue;
 441
 442    t1=sx1->GetSignal("LE",7);
 443    t2=sx2->GetSignal("LE",7);
 444    dt=fabs(t2-t1);
 445    if (dt>fDtmaxA) continue;
 446
 447    if (omx1 && omx2)
 448    {
 449     r1=omx1->GetPosition();
 450     r2=omx2->GetPosition();
 451     dr=r1.GetDistance(r2);
 452     if (dr>fRmaxA) continue;
 453     iso=0;
 454    }
 455   }
 456   if (iso) sx1->SetDead("LE");
 457  }
 458 }
 459 ///////////////////////////////////////////////////////////////////////////
 460 void IceCleanHits::TWRDaq()
 461 {
 462 // Hit cleaning for Amanda TWRDaq data.
 463
 464  // Trigger time setting according to the user selection.
 465  // By default the "main" trigger of the event will be used.
 466  if (fTnamAT != "user")
 467  {
 468   AliDevice* tdev=(AliDevice*)fEvt->GetDevice("Trigger");
 469   if (tdev)
 470   {
 471    AliSignal* trig=tdev->GetHit(fTnamAT);
 472    if (trig) fTtimAT=trig->GetSignal("trig_pulse_le");
 473   }
 474  }
 475
 476  // All Amanda OMs with a signal
 477  TObjArray* aoms=fEvt->GetDevices("IceAOM");
 478  if (!aoms) return;
 479
 480  // Local OM array with bad/dead OMs (as indicated via IceCalibrate) discarded
 481  TObjArray oms;
 482  IceAOM* omx=0;
 483  for (Int_t i=0; i<aoms->GetEntries(); i++)
 484  {
 485   omx=(IceAOM*)aoms->At(i);
 486   if (!omx) continue;
 487   if (omx->GetDeadValue("ADC") || omx->GetDeadValue("LE") || omx->GetDeadValue("TOT")) continue;
 488   oms.Add(omx);
 489  }
 490
 491  // Local array with clean hits
 492  TObjArray hits;
 493  Int_t omid=0;
 494  Int_t clean=1;
 495  AliSignal* sx=0;
 496  Float_t adc,le,tot;
 497  Int_t readout;
 498  for (Int_t iom=0; iom<oms.GetEntries(); iom++)
 499  {
 500   omx=(IceAOM*)oms.At(iom);
 501   if (!omx) continue;
 502   omid=omx->GetUniqueID();
 503   readout=int(omx->GetSignal("READOUT"));
 504   // General readout setting in case info was missing
 505   if (!readout)
 506   {
 507    readout=1;
 508    if (omid>=303) readout=2; // Optical OMs
 509   }
 510   for (Int_t ih=1; ih<=omx->GetNhits(); ih++)
 511   {
 512    sx=omx->GetHit(ih);
 513    if (!sx) continue;
 514    adc=sx->GetSignal("ADC",7);
 515    le=sx->GetSignal("LE",7);
 516    tot=sx->GetSignal("TOT",7);
 517
 518    // Remove hits with an ADC value outside the range
 519    if (adc<fAdcminAT || adc>fAdcmaxAT)
 520    {
 521     sx->SetDead("ADC");
 522     clean=0;
 523    }
 524    // Remove hits with a TOT value outside the range
 525    // Note : Different ranges for electrical and optical modules
 526    if (readout==1) // Electrical OMs
 527    {
 528     if (tot<fTotminAET || tot>fTotmaxAET)
 529     {
 530      sx->SetDead("TOT");
 531      clean=0;
 532     }
 533    }
 534    else // Optical OMs
 535    {
 536     if (tot<fTotminAOT || tot>fTotmaxAOT)
 537     {
 538      sx->SetDead("TOT");
 539      clean=0;
 540     }
 541    }
 542    // Remove hits that are outside the trigger time window.
 543    // Since the trigger time was determined from uncalibrated LE's
 544    // (to include cable length effects) the uncalibrated LE of each
 545    // hit should be used here as well.
 546    le=sx->GetSignal("LE",-7);
 547    if (fabs(le-fTtimAT)>fTwinAT)
 548    {
 549      sx->SetDead("LE");
 550      clean=0;
 551    }
 552    // Store only the current clean hits in our local hit array
 553    // This will save CPU time for the isolation criterion
 554    if (clean) hits.Add(sx);
 555   }
 556  }
 557
 558  // Isolation cut
 559  // Only retain hits that have at least one hit of another OM within a certain
 560  // radius and within a certain time window
 561  Int_t nhits=hits.GetEntries();
 562  AliSignal* sx1=0;
 563  AliSignal* sx2=0;
 564  Float_t t1,t2;
 565  IceAOM* omx1=0;
 566  IceAOM* omx2=0;
 567  AliPosition r1;
 568  AliPosition r2;
 569  Float_t dt,dr;
 570  Int_t iso;
 571  for (Int_t jh1=0; jh1<nhits; jh1++)
 572  {
 573   sx1=(AliSignal*)hits.At(jh1);
 574   if (!sx1) continue;
 575
 576   iso=1;
 577   for (Int_t jh2=0; jh2<nhits; jh2++)
 578   {
 579    sx2=(AliSignal*)hits.At(jh2);
 580    if (!sx2) continue;
 581
 582    omx1=(IceAOM*)sx1->GetDevice();
 583    omx2=(IceAOM*)sx2->GetDevice();
 584    if (omx1==omx2) continue;
 585
 586    t1=sx1->GetSignal("LE",7);
 587    t2=sx2->GetSignal("LE",7);
 588    dt=fabs(t2-t1);
 589    if (dt>fDtmaxA) continue;
 590
 591    if (omx1 && omx2)
 592    {
 593     r1=omx1->GetPosition();
 594     r2=omx2->GetPosition();
 595     dr=r1.GetDistance(r2);
 596     if (dr>fRmaxA) continue;
 597     iso=0;
 598    }
 599   }
 600   if (iso) sx1->SetDead("LE");
 601  }
 602 }
 603 ///////////////////////////////////////////////////////////////////////////
 604 void IceCleanHits::InIce()
 605 {
 606 // Hit cleaning for IceCube InIce DOMs.
 607 }
 608 ///////////////////////////////////////////////////////////////////////////
 609 void IceCleanHits::IceTop()
 610 {
 611 // Hit cleaning for IceTop DOMs.
 612 }
 613 ///////////////////////////////////////////////////////////////////////////