VZERO/AliVZERO.cxx

   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 /* $Id$ */
  17
  18 ///////////////////////////////////////////////////////////////////////////
  19 //                                                                       //
  20 //                          V-Zero   Detector                            //
  21 //  This class contains the base procedures for the VZERO  detector      //
  22 //  Default geometry of November 2003 :   V0R box is 4.4 cm thick        //
  23 //                                  scintillators are 2 cm thick         //
  24 //  All comments should be sent to Brigitte CHEYNIS :                    //
  25 //                                 b.cheynis@ipnl.in2p3.fr               //
  26 //                                                                       //
  27 //                                                                       //
  28 ///////////////////////////////////////////////////////////////////////////
  29
  30
  31 // --- Standard libraries ---
  32 #include <Riostream.h>
  33 #include <stdlib.h>
  34
  35 // --- ROOT libraries ---
  36 #include <TNamed.h>
  37 #include "TROOT.h"
  38 #include "TFile.h"
  39 #include "TNetFile.h"
  40 #include "TRandom.h"
  41 #include "TTree.h"
  42 #include "TBranch.h"
  43 #include "TClonesArray.h"
  44 #include "TStopwatch.h"
  45 #include "TParameter.h"
  46 #include "TF1.h"
  47
  48 // --- AliRoot header files ---
  49 #include "AliRun.h"
  50 #include "AliMC.h"
  51 #include "AliVZERO.h"
  52 #include "AliVZEROLoader.h"
  53 #include "AliVZERODigitizer.h"
  54 #include "AliVZEROBuffer.h"
  55 #include "AliDigitizationInput.h"
  56 #include "AliVZEROdigit.h"
  57 #include "AliVZEROSDigit.h"
  58 #include "AliDAQ.h"
  59 #include "AliRawReader.h"
  60 #include "AliCDBManager.h"
  61 #include "AliCDBEntry.h"
  62 #include "AliVZERORawStream.h"
  63 #include "AliVZEROCalibData.h"
  64 #include "AliVZERORecoParam.h"
  65 #include "AliVZEROReconstructor.h"
  66
  67 ClassImp(AliVZERO)
  68  //__________________________________________________________________
  69 AliVZERO::AliVZERO(): AliDetector(),
  70           fIdSens1(0),
  71           fThickness(0.),
  72           fThickness1(0.),
  73           fMaxStepQua(0.),
  74           fMaxStepAlu(0.),
  75           fMaxDestepQua(0.),
  76           fMaxDestepAlu(0.),
  77           fCalibData(NULL),
  78           fTimeSlewing(NULL),
  79           fSignalShape(NULL),
  80           fRecoParam(NULL)
  81 {
  82 /// Default Constructor
  83
  84     AliDebug(1,Form("default (empty) ctor this=%p",this));
  85     fIshunt          = 0;
  86
  87     for(Int_t i = 0 ; i < 64; ++i) {
  88       fNBins[i] = 0;
  89       fBinSize[i] = 0;
  90     }
  91 }
  92
  93 //_____________________________________________________________________________
  94 AliVZERO::AliVZERO(const char *name, const char *title)
  95        : AliDetector(name,title),
  96          fIdSens1(0),
  97          fThickness(4.4),
  98          fThickness1(2.0),
  99          fMaxStepQua(0.05),
 100          fMaxStepAlu(0.01),
 101          fMaxDestepQua(-1.0),
 102          fMaxDestepAlu(-1.0),
 103          fCalibData(NULL),
 104          fTimeSlewing(NULL),
 105          fSignalShape(NULL),
 106          fRecoParam(NULL)
 107 {
 108
 109   // Standard constructor for VZERO Detector
 110
 111   AliDebug(1,Form("ctor this=%p",this));
 112
 113   //  fIshunt       =  1;  // All hits are associated with primary particles
 114
 115   fHits         =  new TClonesArray("AliVZEROhit", 400);
 116   fDigits       =  new TClonesArray("AliVZEROdigit",400);
 117
 118   gAlice->GetMCApp()->AddHitList(fHits);
 119
 120 //   fThickness    =  4.4;   // total thickness of the V0R box in cm
 121 //   fThickness1   =  2.0;   // thickness of scintillating cells in cm
 122 //
 123 //   fMaxStepQua   =  0.05;
 124 //   fMaxStepAlu   =  0.01;
 125 //
 126 //   fMaxDestepQua =  -1.0;
 127 //   fMaxDestepAlu =  -1.0;
 128
 129   for(Int_t i = 0 ; i < 64; ++i) {
 130     fNBins[i] = 0;
 131     fBinSize[i] = 0;
 132   }
 133 }
 134
 135 //_____________________________________________________________________________
 136 AliVZERO::~AliVZERO()
 137 {
 138   //
 139   // Default destructor for VZERO Detector
 140   //
 141
 142     if (fHits) {
 143         fHits->Delete();
 144         delete fHits;
 145         fHits=0; }
 146
 147     if (fDigits) {
 148         fDigits->Delete();
 149         delete fDigits;
 150         fDigits=0; }
 151     if (fSignalShape) {
 152       delete fSignalShape;
 153       fSignalShape = NULL;
 154     }
 155     if (fRecoParam) {
 156       delete fRecoParam;
 157       fRecoParam = NULL;
 158     }
 159 }
 160
 161 //_____________________________________________________________________________
 162 void AliVZERO::CreateGeometry()
 163 {
 164   //
 165   // Builds simple Geant3 geometry
 166   //
 167 }
 168 //_____________________________________________________________________________
 169 void AliVZERO::CreateMaterials()
 170 {
 171   //
 172   // Creates materials used for Geant3 geometry
 173   //
 174 }
 175
 176 //_____________________________________________________________________________
 177 void AliVZERO::Init()
 178 {
 179   //
 180   // Initialises the VZERO  class after it has been built
 181   //
 182 }
 183
 184
 185 //_____________________________________________________________________________
 186 void AliVZERO::SetMaxStepQua(Float_t p1)
 187 {
 188   //
 189   // Possible parametrisation of steps in active materials
 190   //
 191      fMaxStepQua = p1;
 192 }
 193
 194 //_____________________________________________________________________________
 195 void AliVZERO::SetMaxStepAlu(Float_t p1)
 196 {
 197   //
 198   // Possible parametrisation of steps in Aluminum foils (not used in
 199   // version v2)
 200   //
 201     fMaxStepAlu = p1;
 202 }
 203
 204 //_____________________________________________________________________________
 205 void AliVZERO::SetMaxDestepQua(Float_t p1)
 206 {
 207   //
 208   // Possible parametrisation of steps in active materials (quartz)
 209   //
 210     fMaxDestepQua = p1;
 211 }
 212
 213 //_____________________________________________________________________________
 214 void AliVZERO::SetMaxDestepAlu(Float_t p1)
 215 {
 216   //
 217   // Possible parametrisation of steps in Aluminum (not used in
 218   // version v2)
 219   //
 220     fMaxDestepAlu = p1;
 221 }
 222
 223 //_____________________________________________________________________________
 224 AliLoader* AliVZERO::MakeLoader(const char* topfoldername)
 225 {
 226   //
 227   // Builds VZEROgetter (AliLoader type)
 228   // if detector wants to use customized getter, it must overload this method
 229   //
 230 //  Info("MakeLoader","Creating AliVZEROLoader. Top folder is %s.",topfoldername);
 231
 232   AliDebug(1,Form("Creating AliVZEROLoader, Top folder is %s ",topfoldername));
 233   fLoader = new AliVZEROLoader(GetName(),topfoldername);
 234   return fLoader;
 235 }
 236
 237 //_____________________________________________________________________________
 238 void AliVZERO::SetTreeAddress()
 239 {
 240   //
 241   // Sets tree address for hits.
 242   //
 243   if (fLoader->TreeH() && (fHits == 0x0))
 244     fHits = new  TClonesArray("AliVZEROhit", 400);
 245
 246   AliDetector::SetTreeAddress();
 247 }
 248
 249 //_____________________________________________________________________________
 250 AliDigitizer* AliVZERO::CreateDigitizer(AliDigitizationInput* digInput) const
 251 {
 252   //
 253   // Creates a digitizer for VZERO
 254   //
 255   return new AliVZERODigitizer(digInput);
 256 }
 257
 258 //_____________________________________________________________________________
 259 void AliVZERO::Hits2Digits(){
 260   //
 261   // Converts hits to digits
 262   //
 263   // Creates the VZERO digitizer
 264   AliVZERODigitizer* dig = new AliVZERODigitizer(this,AliVZERODigitizer::kHits2Digits);
 265
 266   // Creates the digits
 267   dig->Digitize("");
 268
 269   // deletes the digitizer
 270   delete dig;
 271 }
 272
 273 //_____________________________________________________________________________
 274 void AliVZERO::Hits2SDigits(){
 275   //
 276   // Converts hits to summable digits
 277   //
 278   // Creates the VZERO digitizer
 279   AliVZERODigitizer* dig = new AliVZERODigitizer(this,AliVZERODigitizer::kHits2SDigits);
 280
 281   // Creates the sdigits
 282   dig->Digitize("");
 283
 284   // deletes the digitizer
 285   delete dig;
 286 }
 287
 288 //_____________________________________________________________________________
 289 void AliVZERO::Digits2Raw()
 290 {
 291    //
 292    //  Converts digits of the current event to raw data
 293    //
 294
 295    AliVZERO *fVZERO = (AliVZERO*)gAlice->GetDetector("VZERO");
 296    fLoader->LoadDigits();
 297    TTree* digits = fLoader->TreeD();
 298    if (!digits) {
 299       Error("Digits2Raw", "no digits tree");
 300       return;
 301    }
 302    TClonesArray * VZEROdigits = new TClonesArray("AliVZEROdigit",1000);
 303    fVZERO->SetTreeAddress();
 304    digits->GetBranch("VZERODigit")->SetAddress(&VZEROdigits);
 305
 306    const char *fileName    = AliDAQ::DdlFileName("VZERO",0);
 307    AliVZEROBuffer* buffer  = new AliVZEROBuffer(fileName);
 308
 309    // Get Trigger information first
 310    // Read trigger inputs from trigger-detector object
 311    AliDataLoader * dataLoader = fLoader->GetDigitsDataLoader();
 312    if( !dataLoader->IsFileOpen() )
 313         dataLoader->OpenFile( "READ" );
 314    AliTriggerDetector* trgdet = (AliTriggerDetector*)dataLoader->GetDirectory()->Get( "Trigger" );
 315    UInt_t triggerInfo = 0;
 316    if(trgdet) {
 317       triggerInfo = trgdet->GetMask() & 0xffff;
 318    }
 319    else {
 320       AliError(Form("There is no trigger object for %s",fLoader->GetName()));
 321    }
 322
 323    buffer->WriteTriggerInfo((UInt_t)triggerInfo);
 324    buffer->WriteTriggerScalers();
 325    buffer->WriteBunchNumbers();
 326
 327    Int_t aBBflagsV0A = 0;
 328    Int_t aBBflagsV0C = 0;
 329    Int_t aBGflagsV0A = 0;
 330    Int_t aBGflagsV0C = 0;
 331
 332    if (digits->GetUserInfo()->FindObject("BBflagsV0A")) {
 333      aBBflagsV0A = ((TParameter<int>*)digits->GetUserInfo()->FindObject("BBflagsV0A"))->GetVal();
 334    }
 335    else
 336      AliWarning("V0A beam-beam flags not found in digits tree UserInfo! The flags will not be written to the raw-data stream!");
 337
 338    if (digits->GetUserInfo()->FindObject("BBflagsV0C")) {
 339      aBBflagsV0C = ((TParameter<int>*)digits->GetUserInfo()->FindObject("BBflagsV0C"))->GetVal();
 340    }
 341    else
 342      AliWarning("V0C beam-beam flags not found in digits tree UserInfo! The flags will not be written to the raw-data stream!");
 343
 344    if (digits->GetUserInfo()->FindObject("BGflagsV0A")) {
 345      aBGflagsV0A = ((TParameter<int>*)digits->GetUserInfo()->FindObject("BGflagsV0A"))->GetVal();
 346    }
 347    else
 348      AliWarning("V0A beam-gas flags not found in digits tree UserInfo! The flags will not be written to the raw-data stream!");
 349
 350    if (digits->GetUserInfo()->FindObject("BGflagsV0C")) {
 351      aBGflagsV0C = ((TParameter<int>*)digits->GetUserInfo()->FindObject("BGflagsV0C"))->GetVal();
 352    }
 353    else
 354      AliWarning("V0C beam-gas flags not found in digits tree UserInfo! The flags will not be written to the raw-data stream!");
 355
 356    // Now retrieve the channel information: charge smaples + time and
 357    // dump it into ADC and Time arrays
 358    Int_t nEntries = Int_t(digits->GetEntries());
 359    Short_t aADC[64][AliVZEROdigit::kNClocks];
 360    Float_t aTime[64];
 361    Float_t aWidth[64];
 362    Bool_t  aIntegrator[64];
 363    Bool_t  aBBflag[64];
 364    Bool_t  aBGflag[64];
 365
 366    for (Int_t i = 0; i < nEntries; i++) {
 367      fVZERO->ResetDigits();
 368      digits->GetEvent(i);
 369      Int_t ndig = VZEROdigits->GetEntriesFast();
 370
 371      if(ndig == 0) continue;
 372      for(Int_t k=0; k<ndig; k++){
 373          AliVZEROdigit* fVZERODigit = (AliVZEROdigit*) VZEROdigits->At(k);
 374          // Convert aliroot channel k into FEE channel iChannel before writing data
 375          Int_t iChannel       = AliVZEROCalibData::GetBoardNumber(fVZERODigit->PMNumber()) * 8 +
 376            AliVZEROCalibData::GetFEEChannelNumber(fVZERODigit->PMNumber());
 377          for(Int_t iClock = 0; iClock < AliVZEROdigit::kNClocks; ++iClock) aADC[iChannel][iClock] = fVZERODigit->ChargeADC(iClock);
 378          aTime[iChannel]      = fVZERODigit->Time();
 379          aWidth[iChannel]     = fVZERODigit->Width();
 380          aIntegrator[iChannel]= fVZERODigit->Integrator();
 381          if(fVZERODigit->PMNumber() < 32) {
 382            aBBflag[iChannel] = (aBBflagsV0C >> fVZERODigit->PMNumber()) & 0x1;
 383            aBGflag[iChannel] = (aBGflagsV0C >> fVZERODigit->PMNumber()) & 0x1;
 384          }
 385          else {
 386            aBBflag[iChannel] = (aBBflagsV0A >> (fVZERODigit->PMNumber()-32)) & 0x1;
 387            aBGflag[iChannel] = (aBGflagsV0A >> (fVZERODigit->PMNumber()-32)) & 0x1;
 388          }
 389          AliDebug(1,Form("DDL: %s\tdigit number: %d\tPM number: %d\tADC: %d\tTime: %f",
 390                          fileName,k,fVZERODigit->PMNumber(),aADC[iChannel][AliVZEROdigit::kNClocks/2],aTime[iChannel]));
 391      }
 392    }
 393
 394    // Now fill raw data
 395
 396    for (Int_t  iCIU = 0; iCIU < 8; iCIU++) {
 397
 398    // decoding of one Channel Interface Unit numbered iCIU - there are 8 channels per CIU (and 8 CIUs) :
 399
 400       for(Int_t iChannel_Offset = iCIU*8; iChannel_Offset < (iCIU*8)+8; iChannel_Offset=iChannel_Offset+4) {
 401          for(Int_t iChannel = iChannel_Offset; iChannel < iChannel_Offset+4; iChannel++) {
 402              buffer->WriteChannel(iChannel, aADC[iChannel], aIntegrator[iChannel]);
 403          }
 404          buffer->WriteBeamFlags(&aBBflag[iChannel_Offset],&aBGflag[iChannel_Offset]);
 405          buffer->WriteMBInfo();
 406          buffer->WriteMBFlags();
 407          buffer->WriteBeamScalers();
 408       }
 409
 410       for(Int_t iChannel = iCIU*8 + 7; iChannel >= iCIU*8; iChannel--) {
 411           buffer->WriteTiming(aTime[iChannel], aWidth[iChannel]);
 412       }
 413
 414     // End of decoding of one CIU card
 415
 416   } // end of decoding the eight CIUs
 417
 418   delete buffer;
 419   fLoader->UnloadDigits();
 420 }
 421
 422 //_____________________________________________________________________________
 423 Bool_t AliVZERO::Raw2SDigits(AliRawReader* rawReader){
 424   // Converts the VZERO raw data into digits
 425   // The method is used for merging simulated and
 426   // real data events
 427   TStopwatch timer;
 428   timer.Start();
 429
 430   if(!fLoader) {
 431     AliError("no VZERO loader found");
 432     return kFALSE;
 433   }
 434   fLoader->LoadSDigits("UPDATE");
 435
 436   if (!fLoader->TreeS()) fLoader->MakeTree("S");
 437   fLoader->MakeSDigitsContainer();
 438   TTree* treeS  = fLoader->TreeS();
 439
 440   TClonesArray *sdigits = new TClonesArray("AliVZEROSDigit", 64);
 441   treeS->Branch("VZEROSDigit", &sdigits);
 442
 443   {
 444     rawReader->Reset();
 445     AliVZERORawStream rawStream(rawReader);
 446
 447     if (!rawStream.Next()) return kFALSE; // No VZERO data found
 448
 449     GetCalibData();
 450
 451     Int_t nSDigits = 0;
 452     Float_t *charges = NULL;
 453     Int_t nbins = 0;
 454     for(Int_t iChannel=0; iChannel < 64; ++iChannel) {
 455       Int_t offlineCh = rawStream.GetOfflineChannel(iChannel);
 456       Short_t chargeADC[AliVZEROdigit::kNClocks];
 457       for(Int_t iClock=0; iClock < AliVZEROdigit::kNClocks; ++iClock) {
 458         chargeADC[iClock] = rawStream.GetPedestal(iChannel,iClock);
 459       }
 460       // Integrator flag
 461       Bool_t integrator = rawStream.GetIntegratorFlag(iChannel,AliVZEROdigit::kNClocks/2);
 462       // HPTDC data (leading time and width)
 463       Int_t board = AliVZEROCalibData::GetBoardNumber(offlineCh);
 464       Float_t time = rawStream.GetTime(iChannel)*fCalibData->GetTimeResolution(board);
 465       //      Float_t width = rawStream.GetWidth(iChannel)*fCalibData->GetWidthResolution(board);
 466       Float_t adc = 0;
 467
 468       // Pedestal retrieval and suppression
 469       Float_t maxadc = 0;
 470       Int_t imax = -1;
 471       Float_t adcPedSub[AliVZEROdigit::kNClocks];
 472       Float_t integral = fSignalShape->Integral(0,200);
 473       for(Int_t iClock=0; iClock < AliVZEROdigit::kNClocks; ++iClock) {
 474         Bool_t iIntegrator = (iClock%2 == 0) ? integrator : !integrator;
 475         Int_t k = offlineCh + 64*iIntegrator;
 476         adcPedSub[iClock] = (Float_t)chargeADC[iClock] - fCalibData->GetPedestal(k);
 477         if(adcPedSub[iClock] <= fRecoParam->GetNSigmaPed()*fCalibData->GetSigma(k)) {
 478           adcPedSub[iClock] = 0;
 479           continue;
 480         }
 481         if(iClock < fRecoParam->GetStartClock() || iClock > fRecoParam->GetEndClock()) continue;
 482         if(adcPedSub[iClock] > maxadc) {
 483           maxadc = adcPedSub[iClock];
 484           imax   = iClock;
 485         }
 486       }
 487       if (imax != -1) {
 488         Int_t start = imax - fRecoParam->GetNPreClocks();
 489         if (start < 0) start = 0;
 490         Int_t end = imax + fRecoParam->GetNPostClocks();
 491         if (end > 20) end = 20;
 492         for(Int_t iClock = start; iClock <= end; iClock++) {
 493           adc += adcPedSub[iClock];
 494         }
 495       }
 496       Float_t correctedTime = CorrectLeadingTime(offlineCh,time,adc);
 497
 498       if (!charges) {
 499         nbins = fNBins[offlineCh];
 500         charges = new Float_t[nbins];
 501       }
 502       else if (nbins != fNBins[offlineCh]) {
 503         delete [] charges;
 504         nbins = fNBins[offlineCh];
 505         charges = new Float_t[nbins];
 506       }
 507       memset(charges,0,nbins*sizeof(Float_t));
 508
 509       // Now lets produce SDigit
 510       if ((correctedTime > (AliVZEROReconstructor::kInvalidTime + 1e-6)) &&
 511           (adc > 1e-6)) {
 512         for(Int_t iBin = 0; iBin < nbins; ++iBin) {
 513           Float_t t = fBinSize[offlineCh]*Float_t(iBin);
 514           if ((t < correctedTime) ||
 515               (t > (correctedTime+200.))) continue;
 516           charges[iBin] = kChargePerADC*adc*(fSignalShape->Eval(t-correctedTime)*fBinSize[offlineCh]/integral);
 517         }
 518       }
 519
 520       TClonesArray &sdigitsref = *sdigits;
 521       new (sdigitsref[nSDigits++]) AliVZEROSDigit(offlineCh,fNBins[offlineCh],charges);
 522     }
 523     if (charges) delete [] charges;
 524   }
 525
 526   treeS->Fill();
 527   fLoader->WriteSDigits("OVERWRITE");
 528   fLoader->UnloadSDigits();
 529
 530   timer.Stop();
 531   timer.Print();
 532   return kTRUE;
 533 }
 534
 535 //_____________________________________________________________________________
 536 void AliVZERO::GetCalibData()
 537 {
 538   // Gets calibration object for VZERO set
 539   // Do nothing in case it is already loaded
 540   if (fCalibData) return;
 541
 542   AliCDBEntry *entry = AliCDBManager::Instance()->Get("VZERO/Calib/Data");
 543   if (entry) fCalibData = (AliVZEROCalibData*) entry->GetObject();
 544   if (!fCalibData)  AliFatal("No calibration data from calibration database !");
 545
 546   AliCDBEntry *entry2 = AliCDBManager::Instance()->Get("VZERO/Calib/TimeSlewing");
 547   if (!entry2) AliFatal("VZERO time slewing function is not found in OCDB !");
 548   fTimeSlewing = (TF1*)entry2->GetObject();
 549
 550   for(Int_t i = 0 ; i < 64; ++i) {
 551     Int_t board = AliVZEROCalibData::GetBoardNumber(i);
 552     fNBins[i] = TMath::Nint(((Float_t)(fCalibData->GetMatchWindow(board)+1)*25.0+
 553                              (Float_t)kMaxTDCWidth*fCalibData->GetWidthResolution(board))/
 554                             fCalibData->GetTimeResolution(board));
 555     fBinSize[i] = fCalibData->GetTimeResolution(board);
 556   }
 557
 558   fSignalShape = new TF1("VZEROSDigitSignalShape",this,&AliVZERO::SignalShape,0,200,6,"AliVZERO","SignalShape");
 559   fSignalShape->SetParameters(0,1.57345e1,-4.25603e-1,
 560                               2.9,6.40982,3.69339e-01);
 561
 562   fRecoParam = new AliVZERORecoParam;
 563
 564   return;
 565 }
 566
 567 Float_t AliVZERO::CorrectLeadingTime(Int_t i, Float_t time, Float_t adc) const
 568 {
 569   // Correct the leading time
 570   // for slewing effect and
 571   // misalignment of the channels
 572   if (time < 1e-6) return -1024;
 573
 574   // In case of pathological signals
 575   if (adc < 1e-6) return time;
 576
 577   // Slewing correction
 578   Float_t thr = fCalibData->GetCalibDiscriThr(i,kTRUE);
 579   time -= fTimeSlewing->Eval(adc/thr);
 580
 581   return time;
 582 }
 583
 584 double AliVZERO::SignalShape(double *x, double *par)
 585 {
 586   // this function simulates the signal
 587   // shape used in Raw->SDigits method
 588
 589   Double_t xx = x[0];
 590   if (xx <= par[0]) return 0;
 591   Double_t a = 1./TMath::Power((xx-par[0])/par[1],1./par[2]);
 592   if (xx <= par[3]) return a;
 593   Double_t b = 1./TMath::Power((xx-par[3])/par[4],1./par[5]);
 594   Double_t f = a*b/(a+b);
 595   AliDebug(100,Form("x=%f func=%f",xx,f));
 596   return f;
 597 }