HLT/BASE/AliHLTSpacePointContainer.cxx

   1 // $Id$
   2
   3 //**************************************************************************
   4 //* This file is property of and copyright by the ALICE HLT Project        *
   5 //* ALICE Experiment at CERN, All rights reserved.                         *
   6 //*                                                                        *
   7 //* Primary Authors: Matthias Richter <Matthias.Richter@ift.uib.no>        *
   8 //*                  for The ALICE HLT Project.                            *
   9 //*                                                                        *
  10 //* Permission to use, copy, modify and distribute this software and its   *
  11 //* documentation strictly for non-commercial purposes is hereby granted   *
  12 //* without fee, provided that the above copyright notice appears in all   *
  13 //* copies and that both the copyright notice and this permission notice   *
  14 //* appear in the supporting documentation. The authors make no claims     *
  15 //* about the suitability of this software for any purpose. It is          *
  16 //* provided "as is" without express or implied warranty.                  *
  17 //**************************************************************************
  18
  19 /// @file   AliHLTSpacePointContainer.cxx
  20 /// @author Matthias Richter
  21 /// @date   2011-04-29
  22 /// @brief  Base helper class for handling of HLT space point data blocks
  23 ///
  24
  25 #include "AliHLTSpacePointContainer.h"
  26 #include "AliHLTComponent.h"
  27 #include "TFile.h"
  28 #include "TString.h"
  29 #include "TArrayC.h"
  30 #include "TObjArray.h"
  31 #include "TH2F.h"
  32 #include "TMath.h"
  33 #include "TMarker.h"
  34 #include "TTree.h"
  35 #include <memory>
  36 #include <algorithm>
  37 #include <iostream>
  38 #include <iomanip>
  39
  40 /** ROOT macro for the implementation of ROOT specific class methods */
  41 ClassImp(AliHLTSpacePointContainer)
  42
  43 AliHLTSpacePointContainer::AliHLTSpacePointContainer()
  44   : TObject(), AliHLTLogging()
  45   , fBuffers()
  46 {
  47   // see header file for class documentation
  48   // or
  49   // refer to README to build package
  50   // or
  51   // visit http://web.ift.uib.no/~kjeks/doc/alice-hlt
  52 }
  53
  54 AliHLTSpacePointContainer::AliHLTSpacePointContainer(const AliHLTSpacePointContainer& c)
  55   : TObject(c), AliHLTLogging()
  56   , fBuffers()
  57 {
  58   /// copy constructor
  59
  60   // currently no need to copy or allocate, the fBuffers of the source object are expected
  61   // to be persistent throughout the lifetime of the new object. Only pointers are copied.
  62 }
  63
  64 AliHLTSpacePointContainer& AliHLTSpacePointContainer::operator=(const AliHLTSpacePointContainer& c)
  65 {
  66   /// assignment operator
  67   if (&c==this) return *this;
  68
  69   // currently no need to copy or allocate
  70   return *this;
  71 }
  72
  73 AliHLTSpacePointContainer::~AliHLTSpacePointContainer()
  74 {
  75   // destructor
  76   Clear();
  77 }
  78
  79 void AliHLTSpacePointContainer::Clear(Option_t * /*option*/)
  80 {
  81   // internal cleanup
  82   vector<TArrayC*>::iterator element=fBuffers.begin();
  83   while (element!=fBuffers.end()) {
  84     TArrayC* pBuffer=*element;
  85     element=fBuffers.erase(element);
  86     delete pBuffer;
  87   }
  88 }
  89
  90 void AliHLTSpacePointContainer::Print(Option_t *option) const
  91 {
  92   // print info
  93   Print(cout, option);
  94 }
  95
  96 void AliHLTSpacePointContainer::Print(ostream& out, Option_t */*option*/) const
  97 {
  98   // print to stream
  99   out << "AliHLTSpacePointContainer::Print" << endl;
 100 }
 101
 102 int AliHLTSpacePointContainer::AddInputBlock(const char* filename, AliHLTComponentDataType dt, unsigned specification)
 103 {
 104   // open block from file and add to collection
 105   if (!filename) return -EINVAL;
 106
 107   TString input=filename;
 108   input+="?filetype=raw";
 109   std::auto_ptr<TFile> pFile(new TFile(input));
 110   if (!pFile.get()) return -ENOMEM;
 111   if (pFile->IsZombie()) return -ENOENT;
 112
 113   int iResult=0;
 114   pFile->Seek(0);
 115   std::auto_ptr<TArrayC> buffer(new TArrayC);
 116   if (!buffer.get()) return -ENOMEM;
 117
 118   buffer->Set(pFile->GetSize());
 119   if (pFile->ReadBuffer(buffer->GetArray(), buffer->GetSize())==0) {
 120     AliHLTComponentBlockData bd;
 121     AliHLTComponent::FillBlockData(bd);
 122     bd.fPtr=buffer->GetArray();
 123     bd.fSize=buffer->GetSize();
 124     bd.fDataType=dt;
 125     bd.fSpecification=specification;
 126     HLTDebug("adding data block %d byte(s) from file %s", pFile->GetSize(), filename);
 127     iResult=AddInputBlock(&bd);
 128   } else {
 129     HLTError("failed reading %d byte(s) from file %s", pFile->GetSize(), filename);
 130     iResult=-ENODATA;
 131   }
 132
 133   fBuffers.push_back(buffer.release());
 134   return iResult;
 135 }
 136
 137 int AliHLTSpacePointContainer::AddInputBlocks(const char* listfile, AliHLTComponentDataType dt)
 138 {
 139   // open blank separated list of files and add data
 140   ifstream list(listfile);
 141   if (!list.good()) return -ENOENT;
 142
 143   int count=0;
 144   TString file;
 145   while (file.ReadLine(list)) {
 146     HLTInfo("adding data from file %s", file.Data());
 147     int iResult=AddInputBlock(file.Data(), dt, kAliHLTVoidDataSpec);
 148     if (iResult<0) {
 149       HLTInfo("failed to add data from file %s: error %d", file.Data(), iResult);
 150       return iResult;
 151     }
 152     count+=iResult;
 153   }
 154
 155   // std::auto_ptr<TObjArray> tokens(files.Tokenize(" "));
 156   // if (!tokens.get()) return 0;
 157   // for (int i=0; i<tokens->GetEntriesFast(); i++) {
 158   //   if (!tokens->At(i)) continue;
 159   //   cout << "adding data from file " << tokens->At(i)->GetName() << endl;
 160   //   AddInputBlock(tokens->At(i)->GetName(), dt, kAliHLTVoidDataSpec);
 161   // }
 162
 163   return count;
 164 }
 165
 166 AliHLTUInt8_t* AliHLTSpacePointContainer::Alloc(int size)
 167 {
 168   /// alloc memory for a space point data block
 169   TArrayC* buffer=new TArrayC;
 170   if (!buffer) return NULL;
 171
 172   buffer->Set(size);
 173   fBuffers.push_back(buffer);
 174   return reinterpret_cast<AliHLTUInt8_t*>(buffer->GetArray());
 175 }
 176
 177 AliHLTSpacePointContainer* AliHLTSpacePointContainer::SelectByMask(AliHLTUInt32_t /*mask*/, bool /*bAlloc*/) const
 178 {
 179   /// create a collection of clusters for a space point mask
 180   /// default implementation, nothing to do
 181   return NULL;
 182 }
 183
 184 AliHLTSpacePointContainer* AliHLTSpacePointContainer::SelectByTrack(int /*trackId*/, bool /*bAlloc*/) const
 185 {
 186   /// create a collection of clusters for a specific track
 187   /// default implementation, nothing to do
 188   return NULL;
 189 }
 190
 191 AliHLTSpacePointContainer* AliHLTSpacePointContainer::SelectByMC(int /*mcId*/, bool /*bAlloc*/) const
 192 {
 193   /// create a collection of clusters for a specific MC track
 194   /// default implementation, nothing to do
 195   return NULL;
 196 }
 197
 198 AliHLTSpacePointContainer* AliHLTSpacePointContainer::UsedClusters(bool /*bAlloc*/) const
 199 {
 200   /// create a collection of all used clusters
 201   /// default implementation, nothing to do
 202   return NULL;
 203 }
 204
 205 AliHLTSpacePointContainer* AliHLTSpacePointContainer::UnusedClusters(bool /*bAlloc*/) const
 206 {
 207   /// create a collection of all unused clusters
 208   /// default implementation, nothing to do
 209   return NULL;
 210 }
 211
 212 int AliHLTSpacePointContainer::MarkUsed(const AliHLTUInt32_t* /*clusterIDs*/, int /*arraySize*/)
 213 {
 214   /// default implementation, nothing to do
 215   return -ENOSYS;
 216 }
 217
 218 int AliHLTSpacePointContainer::SetTrackID(int /*trackID*/, const AliHLTUInt32_t* /*clusterIDs*/, int /*arraySize*/)
 219 {
 220   /// default implementation, nothing to do
 221   return -ENOSYS;
 222 }
 223
 224 int AliHLTSpacePointContainer::SetMCID(int /*mcID*/, const AliHLTUInt32_t* /*clusterIDs*/, int /*arraySize*/)
 225 {
 226   /// default implementation, nothing to do
 227   return -ENOSYS;
 228 }
 229
 230 int AliHLTSpacePointContainer::GetNumberOfSpacePoints() const
 231 {
 232   // get number of space points
 233   vector<AliHLTUInt32_t> clusterIDs;
 234   if (GetClusterIDs(clusterIDs)<0) return 0;
 235   return clusterIDs.size();
 236 }
 237
 238 bool AliHLTSpacePointContainer::Check(AliHLTUInt32_t id) const
 239 {
 240   // check if space point exists
 241   vector<AliHLTUInt32_t> clusterIDs;
 242   if (GetClusterIDs(clusterIDs)<0) return false;
 243   return find(clusterIDs.begin(), clusterIDs.end(), id)!=clusterIDs.end();
 244 }
 245
 246 TH1* AliHLTSpacePointContainer::DrawProjection(const char* plane, const vector<AliHLTUInt32_t>& selection) const
 247 {
 248   // draw the projection of space points in a specified plane
 249   // "xy", "xz", "yz"
 250
 251   int mode=0;
 252   if (strcmp(plane, "xy")==0) mode=0;
 253   else if (strcmp(plane, "xz")==0) mode=1;
 254   else if (strcmp(plane, "yz")==0) mode=2;
 255   else {
 256     HLTError("invalid plane specification %s, allowed 'xy', 'xz', 'yz'", plane);
 257     return NULL;
 258   }
 259
 260   float maxXAxis=100.0;
 261   float maxYAxis=100.0;
 262   vector<AliHLTUInt32_t> clusterIDs;
 263   GetClusterIDs(clusterIDs);
 264   vector<AliHLTUInt32_t>::iterator clusterID=clusterIDs.begin();
 265   while (clusterID!=clusterIDs.end()) {
 266     vector<AliHLTUInt32_t>::const_iterator element=selection.begin();
 267     for (; element!=selection.end(); element++) {
 268       if (((*element)&0xffff0000)==((*clusterID)&0xffff0000)) break;
 269     }
 270     if (element==selection.end() && selection.size()>0) {
 271       clusterID=clusterIDs.erase(clusterID);
 272       continue;
 273     }
 274     float XValue=0.0;
 275     float YValue=0.0;
 276     if (mode==0) {
 277       XValue=GetX(*clusterID);
 278       YValue=GetY(*clusterID);
 279     } else if (mode==1) {
 280       XValue=GetX(*clusterID);
 281       YValue=GetZ(*clusterID);
 282     } else if (mode==2) {
 283       XValue=GetY(*clusterID);
 284       YValue=GetZ(*clusterID);
 285     }
 286     if (maxXAxis<XValue) maxXAxis=XValue;
 287     if (maxYAxis<YValue) maxYAxis=YValue;
 288     clusterID++;
 289   }
 290   if (maxXAxis<maxYAxis) maxXAxis=maxYAxis;
 291   else maxYAxis=maxXAxis;
 292
 293   TH2F* histogram=new TH2F(plane, plane, 1000, -maxXAxis, maxXAxis, 1000, -maxYAxis, maxYAxis);
 294   if (!histogram) return NULL;
 295     if (mode==0) {
 296       histogram->GetXaxis()->SetTitle("X [cm]");
 297       histogram->GetYaxis()->SetTitle("Y [cm]");
 298     } else if (mode==1) {
 299       histogram->GetXaxis()->SetTitle("Z [cm]");
 300       histogram->GetYaxis()->SetTitle("X [cm]");
 301     } else if (mode==2) {
 302       histogram->GetXaxis()->SetTitle("Z [cm]");
 303       histogram->GetYaxis()->SetTitle("Y [cm]");
 304     }
 305
 306   for (clusterID=clusterIDs.begin(); clusterID!=clusterIDs.end(); clusterID++) {
 307     float XValue=GetX(*clusterID);
 308     float YValue=GetY(*clusterID);
 309     float ZValue=GetZ(*clusterID);
 310     Float_t phi     = GetPhi(*clusterID);
 311     Float_t cos     = TMath::Cos( phi );
 312     Float_t sin     = TMath::Sin( phi );
 313     if (mode==0) {
 314       histogram->SetTitle("XY projection");
 315       histogram->Fill(cos*XValue - sin*YValue, sin*XValue + cos*YValue);
 316     } else if (mode==1) {
 317       // should be maybe 'ZX' but somehow 'XZ' is more 'lingual convenient'
 318       histogram->SetTitle("XZ projection");
 319       histogram->Fill(ZValue, cos*XValue - sin*YValue);
 320     } else if (mode==2) {
 321       // same comment
 322       histogram->SetTitle("YZ projection");
 323       histogram->Fill(ZValue, sin*XValue + cos*YValue);
 324     }
 325   }
 326
 327   return histogram;
 328 }
 329
 330 void AliHLTSpacePointContainer::Draw(Option_t *option)
 331 {
 332   /// Inherited from TObject, draw clusters
 333   float scale=250;
 334   float center[2]={0.5,0.5};
 335   int markersize=1;
 336   int markercolor=5;
 337
 338   TString strOption(option);
 339   std::auto_ptr<TObjArray> tokens(strOption.Tokenize(" "));
 340   if (!tokens.get()) return;
 341   for (int i=0; i<tokens->GetEntriesFast(); i++) {
 342     if (!tokens->At(i)) continue;
 343     const char* key="";
 344     TString arg=tokens->At(i)->GetName();
 345
 346     key="scale=";
 347     if (arg.BeginsWith(key)) {
 348       arg.ReplaceAll(key, "");
 349       scale=arg.Atof();
 350     }
 351     key="centerx=";
 352     if (arg.BeginsWith(key)) {
 353       arg.ReplaceAll(key, "");
 354       center[0]=arg.Atof();
 355     }
 356     key="centery=";
 357     if (arg.BeginsWith(key)) {
 358       arg.ReplaceAll(key, "");
 359       center[1]=arg.Atof();
 360     }
 361     key="markersize=";
 362     if (arg.BeginsWith(key)) {
 363       arg.ReplaceAll(key, "");
 364       markersize=arg.Atoi();
 365     }
 366     key="markercolor=";
 367     if (arg.BeginsWith(key)) {
 368       arg.ReplaceAll(key, "");
 369       markercolor=arg.Atoi();
 370     }
 371   }
 372
 373   vector<AliHLTUInt32_t> clusterIDs;
 374   GetClusterIDs(clusterIDs);
 375   for (vector<AliHLTUInt32_t>::iterator clusterID=clusterIDs.begin();
 376        clusterID!=clusterIDs.end();
 377        clusterID++) {
 378     float clusterphi=GetPhi(*clusterID);
 379     float cosphi=TMath::Cos(clusterphi);
 380     float sinphi=TMath::Sin(clusterphi);
 381     float clusterx=GetX(*clusterID);
 382     float clustery=GetY(*clusterID);
 383     //cout << " " << *clusterID << " " << clusterx << " " << clustery << " " << clusterphi << endl;
 384     // rotate
 385     float pointx= clusterx*sinphi + clustery*cosphi;
 386     float pointy=-clusterx*cosphi + clustery*sinphi;
 387
 388     TMarker* m=new TMarker(pointx/(2*scale)+center[0], pointy/(2*scale)+center[1], 3);
 389     m->SetMarkerSize(markersize);
 390     m->SetMarkerColor(markercolor);
 391     m->Draw("same");
 392   }
 393 }
 394
 395 TTree* AliHLTSpacePointContainer::FillTree(const char* name, const char* title)
 396 {
 397   // create tree and fill the space point coordinates
 398   TString treename=name;
 399   TString treetitle=title;
 400   if (treename.IsNull()) treename="spacepoints";
 401   if (treetitle.IsNull()) treetitle="HLT space point coordinates";
 402   std::auto_ptr<TTree> tree(new TTree(treename, treetitle));
 403   if (!tree.get()) return NULL;
 404
 405   const unsigned dimension=8;
 406   float values[dimension];
 407   memset(values, 0, sizeof(values));
 408   const char* names[dimension]={
 409     "x", "y", "z", "sigmaY2", "sigmaZ2", "charge", "qmax", "alpha"
 410   };
 411
 412   for (unsigned i=0; i<dimension; i++) {
 413     TString type=names[i]; type+="/F";
 414     tree->Branch(names[i], &values[i], type);
 415   }
 416
 417   const vector<AliHLTUInt32_t>* clusterIDs=GetClusterIDs(0xffffffff);
 418   for (vector<AliHLTUInt32_t>::const_iterator clusterID=clusterIDs->begin();
 419        clusterID!=clusterIDs->end();
 420        clusterID++) {
 421     unsigned pos=0;
 422     values[pos++]=GetX(*clusterID);
 423     values[pos++]=GetY(*clusterID);
 424     values[pos++]=GetZ(*clusterID);
 425     values[pos++]=GetYWidth(*clusterID);
 426     values[pos++]=GetZWidth(*clusterID);
 427     values[pos++]=GetCharge(*clusterID);
 428     values[pos++]=GetMaxSignal(*clusterID);
 429     values[pos++]=GetPhi(*clusterID);
 430
 431     tree->Fill();
 432   }
 433
 434   return tree.release();
 435 }
 436
 437 ostream& operator<<(ostream &out, const AliHLTSpacePointContainer& c)
 438 {
 439   c.Print(out);
 440   return out;
 441 }
 442
 443 ostream& operator<<(ostream &out, const AliHLTSpacePointContainer::AliHLTSpacePointProperties& p)
 444 {
 445   // print
 446   cout << p.fId;
 447   return out;
 448 }
 449
 450 bool operator==(const AliHLTSpacePointContainer::AliHLTSpacePointProperties& a,
 451                 const AliHLTSpacePointContainer::AliHLTSpacePointProperties& b) {
 452   return a.fId==b.fId;
 453 }