| 930d0706 |
1 | /************************************************************************* |
| 2 | * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * |
| 3 | * * |
| 4 | * Author: The ALICE Off-line 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. The authors make no claims * |
| 12 | * about the suitability of this software for any purpose. It is * |
| 13 | * provided "as is" without express or implied warranty. * |
| 14 | **************************************************************************/ |
| 15 | |
| 16 | // |
| 17 | // Implementation of a branch replicator |
| 18 | // to produce aods with only few branches. |
| 19 | // |
| 20 | // This replicator is in charge of replicating the nuclei primary vertices |
| 21 | // tracks identified as nuclei with Z>=2, secondary vertices in form of |
| 22 | // AliAODRecoDecayLF2Prong and their daughter tracks. |
| 23 | // These informations are stored into a reduced AODs (AliAOD.NuclEx.root) |
| 24 | // |
| 25 | // The vertices are filtered so that only the primary vertices make it |
| 26 | // to the output aods. |
| 27 | // |
| 28 | // The secondary vertices are recreated here, as a AliAODRecoDecayLF2Prong |
| 29 | // plus cuts that select secondary vericesoutside the primary vertex |
| 30 | |
| 31 | // Authors: S. Bufalino (stefania.bufalino@cern.ch) |
| 32 | // R. Lea (ramona.lea@cern.ch) |
| 33 | // Based on AliAODMuonReplicator.cxx |
| 34 | |
| 35 | //NOTE : This is a test on MC : no PID response only MC truth + select only 3LH |
| 36 | // daughters : NOT INTENDED for any efficiency!! |
| 37 | |
| 38 | class AliESDv0; |
| 39 | class AliESDVertex; |
| 40 | class AliAODVertex; |
| 41 | class AliAODRecoDecay; |
| 42 | |
| 43 | #include "AliStack.h" |
| 44 | #include "AliAODEvent.h" |
| 45 | #include "AliAODMCHeader.h" |
| 46 | #include "AliAODMCParticle.h" |
| 47 | #include "AliAODTZERO.h" |
| 48 | #include "AliAODTrack.h" |
| 49 | #include "AliAODVZERO.h" |
| 50 | //#include "AliAnalysisCuts.h" |
| 51 | #include "TF1.h" |
| 52 | #include "AliExternalTrackParam.h" |
| 53 | #include "AliESDv0.h" |
| 54 | #include "AliAODv0.h" |
| 55 | //#include "AliPIDResponse.h" |
| 56 | #include <iostream> |
| 57 | #include <cassert> |
| 58 | #include "AliESDtrack.h" |
| 59 | #include "TObjArray.h" |
| 60 | #include "AliAnalysisFilter.h" |
| 61 | #include "AliAODRecoDecay.h" |
| 62 | #include "AliAODRecoDecayLF.h" |
| 63 | #include "AliAODRecoDecayLF2Prong.h" |
| 64 | |
| 65 | #include <TFile.h> |
| 66 | #include <TDatabasePDG.h> |
| 67 | #include <TString.h> |
| 68 | #include <TList.h> |
| 69 | #include "AliLog.h" |
| 70 | #include "AliVEvent.h" |
| 71 | #include "AliVVertex.h" |
| 72 | #include "AliVTrack.h" |
| 73 | #include "AliVertexerTracks.h" |
| 74 | #include "AliKFVertex.h" |
| 75 | #include "AliESDEvent.h" |
| 76 | #include "AliESDVertex.h" |
| 77 | #include "AliESDtrackCuts.h" |
| 78 | #include "AliAODEvent.h" |
| 79 | #include "AliAnalysisFilter.h" |
| 80 | //#include "AliAnalysisVertexingLF.h" |
| 81 | #include "AliAnalysisManager.h" |
| 82 | #include "AliAODMCNuclExReplicator.h" |
| 83 | #include "TH1.h" |
| 84 | #include "TCanvas.h" |
| 85 | #include "AliInputEventHandler.h" |
| 86 | |
| 87 | using std::cout; |
| 88 | using std::endl; |
| 89 | |
| 90 | ClassImp(AliAODMCNuclExReplicator) |
| 91 | |
| 92 | //_____________________________________________________________________________ |
| 93 | AliAODMCNuclExReplicator::AliAODMCNuclExReplicator(const char* name, const char* title, |
| 94 | Int_t mcMode, |
| 95 | Int_t nsigmaTrk1, Int_t partType1, |
| 96 | Int_t nsigmaTrk2, Int_t partType2 |
| 97 | ) |
| 98 | :AliAODBranchReplicator(name,title), |
| 99 | fBzkG(0.), |
| 100 | fCosMin(), |
| 101 | fDCAtracksMin(), |
| 102 | fRmax(), |
| 103 | fRmin(), |
| 104 | fDNmin(), |
| 105 | fDPmin(), |
| 106 | fHeader(0x0), |
| 107 | fVertices(0x0), |
| 108 | fNuclei(0x0), |
| 109 | fSecondaryVerices(0x0), |
| 110 | fDaughterTracks(0x0), |
| 111 | fList(0x0), |
| 112 | fMCParticles(0x0), |
| 113 | fMCHeader(0x0), |
| 114 | fMCMode(mcMode), |
| 115 | fLabelMap(), |
| 116 | fParticleSelected(), |
| 117 | fReplicateHeader(kTRUE), //replicateHeader //to be fixed |
| 118 | fnSigmaTrk1(nsigmaTrk1), |
| 119 | fnSigmaTrk2(nsigmaTrk2), |
| 120 | fpartType1(partType1), |
| 121 | fpartType2(partType2), |
| 122 | fSecVtxWithKF(kFALSE), |
| 123 | fVertexerTracks(0x0), |
| 124 | fV1(0x0), |
| 125 | fAODMapSize(0), |
| 126 | fAODMap(0) |
| 127 | |
| 128 | { |
| 129 | // default ctor |
| 130 | } |
| 131 | |
| 132 | //_____________________________________________________________________________ |
| 133 | AliAODMCNuclExReplicator::~AliAODMCNuclExReplicator() |
| 134 | { |
| 135 | // destructor |
| 136 | // delete fTrackCut; |
| 137 | // delete fVertexCut; |
| 138 | delete fList; |
| 139 | } |
| 140 | |
| 141 | //_____________________________________________________________________________ |
| 142 | void AliAODMCNuclExReplicator::SelectParticle(Int_t i) |
| 143 | { |
| 144 | // taking the absolute values here, need to take care |
| 145 | // of negative daughter and mother |
| 146 | // IDs when setting! |
| 147 | |
| 148 | if (!IsParticleSelected(TMath::Abs(i))) |
| 149 | { |
| 150 | fParticleSelected.Add(TMath::Abs(i),1); |
| 151 | } |
| 152 | } |
| 153 | |
| 154 | //_____________________________________________________________________________ |
| 155 | Bool_t AliAODMCNuclExReplicator::IsParticleSelected(Int_t i) |
| 156 | { |
| 157 | // taking the absolute values here, need to take |
| 158 | // care with negative daughter and mother |
| 159 | // IDs when setting! |
| 160 | return (fParticleSelected.GetValue(TMath::Abs(i))==1); |
| 161 | } |
| 162 | |
| 163 | |
| 164 | //_____________________________________________________________________________ |
| 165 | void AliAODMCNuclExReplicator::CreateLabelMap(const AliAODEvent& source) |
| 166 | { |
| 167 | // |
| 168 | // this should be called once all selections are done |
| 169 | // |
| 170 | |
| 171 | fLabelMap.Delete(); |
| 172 | |
| 173 | TClonesArray* mcParticles = static_cast<TClonesArray*>(source.FindListObject(AliAODMCParticle::StdBranchName())); |
| 174 | |
| 175 | Int_t i(0); |
| 176 | Int_t j(0); |
| 177 | |
| 178 | TIter next(mcParticles); |
| 179 | |
| 180 | while ( next() ) |
| 181 | { |
| 182 | if (IsParticleSelected(i)) |
| 183 | { |
| 184 | fLabelMap.Add(i,j++); |
| 185 | } |
| 186 | ++i; |
| 187 | } |
| 188 | } |
| 189 | |
| 190 | //_____________________________________________________________________________ |
| 191 | Int_t AliAODMCNuclExReplicator::GetNewLabel(Int_t i) |
| 192 | { |
| 193 | // Gets the label from the new created Map |
| 194 | // Call CreatLabelMap before |
| 195 | // otherwise only 0 returned |
| 196 | return fLabelMap.GetValue(TMath::Abs(i)); |
| 197 | } |
| 198 | |
| 199 | |
| 200 | //_____________________________________________________________________________ |
| 201 | TList* AliAODMCNuclExReplicator::GetList() const |
| 202 | { |
| 203 | // return (and build if not already done) our internal list of managed objects |
| 204 | if (!fList) |
| 205 | { |
| 206 | |
| 207 | if ( fReplicateHeader ) |
| 208 | { |
| 209 | fHeader = new AliAODHeader; |
| 210 | } |
| 211 | |
| 212 | |
| 213 | fSecondaryVerices = new TClonesArray("AliAODRecoDecayLF2Prong",30); |
| 214 | fSecondaryVerices->SetName("SecondaryVertices"); |
| 215 | |
| 216 | fVertices = new TClonesArray("AliAODVertex",2); |
| 217 | fVertices->SetName("vertices"); |
| 218 | |
| 219 | fNuclei = new TClonesArray("AliAODTrack",30); |
| 220 | fNuclei->SetName("Nuclei"); |
| 221 | |
| 222 | fDaughterTracks = new TClonesArray("AliAODTrack",30); |
| 223 | fDaughterTracks->SetName("DaughterTracks"); |
| 224 | |
| 225 | |
| 226 | fList = new TList; |
| 227 | fList->SetOwner(kTRUE); |
| 228 | |
| 229 | fList->Add(fHeader); |
| 230 | fList->Add(fVertices); |
| 231 | fList->Add(fNuclei); |
| 232 | fList->Add(fSecondaryVerices); |
| 233 | fList->Add(fDaughterTracks); |
| 234 | |
| 235 | |
| 236 | if ( fMCMode > 0 ) |
| 237 | { |
| 238 | fMCHeader = new AliAODMCHeader; |
| 239 | fMCParticles = new TClonesArray("AliAODMCParticle",1000); |
| 240 | fMCParticles->SetName(AliAODMCParticle::StdBranchName()); |
| 241 | fList->Add(fMCHeader); |
| 242 | fList->Add(fMCParticles); |
| 243 | } |
| 244 | } |
| 245 | return fList; |
| 246 | } |
| 247 | |
| 248 | //_____________________________________________________________________________ |
| 249 | void AliAODMCNuclExReplicator::ReplicateAndFilter(const AliAODEvent& source) |
| 250 | { |
| 251 | // Replicate (and filter if filters are there) the relevant parts we're interested in AODEvent |
| 252 | |
| 253 | // cout<<"-------------------->QUESTO"<<endl; |
| 254 | |
| 255 | //----------------------------------------------- |
| 256 | // AliPIDResponse |
| 257 | |
| 258 | // AliAnalysisManager *man=AliAnalysisManager::GetAnalysisManager(); |
| 259 | // AliInputEventHandler* inputHandler = (AliInputEventHandler*) (man->GetInputEventHandler()); |
| 260 | |
| 261 | //-------------------------------------------------------- |
| 262 | |
| 263 | // printf("Execute NuclEx Replicator\n"); |
| 264 | |
| 265 | //--------------------------------- |
| 266 | |
| 267 | if (fReplicateHeader) |
| 268 | { |
| 0a918d8d |
269 | AliAODHeader * header = dynamic_cast<AliAODHeader*>(source.GetHeader()); |
| 270 | if(!header) AliFatal("Not a standard AOD"); |
| 271 | *fHeader = *(header); |
| 930d0706 |
272 | } |
| 273 | |
| 274 | fVertices->Clear("C"); |
| 275 | |
| 276 | fNuclei->Clear("C"); |
| 277 | |
| 278 | fSecondaryVerices->Clear("C"); |
| 279 | |
| 280 | fDaughterTracks->Clear("C"); |
| 281 | |
| 282 | //---------------------------------- |
| 283 | |
| 284 | //retrive MC infos |
| 285 | |
| 286 | TClonesArray *arrayMC = 0; |
| 287 | AliAODMCHeader *mcHeader=0; |
| 288 | Int_t mumpdg=-100; |
| 7b23fc90 |
289 | // Int_t mumpdgNeg=-100; |
| be4df180 |
290 | |
| 291 | |
| 930d0706 |
292 | arrayMC = (TClonesArray*) source.GetList()->FindObject(AliAODMCParticle::StdBranchName()); |
| 293 | if (!arrayMC) { |
| 294 | Printf("Error: MC particles branch not found!\n"); |
| 295 | return; |
| 296 | } |
| be4df180 |
297 | |
| 930d0706 |
298 | // if(arrayMC) |
| 299 | // cout<<"Ho caricato array mc"<<endl; |
| 300 | |
| 301 | mcHeader = (AliAODMCHeader*)source.GetList()->FindObject(AliAODMCHeader::StdBranchName()); |
| 302 | if(!mcHeader) { |
| 303 | printf("AliAnalysisTaskSEDs::UserExec: MC header branch not found!\n"); |
| 304 | return; |
| 305 | } |
| be4df180 |
306 | |
| 930d0706 |
307 | // if(mcHeader) |
| 308 | // cout<<"Ho caricato MC header"<<endl; |
| be4df180 |
309 | |
| 930d0706 |
310 | // cout<<"Centrality AOD source: "<<source.GetHeader()->GetCentrality()<<endl; |
| 311 | |
| 312 | Int_t nsv(0); |
| 313 | // Int_t nnuclei(0); |
| 314 | Int_t ntracksD(0); |
| 315 | |
| 316 | Int_t input(0); |
| 317 | Double_t xdummy,ydummy; |
| 318 | |
| 319 | AliAODRecoDecayLF2Prong *io2Prong = 0; |
| 320 | |
| 321 | TObjArray *twoTrackArray = new TObjArray(2); |
| 322 | Double_t dispersion; |
| 323 | |
| 324 | // cout<<"Qui"<<endl; |
| be4df180 |
325 | // cout<<source.GetMagneticField()<<endl; |
| 930d0706 |
326 | |
| 327 | AliAODVertex *vtx = source.GetPrimaryVertex(); |
| 328 | |
| be4df180 |
329 | // cout<<"Source "<<source<<endl; |
| 330 | // cout<<"vtx: "<<vtx<<endl; |
| 930d0706 |
331 | |
| 332 | // A Set of very loose cut for a weak strange decay |
| 333 | |
| 334 | fCosMin = 0.99; |
| 335 | fDCAtracksMin = 1; |
| 336 | fRmax = 200.; |
| 337 | fRmin = 0.1; |
| 338 | fDNmin = 0.05; |
| 339 | fDPmin = 0.05; |
| 340 | |
| 341 | //---------------------------------------------------------- |
| 342 | |
| 343 | // Int_t nindices=0; |
| 344 | UShort_t *indices = 0; |
| 345 | const Int_t entries = source.GetNumberOfTracks(); |
| 346 | |
| 347 | Double_t pos[3],cov[6]; |
| 348 | vtx->GetXYZ(pos); |
| 349 | vtx->GetCovarianceMatrix(cov); |
| 350 | fV1 = new AliESDVertex(pos,cov,100.,100,vtx->GetName()); |
| be4df180 |
351 | cout<<"fV1 pricipal loop: "<<fV1<<endl; |
| 930d0706 |
352 | |
| 353 | if(entries<=0) return; |
| 354 | indices = new UShort_t[entries]; |
| 355 | memset(indices,0,sizeof(UShort_t)*entries); |
| 356 | fAODMapSize = 100000; |
| 357 | fAODMap = new Int_t[fAODMapSize]; |
| 358 | memset(fAODMap,0,sizeof(Int_t)*fAODMapSize); |
| 359 | // cent=((AliAODEvent&)source)->GetCentrality(); |
| 360 | |
| 361 | //------------------------------------------------------------- |
| 362 | |
| 363 | if(vtx->GetNContributors()<1) { |
| 364 | |
| 365 | // SPD vertex cut |
| 366 | vtx =source.GetPrimaryVertexSPD(); |
| 367 | |
| 368 | if(vtx->GetNContributors()<1) { |
| 369 | Info("AliAnalysisTaskHelium3Pi","No good vertex, skip event"); |
| 370 | return; // NO GOOD VERTEX, SKIP EVENT |
| 371 | } |
| 372 | } |
| 373 | |
| 374 | Double_t xPrimaryVertex=0.,yPrimaryVertex=0.; |
| 375 | xPrimaryVertex=vtx->GetX(); |
| 376 | yPrimaryVertex=vtx->GetY(); |
| 377 | |
| 378 | fBzkG=source.GetMagneticField(); |
| 379 | fVertexerTracks=new AliVertexerTracks(fBzkG); |
| 380 | |
| 381 | Double_t TrackNumber = source.GetNumberOfTracks(); |
| 382 | Int_t label =-1; |
| 383 | |
| 384 | //Tracks arrays |
| 385 | |
| 386 | TArrayI Track0(TrackNumber); //Pions |
| 387 | Int_t nTrack0=0; |
| 388 | |
| 389 | TArrayI Track1(TrackNumber); //Helium3 |
| 390 | Int_t nTrack1=0; |
| 391 | |
| 392 | for(Int_t j=0; j<TrackNumber; j++){ |
| 393 | |
| 394 | // cout<<"Inside loop tracks"<<endl; |
| 395 | |
| 396 | |
| 397 | AliVTrack *track = (AliVTrack*)source.GetTrack(j); |
| 398 | |
| 399 | AliAODTrack *aodtrack =(AliAODTrack*)track; |
| 400 | |
| 401 | //----------------------------------------------------------- |
| 402 | //Track cuts |
| 403 | if(aodtrack->GetTPCNcls() < 70 )continue; |
| 404 | if(aodtrack->Chi2perNDF() > 4 )continue; |
| 405 | |
| 406 | if (!aodtrack->IsOn(AliAODTrack::kTPCrefit)) continue; |
| 407 | if (!aodtrack->IsOn(AliAODTrack::kTPCin)) continue; |
| 408 | if (aodtrack->IsOn(AliAODTrack::kITSpureSA)) continue; |
| 409 | |
| 410 | //--------------------------------------------------------------- |
| 411 | |
| 412 | Double_t mom = aodtrack->P(); |
| 413 | |
| 414 | if(mom<0.150)continue; |
| 415 | |
| 416 | label = TMath::Abs(aodtrack->GetLabel()); |
| 417 | |
| 418 | AliAODMCParticle *part = (AliAODMCParticle*) arrayMC->At(label); |
| 419 | |
| 420 | Int_t PDGCode=part->GetPdgCode(); |
| 421 | |
| 422 | Int_t mumid = part->GetMother(); |
| 423 | |
| 424 | if(mumid>-1){ |
| 425 | AliAODMCParticle *mother = (AliAODMCParticle*) arrayMC->At(mumid); |
| 426 | mumpdg = mother->GetPdgCode(); |
| 427 | } |
| 428 | |
| 429 | // if(mumpdg == 1010010030 ||mumpdg == -1010010030 ){ |
| 430 | |
| 2f82005f |
431 | if(mumpdg == -1010010030){ |
| 930d0706 |
432 | |
| be4df180 |
433 | if(PDGCode==-211 || PDGCode==+211){ |
| 434 | // if(PDGCode==-211){ |
| 435 | |
| 436 | // if(aodtrack->Charge()>0)continue; |
| 437 | |
| 930d0706 |
438 | Track0[nTrack0++]=j; |
| 439 | } |
| 440 | |
| be4df180 |
441 | if(PDGCode==1000020030 ||PDGCode==-1000020030 ){ |
| 442 | //if(PDGCode==1000020030){ |
| 443 | |
| 444 | // if(aodtrack->Charge()<0)continue; |
| 445 | |
| 930d0706 |
446 | Track1[nTrack1++]=j; |
| 447 | // new((*fNuclei)[nnuclei++]) AliAODTrack(*aodtrack); |
| 448 | } |
| 449 | } |
| 450 | } |
| 451 | |
| 452 | //Set Track Daughters |
| 453 | Track0.Set(nTrack0); |
| 454 | Track1.Set(nTrack1); |
| 455 | |
| 456 | |
| 457 | // cout<<"Track loops..."<<endl; |
| 458 | // cout<<"npos "<<nTrack1<<endl; |
| 459 | // cout<<"nneg "<<nTrack0<<endl; |
| 460 | |
| 461 | |
| 462 | AliAODTrack *postrackAOD = 0; |
| 463 | AliAODTrack *negtrackAOD = 0; |
| 464 | |
| 465 | AliESDtrack *postrack = 0; |
| 466 | AliESDtrack *negtrack = 0; |
| 467 | |
| 468 | Bool_t isOk=kFALSE; |
| 469 | |
| 470 | for (Int_t i=0; i<nTrack1; i++){ //! He Tracks Loop |
| 471 | |
| 472 | Int_t Track1idx=Track1[i]; |
| 473 | |
| 474 | AliVTrack *trackpos = (AliVTrack*)source.GetTrack(Track1idx); |
| 475 | |
| 476 | postrackAOD = (AliAODTrack*)trackpos; |
| 477 | postrack = new AliESDtrack(trackpos); |
| 478 | |
| 479 | //--- MC infos |
| 480 | |
| 481 | Int_t labelpos = TMath::Abs(postrack->GetLabel()); |
| 482 | AliAODMCParticle *partPos = (AliAODMCParticle*) arrayMC->At(labelpos); |
| 483 | Int_t mumidPos = partPos->GetMother(); |
| 484 | |
| be4df180 |
485 | // cout<<"\n\nmumidPos: "<<mumidPos <<endl; |
| 486 | |
| 930d0706 |
487 | //------------------------------ |
| 488 | |
| 489 | for (Int_t k=0; k <nTrack0 ; k++) { //! Pion Tracks Loop |
| 490 | |
| 491 | Int_t Track0idx=Track0[k]; |
| 492 | |
| 493 | AliVTrack *trackneg = (AliVTrack*)source.GetTrack(Track0idx); |
| 494 | negtrackAOD =(AliAODTrack*)trackneg; |
| 495 | negtrack = new AliESDtrack(trackneg); |
| 496 | |
| 497 | //--- MC infos |
| 498 | |
| 499 | Int_t labelneg = TMath::Abs(negtrack->GetLabel()); |
| 500 | AliAODMCParticle *partNeg = (AliAODMCParticle*) arrayMC->At(labelneg); |
| 501 | Int_t mumidNeg = partNeg->GetMother(); |
| 502 | |
| be4df180 |
503 | // cout<<"mumidNeg: "<<mumidNeg <<endl; |
| 504 | |
| 505 | // AliAODMCParticle *motherNeg = (AliAODMCParticle*) arrayMC->At(mumidNeg); |
| 506 | // mumpdgNeg = motherNeg->GetPdgCode(); |
| 507 | |
| 930d0706 |
508 | //------------------------------ |
| 509 | // if(mumidPos == mumidNeg && mumidNeg > 0){ |
| 510 | isOk=kFALSE; |
| 511 | |
| be4df180 |
512 | // if(mumidPos == mumidNeg && mumidNeg > 0 && mumpdgNeg == 1010010030) |
| 930d0706 |
513 | if(mumidPos == mumidNeg && mumidNeg > 0) |
| 514 | isOk = kTRUE; |
| 515 | |
| 516 | twoTrackArray->AddAt(negtrack,0); |
| 517 | twoTrackArray->AddAt(postrack,1); |
| 518 | |
| 519 | Double_t dcap1n1 = postrack->GetDCA(negtrack,fBzkG,xdummy,ydummy); |
| 520 | |
| 521 | Double_t dcap1toPV = TMath::Abs(postrack->GetD(xPrimaryVertex, yPrimaryVertex,fBzkG)); |
| 522 | Double_t dcan1toPV = TMath::Abs(negtrack->GetD(xPrimaryVertex, yPrimaryVertex,fBzkG)); |
| 523 | |
| 524 | if(dcap1n1>fDCAtracksMin)continue; |
| 525 | if((xdummy+ydummy)>fRmax )continue; |
| 526 | if((xdummy+ydummy)< fRmin)continue; |
| 527 | |
| 528 | if ( dcan1toPV < fDNmin) |
| 529 | if ( dcap1toPV < fDPmin) continue; |
| 530 | |
| 531 | // cout<<"dcap1n1: "<<dcap1n1<<endl; |
| 532 | |
| 533 | AliAODVertex *vertexp1n1 = ReconstructSecondaryVertex(twoTrackArray,dispersion,kTRUE); |
| 534 | |
| 535 | if(!vertexp1n1) { |
| 536 | |
| 537 | twoTrackArray->Clear(); |
| 538 | delete negtrack; |
| 539 | negtrack=NULL; |
| 540 | continue; |
| 541 | } |
| 542 | |
| 543 | io2Prong = Make2Prong(twoTrackArray,source,vertexp1n1,dcap1n1); |
| 544 | |
| 545 | if(io2Prong->CosPointingAngle()<fCosMin)continue; |
| 546 | |
| be4df180 |
547 | // cout<<"CPA: "<<io2Prong->CosPointingAngle()<<endl; |
| 930d0706 |
548 | |
| 549 | // AliAODTrack *trk0 = (AliAODTrack*)io2Prong->GetDaughter(0); |
| 550 | // AliAODTrack *trk1 = (AliAODTrack*)io2Prong->GetDaughter(1); |
| 551 | |
| 552 | // cout<<"**********************************************"<<endl; |
| 553 | // cout<<trk0/*->GetID()*/<<" "<<negtrackAOD->GetID()<<endl; |
| 554 | // cout<<trk1/*->GetID()*/<<" "<<postrackAOD->GetID()<<endl; |
| 555 | // cout<<"d0 io2Prong: "<<io2Prong->GetProngID(1)<<endl; |
| 556 | // cout<<"d1 io2Prong: "<<io2Prong->GetProngID(0)<<endl; |
| 557 | // cout<<"**********************************************"<<endl; |
| 558 | |
| 559 | // rd = new((*fSecondaryVerices)[nsv++]) AliAODRecoDecayLF2Prong(*io2Prong); |
| be4df180 |
560 | |
| 561 | // cout<<"is ok??? "<<isOk<<endl; |
| 930d0706 |
562 | if(isOk){ |
| 563 | new((*fSecondaryVerices)[nsv++]) AliAODRecoDecayLF2Prong(*io2Prong); |
| 564 | new((*fDaughterTracks)[ntracksD++]) AliAODTrack(*negtrackAOD); |
| 565 | new((*fDaughterTracks)[ntracksD++]) AliAODTrack(*postrackAOD); |
| 566 | |
| 567 | } |
| 568 | // rd->SetSecondaryVtx(vertexp1n1); |
| 569 | // vertexp1n1->SetParent(rd); |
| 570 | // AddRefs(vertexp1n1,rd,source,twoTrackArray); |
| 571 | |
| 572 | delete negtrack; |
| 573 | negtrack = NULL; |
| 574 | |
| 575 | delete vertexp1n1; |
| 576 | vertexp1n1= NULL; |
| 577 | continue; |
| 578 | // } |
| 579 | } |
| 580 | |
| 581 | delete postrack; |
| 582 | postrack = NULL; |
| 583 | |
| 584 | } |
| 585 | |
| 586 | //---------------------------------------------------------- |
| 587 | |
| 588 | assert(fVertices!=0x0); |
| 589 | fVertices->Clear("C"); |
| 590 | TIter nextV(source.GetVertices()); |
| 591 | AliAODVertex* v; |
| 592 | Int_t nvertices(0); |
| 593 | |
| 594 | while ( ( v = static_cast<AliAODVertex*>(nextV()) ) ) |
| 595 | { |
| 596 | if ( v->GetType() == AliAODVertex::kPrimary || |
| 597 | v->GetType() == AliAODVertex::kMainSPD || |
| 598 | v->GetType() == AliAODVertex::kPileupSPD || |
| 599 | v->GetType() == AliAODVertex::kPileupTracks|| |
| 600 | v->GetType() == AliAODVertex::kMainTPC ) |
| 601 | { |
| 602 | |
| 603 | AliAODVertex* tmp = v->CloneWithoutRefs(); |
| 604 | //AliAODVertex* copiedVertex = new((*fVertices)[nvertices++]) AliAODVertex(*tmp); |
| 605 | new((*fVertices)[nvertices++]) AliAODVertex(*tmp); |
| 606 | |
| 607 | // to insure the main vertex retains the ncontributors information |
| 608 | // (which is otherwise computed dynamically from |
| 609 | // references to tracks, which we do not keep in muon aods...) |
| 610 | // we set it here |
| 611 | |
| 612 | //copiedVertex->SetNContributors(v->GetNContributors()); |
| 613 | |
| 614 | // fVertices->Delete(); |
| 615 | // delete copiedVertex; |
| 616 | delete tmp; |
| 617 | // printf("....Prendo il vertice primario...\n"); |
| 618 | } |
| 619 | // printf("....Prendo il vertice primario...\n"); |
| 620 | } |
| 621 | |
| 622 | //printf("....Done NuclEx Replicator...\n"); |
| 623 | |
| 624 | AliDebug(1,Form("input mu tracks=%d tracks=%d vertices=%d nnuclei=%d", |
| 625 | input,fSecondaryVerices->GetEntries(),fVertices->GetEntries(),fNuclei->GetEntries())); |
| 626 | |
| 627 | // cout<<"Delete..."<<endl; |
| 628 | // delete foPion; |
| 629 | // foPion = NULL; |
| 630 | //cout<<"Delete 1"<< endl; |
| 631 | |
| 632 | if(io2Prong) {delete io2Prong; io2Prong=NULL;} |
| 633 | //cout<<"Delete 2"<< endl; |
| 634 | twoTrackArray->Delete(); delete twoTrackArray; |
| 635 | // cout<<"Delete 3"<< endl; |
| 636 | // vtx->Delete(); delete vtx; |
| 637 | // cout<<"Delete 4"<< endl; |
| 638 | if(fV1) { delete fV1; fV1=NULL; } |
| 639 | // cout<<"Delete 5"<< endl; |
| 640 | if(fAODMap) { delete [] fAODMap; fAODMap=NULL; } |
| 641 | delete []indices; |
| 642 | // cout<<"Delete 6"<< endl; |
| 643 | delete fVertexerTracks; |
| 644 | |
| 645 | // cout<<"Mi rompo alla fine. Perche???"<<endl; |
| 646 | |
| 647 | } |
| 648 | |
| 649 | //----------------------------------------------------------------------------- |
| 650 | |
| 651 | AliAODVertex *AliAODMCNuclExReplicator::ReconstructSecondaryVertex(TObjArray *trkArray, |
| 652 | Double_t &dispersion,Bool_t useTRefArray) const |
| 653 | |
| 654 | { |
| 655 | // Secondary vertex reconstruction with AliVertexerTracks or AliKFParticle |
| 656 | //AliCodeTimerAuto("",0); |
| 657 | |
| 658 | AliESDVertex *vertexESD = 0; |
| 659 | AliAODVertex *vertexAOD = 0; |
| 660 | |
| 661 | if(!fSecVtxWithKF) { // AliVertexerTracks |
| 662 | |
| 663 | fVertexerTracks->SetVtxStart(fV1); |
| 664 | vertexESD = (AliESDVertex*)fVertexerTracks->VertexForSelectedESDTracks(trkArray); |
| 665 | |
| 666 | if(!vertexESD) return vertexAOD; |
| 667 | |
| 668 | |
| e690d4d0 |
669 | Double_t vertRadius2=vertexESD->GetX()*vertexESD->GetX()+vertexESD->GetY()*vertexESD->GetY(); |
| 930d0706 |
670 | if(vertRadius2>200.){ |
| 671 | // vertex outside beam pipe, reject candidate to avoid propagation through material |
| 672 | delete vertexESD; vertexESD=NULL; |
| 673 | return vertexAOD; |
| 674 | } |
| 675 | |
| 676 | } else { // Kalman Filter vertexer (AliKFParticle) |
| 677 | |
| 678 | AliKFParticle::SetField(fBzkG); |
| 679 | |
| 680 | AliKFVertex vertexKF; |
| 681 | |
| 682 | Int_t nTrks = trkArray->GetEntriesFast(); |
| 683 | for(Int_t i=0; i<nTrks; i++) { |
| 684 | AliESDtrack *esdTrack = (AliESDtrack*)trkArray->At(i); |
| 685 | AliKFParticle daughterKF(*esdTrack,211); |
| 686 | vertexKF.AddDaughter(daughterKF); |
| 687 | } |
| 688 | vertexESD = new AliESDVertex(vertexKF.Parameters(), |
| 689 | vertexKF.CovarianceMatrix(), |
| 690 | vertexKF.GetChi2(), |
| 691 | vertexKF.GetNContributors()); |
| 692 | |
| 693 | } |
| 694 | // convert to AliAODVertex |
| 695 | Double_t pos[3],cov[6],chi2perNDF; |
| 696 | vertexESD->GetXYZ(pos); // position |
| 697 | vertexESD->GetCovMatrix(cov); //covariance matrix |
| 698 | chi2perNDF = vertexESD->GetChi2toNDF(); |
| 699 | dispersion = vertexESD->GetDispersion(); |
| 700 | delete vertexESD; |
| 701 | vertexESD=NULL; |
| 702 | |
| 703 | Int_t nprongs= (useTRefArray ? 0 : trkArray->GetEntriesFast()); |
| 704 | vertexAOD = new AliAODVertex(pos,cov,chi2perNDF,0x0,-1,AliAODVertex::kUndef,nprongs); |
| 705 | |
| 706 | // cout<<"------------------> Reconstruct vertexAOD: "<<vertexAOD<<endl; |
| 707 | |
| 708 | return vertexAOD; |
| 709 | |
| 710 | |
| 711 | } |
| 712 | |
| 713 | //----------------------------------------------------------------------------- |
| 714 | |
| 715 | AliAODRecoDecayLF2Prong* AliAODMCNuclExReplicator::Make2Prong(TObjArray *twoTrackArray,const AliAODEvent &evento, |
| 716 | AliAODVertex *secVert,Double_t dca) |
| 717 | |
| 718 | |
| 719 | { |
| 720 | |
| 721 | //cout<<"Inside make 2 prong"<<endl; |
| 722 | |
| 723 | Int_t charge[2]; |
| 724 | charge[0]=1; //it was -1 //Ramona |
| 725 | charge[1]=2; |
| 726 | |
| 727 | // From AliAnalysisVertexingLF.cxx Method:Make2Prongs |
| 728 | |
| 729 | fBzkG = evento.GetMagneticField(); |
| 730 | |
| 731 | Double_t px[2],py[2],pz[2],d0[2],d0err[2]; |
| 732 | // Also this has been changed //Ramona |
| 733 | AliESDtrack *negtrack = (AliESDtrack*)twoTrackArray->UncheckedAt(0); |
| 734 | AliESDtrack *postrack = (AliESDtrack*)twoTrackArray->UncheckedAt(1); |
| 735 | |
| 736 | //cout<<"negtrack: "<<negtrack<<" postrack: "<<postrack<<endl; |
| 737 | //cout<<"kVeryBig: "<<kVeryBig<<endl; |
| 738 | //cout<<"secVert: "<<secVert<<endl; |
| 739 | |
| 740 | // // propagate tracks to secondary vertex, to compute inv. mass |
| 741 | |
| 742 | postrack->PropagateToDCA(secVert,fBzkG,kVeryBig); |
| 743 | negtrack->PropagateToDCA(secVert,fBzkG,kVeryBig); |
| 744 | |
| 745 | Double_t momentum[3]; |
| 746 | |
| 747 | negtrack->GetPxPyPz(momentum); |
| 748 | px[0] = charge[0]*momentum[0]; |
| 749 | py[0] = charge[0]*momentum[1]; |
| 750 | pz[0] = charge[0]*momentum[2]; |
| 751 | |
| 752 | postrack->GetPxPyPz(momentum); |
| 753 | px[1] = charge[1]*momentum[0]; |
| 754 | py[1] = charge[1]*momentum[1]; |
| 755 | pz[1] = charge[1]*momentum[2]; |
| 756 | |
| 757 | //cout<< px[0] <<" "<< " "<< py[0] << " "<< pz[0]<<endl; |
| 758 | // px[0] = momentum[0]; py[0] = momentum[1]; pz[0] = momentum[2]; |
| 759 | //cout<< px[1] <<" "<< " "<< py[1] << " "<< pz[1]<<endl; |
| 760 | //px[1] = momentum[0]; py[1] = momentum[1]; pz[1] = momentum[2]; |
| 761 | |
| 762 | |
| 763 | // primary vertex to be used by this candidate |
| 764 | AliAODVertex *primVertexAOD = evento.GetPrimaryVertex(); |
| 765 | //cout<<"primVertexAOD "<<primVertexAOD<<endl; |
| 766 | if(!primVertexAOD) return 0x0; |
| 767 | |
| 768 | Double_t d0z0[2],covd0z0[3]; |
| 769 | |
| 770 | d0z0[0] = -999.; |
| 771 | d0z0[1] = -999.; |
| 772 | covd0z0[0] = -999.; |
| 773 | covd0z0[1] = -999.; |
| 774 | covd0z0[2] = -999.; |
| 775 | |
| 776 | d0[0] = d0z0[0]; |
| 777 | d0err[0] = TMath::Sqrt(TMath::Abs(covd0z0[0])); |
| 778 | |
| 779 | d0[1] = d0z0[0]; |
| 780 | d0err[1] = TMath::Sqrt(TMath::Abs(covd0z0[0])); |
| 781 | |
| 782 | // create the object AliAODRecoDecayLF2Prong |
| 783 | // AliAODRecoDecayLF2Prong *the2Prong = new AliAODRecoDecayLF2Prong(secVert,px,py,pz,d0,d0err,dcap1n1); |
| 784 | AliAODRecoDecayLF2Prong *the2Prong = new AliAODRecoDecayLF2Prong(secVert,px,py,pz,d0,d0err,dca); |
| 785 | the2Prong->SetOwnPrimaryVtx(primVertexAOD); |
| 786 | |
| 787 | // the2Prong->SetProngIDs(2,{-999,999}); |
| 788 | // UShort_t id0[2]={99999,999999}; |
| 789 | // the2Prong->SetProngIDs(2,id0); |
| 790 | |
| 791 | UShort_t id[2]={(UShort_t)postrack->GetID(),(UShort_t)negtrack->GetID()}; |
| 792 | the2Prong->SetProngIDs(2,id); |
| 793 | |
| 794 | // cout<<"\n\n\nMake 2 Prong: id[0]"<<id[0]<<" id[1]: "<<id[1]<<endl; |
| 795 | // cout<<"Get: 1 "<<the2Prong->GetProngID(0)<<" 2 "<<the2Prong->GetProngID(1)<<endl; |
| 796 | // cout<<"Get: 3 "<<the2Prong->GetProngID(2)<<" 4 "<<the2Prong->GetProngID(3)<<endl<<endl<<endl; |
| 797 | //delete primVertexAOD; primVertexAOD=NULL; |
| 798 | |
| 799 | if(postrack->Charge()!=0 && negtrack->Charge()!=0) { |
| 800 | |
| 801 | AddDaughterRefs(secVert,(AliAODEvent&)evento,twoTrackArray); |
| 802 | // AddDaughterRefs(secVert,(AliAODEvent*)evento,twoTrackArray); |
| 803 | |
| 804 | } |
| 805 | |
| 806 | return the2Prong; |
| 807 | |
| 808 | delete the2Prong; |
| 809 | } |
| 810 | |
| 811 | //---------------------------------------------------------------------------- |
| 812 | void AliAODMCNuclExReplicator::AddDaughterRefs(AliAODVertex *v, |
| 813 | const AliAODEvent &event, |
| 814 | const TObjArray *trkArray) const |
| 815 | |
| 816 | { |
| 817 | // Add the AOD tracks as daughters of the vertex (TRef) |
| 818 | // AliCodeTimerAuto("",0); |
| 819 | // cout<<"Inside AddDaughterRefs "<<endl; |
| 820 | |
| 821 | Int_t nDg = v->GetNDaughters(); |
| 822 | |
| 823 | TObject *dg = 0; |
| 824 | if(nDg) dg = v->GetDaughter(0); |
| 825 | if(dg) return; // daughters already added |
| 826 | |
| 827 | Int_t nTrks = trkArray->GetEntriesFast(); |
| 828 | |
| 829 | AliExternalTrackParam *track = 0; |
| 830 | AliAODTrack *aodTrack = 0; |
| 831 | Int_t id; |
| 832 | |
| 833 | for(Int_t i=0; i<nTrks; i++) { |
| 834 | track = (AliExternalTrackParam*)trkArray->UncheckedAt(i); |
| 835 | |
| 836 | id = (Int_t)track->GetID(); |
| 837 | // printf("---> %d\n",id); |
| 838 | if(id<0) continue; // this track is a AliAODRecoDecay |
| 839 | |
| 840 | aodTrack = (AliAODTrack*)event.GetTrack(fAODMap[id]); |
| 841 | v->AddDaughter(aodTrack); |
| 842 | } |
| 843 | return; |
| 844 | } |
| 845 | //---------------------------------------------------------------------------- |
| 846 | |
| 847 | void AliAODMCNuclExReplicator::AddRefs(AliAODVertex *v,AliAODRecoDecayLF *rd, |
| 848 | const AliAODEvent &event, |
| 849 | const TObjArray *trkArray) const |
| 850 | |
| 851 | { |
| 852 | // Add the AOD tracks as daughters of the vertex (TRef) |
| 853 | // Also add the references to the primary vertex and to the cuts |
| 854 | //AliCodeTimerAuto("",0); |
| 855 | |
| 856 | AddDaughterRefs(v,event,trkArray); |
| 857 | rd->SetPrimaryVtxRef((AliAODVertex*)event.GetPrimaryVertex()); |
| 858 | return; |
| 859 | } |
| 860 | |
| 861 | //--------------------------------------------------------------------------- |
| 862 | |
| 863 | void AliAODMCNuclExReplicator::Terminate(){ |
| 864 | |
| 865 | } |
git://git.uio.no / u / mrichter / AliRoot.git / blame