[u/mrichter/AliRoot.git] / PWGCF / Correlations / macros / dptdptcorrelations / AddTaskCorrMC.C

//
// Macro designed for use with the AliDptDptInMC task.
//
// Author: Prabhat Pujahari & Claude Pruneau, Wayne State
// 
//           system:  0: PbPb                 1: pPb
//      singlesOnly:  0: full correlations    1: singles only
//       useWeights:  0: no                   1: yes
// centralityMethod:  3: track count  4: V0 centrality  7: V0A centrality for pPb
//        chargeSet:  0: ++    1: +-    2: -+    3: --
/////////////////////////////////////////////////////////////////////////////////

AliDptDptInMC *AddTaskCorrMC
(int    singlesOnly            = 1,
 int    useWeights             = 0,
 int    centralityMethod       = 4,
 int    chargeSet              = 1,
 int    trackFilterBit         = 128,
 double etaMin                 = -0.8,
 double etaMax                 =  0.8,
 double dcaZMin                = -3.0,
 double dcaZMax                =  3.0,
 double dcaXYMin               = -2.4,
 double dcaXYMax               =  2.4,
 int nCentrality               =  5, 
 TString anadata               = "MCAODreco",
 Bool_t NoResonances           =  kTRUE,
 Bool_t NoElectron             =  kTRUE,
 const char* taskname          = "MCTruth",
 char *inputHistogramFileName  = "alien:///alice/cern.ch/user/p/prabhat/MCCalib/MCTruthCaalib.root")
 
{
 
  // Set Default Configuration of this analysis
  // ==========================================
  int    debugLevel             = 0;
  int    rejectPileup           = 1;
  int    rejectPairConversion   = 1;
  int    sameFilter             = 1;

  
  //int    nCentrality;
  double minCentrality[10];
  double maxCentrality[10];
  
      
      if (centralityMethod == 4)
	{
	  
	  minCentrality[0] = 0.0; maxCentrality[0] = 10.0; 
	  minCentrality[1] = 10.0; maxCentrality[1] = 20.0; 
	  minCentrality[2] = 20.0; maxCentrality[2] = 30.0; 
	  minCentrality[3] = 30.0; maxCentrality[3] = 40.0;
	  minCentrality[4] = 60.0; maxCentrality[4] = 70.0;
	  
	}
      else
	{
	  cout << "-F- AddTaskDptDptCorrelations() system:" << system << ". centralityMethod:" << centralityMethod << " Option NOT AVAILABLE. ABORT."
	    return 0;
	}
    
  
  double zMin                   = -10.;
  double zMax                   =  10.;
  double ptMin                  =  0.2;
  double ptMax                  =  2.0;
  double dedxMin                =  0.0;
  double dedxMax                =  20000.0;
  int    requestedCharge1       =  1; //default
  int    requestedCharge2       = -1; //default
  
  
  // Get the pointer to the existing analysis manager via the static access method.
  // ==============================================================================
  AliAnalysisManager *analysisManager = AliAnalysisManager::GetAnalysisManager();
  
  if (!analysisManager) 
    {
    ::Error("AddTaskDptDptCorrelations", "No analysis manager to connect to.");
    return NULL;
    }  
  
  TString part1Name;
  TString part2Name;
  TString eventName;
  TString prefixName        = "Corr_";
  TString pileupRejecSuffix = "_PileupRejec";
  TString pairRejecSuffix   = "_PairRejec";
  TString calibSuffix       = "_calib";
  TString singlesOnlySuffix = "_SO";
  TString suffix;
  
  TString inputPath         = ".";
  TString outputPath        = ".";
  TString baseName;
  TString listName;
  TString taskName;
  //TString inputHistogramFileName;
  TString outputHistogramFileName;
  
  // Create the task and add subtask.
  // ===========================================================================
  int iTask = 0; // task counter
  AliAnalysisDataContainer *taskInputContainer;
  AliAnalysisDataContainer *taskOutputContainer;
  AliDptDptInMC* task;
  
  for (int iCentrality=0; iCentrality < nCentrality; ++iCentrality)
    {
      switch (chargeSet)
        {
          case 0: part1Name = "P_"; part2Name = "P_"; requestedCharge1 =  1; requestedCharge2 =  1; sameFilter = 1; break;
          case 1: part1Name = "P_"; part2Name = "M_"; requestedCharge1 =  1; requestedCharge2 = -1; sameFilter = 0;   break;
          case 2: part1Name = "M_"; part2Name = "P_"; requestedCharge1 = -1; requestedCharge2 =  1; sameFilter = 0;   break;
          case 3: part1Name = "M_"; part2Name = "M_"; requestedCharge1 = -1; requestedCharge2 = -1; sameFilter = 1;   break;
        }
      //part1Name += int(1000*etaMin);                                                                                        
      part1Name += "eta";
      part1Name += int(1000*etaMax);
      part1Name += "_";
      part1Name += int(1000*ptMin);
      part1Name += "pt";
      part1Name += int(1000*ptMax);
      part1Name += "_";
      part1Name += int(1000*dcaZMin);
      part1Name += "DCA";
      part1Name += int(1000*dcaZMax);
      part1Name += "_";


      //part2Name += int(1000*etaMin);                                                                                        
      part2Name += "eta";
      part2Name += int(1000*etaMax);
      part2Name += "_";
      part2Name += int(1000*ptMin);
      part2Name += "pt";
      part2Name += int(1000*ptMax);
      part2Name += "_";
      part2Name += int(1000*dcaZMin);
      part2Name += "DCA";
      part2Name += int(1000*dcaZMax);
      part2Name += "_";

      eventName =  "";
      eventName += int(10.*minCentrality[iCentrality] );
      eventName += "Vo";
      eventName += int(10.*maxCentrality[iCentrality] );
      

      baseName     =   prefixName;
      baseName     +=  part1Name;
      baseName     +=  part2Name;
      baseName     +=  eventName;
      listName     =   baseName;
      taskName     =   baseName;
      
      
      outputHistogramFileName = baseName;
      if (singlesOnly) outputHistogramFileName += singlesOnlySuffix;
      outputHistogramFileName += ".root";
      
    
      TFile  * inputFile  = 0;
      TList  * histoList  = 0;
      TH3F   * weight_1   = 0;
      TH3F   * weight_2   = 0;
      if (useWeights)
        {
        TGrid::Connect("alien:");
        inputFile = TFile::Open(inputHistogramFileName,"OLD");
        if (!inputFile)
          {
          cout << "Requested file:" << inputHistogramFileName << " was not opened. ABORT." << endl;
          return;
          }
        TString nameHistoBase = "correction_";
        TString nameHisto;
        nameHistoBase += eventName;
        if (requestedCharge1 == 1)
          {
          nameHisto = nameHistoBase + "_p";
          cout << "Input Histogram named: " << nameHisto << endl;
          weight_1 = (TH3F *) inputFile->Get(nameHisto);
          }
        else
          {
          nameHisto = nameHistoBase + "_m";
          cout << "Input Histogram named: " << nameHisto << endl;
          weight_1 = (TH3F *) inputFile->Get(nameHisto);
          }
        if (!weight_1) 
          {
          cout << "Requested histogram 'correction_p/m' was not found. ABORT." << endl;
          return 0;
          }
        
        if (!sameFilter)
          {
          weight_2 = 0;
          if (requestedCharge2 == 1)
            {
            nameHisto = nameHistoBase + "_p";
            cout << "Input Histogram named: " << nameHisto << endl;
            weight_2 = (TH3F *) inputFile->Get(nameHisto);
            }
          else
            {
            nameHisto = nameHistoBase + "_m";
            cout << "Input Histogram named: " << nameHisto << endl;
            weight_2 = (TH3F *) inputFile->Get(nameHisto);
            }
          if (!weight_2) 
            {
            cout << "Requested histogram 'correction_p/m' was not found. ABORT." << endl;
            return 0;
            }
          }  
        }

      task = new  AliDptDptInMC(taskName);
      //configure my task
      task->SetDebugLevel(          debugLevel      ); 
      task->SetSameFilter(          sameFilter      );
      task->SetSinglesOnly(         singlesOnly     ); 
      task->SetUseWeights(          useWeights      ); 
      task->SetRejectPileup(        rejectPileup    ); 
      task->SetRejectPairConversion(rejectPairConversion); 
      task->SetVertexZMin(          zMin            ); 
      task->SetVertexZMax(          zMax            ); 
      task->SetVertexXYMin(         -1.            ); 
      task->SetVertexXYMax(          1.            ); 
      task->SetCentralityMethod(    centralityMethod);
      task->SetCentrality(          minCentrality[iCentrality], maxCentrality[iCentrality]);
      task->SetPtMin1(              ptMin           ); 
      task->SetPtMax1(              ptMax           ); 
      task->SetEtaMin1(             etaMin          ); 
      task->SetEtaMax1(             etaMax          ); 
      task->SetPtMin2(              ptMin           ); 
      task->SetPtMax2(              ptMax           ); 
      task->SetEtaMin2(             etaMin          ); 
      task->SetEtaMax2(             etaMax          ); 
      task->SetDcaZMin(             dcaZMin         ); 
      task->SetDcaZMax(             dcaZMax         ); 
      task->SetDcaXYMin(            dcaXYMin        ); 
      task->SetDcaXYMax(            dcaXYMax        ); 
      task->SetDedxMin(             dedxMin         ); 
      task->SetDedxMax(             dedxMax         ); 
      //task->SetNClusterMin(         nClusterMin     ); 
      task->SetTrackFilterBit(      trackFilterBit  );
      task->SetRequestedCharge_1(   requestedCharge1); 
      task->SetRequestedCharge_2(   requestedCharge2); 
      task->SetWeigth_1(            weight_1        );
      task->SetWeigth_2(            weight_2        );
      task->SetAnalysisType(anadata);
      
      task->SelectCollisionCandidates(AliVEvent::kMB);

      task->SetResonancesCut(NoResonances);
      task->SetElectronCut(NoElectron);
            
      cout << "Creating task output container" << endl;
      
      taskOutputContainer = analysisManager->CreateContainer(listName, 
                                                             TList::Class(),    
                                                             AliAnalysisManager::kOutputContainer, 
                                                             Form("%s:%s", AliAnalysisManager::GetCommonFileName(),taskname));
      cout << "Add task to analysis manager and connect it to input and output containers" << endl;
           

      analysisManager->AddTask(task);
      analysisManager->ConnectInput( task,  0, analysisManager->GetCommonInputContainer());
      analysisManager->ConnectOutput(task,  0, taskOutputContainer );
      cout << "Task added ...." << endl;
      
      iTask++;
    
    }
  
  return task;
}
Commit	Line	Data
1486c432 MW	1	//
	2	// Macro designed for use with the AliDptDptInMC task.
	3	//
	4	// Author: Prabhat Pujahari & Claude Pruneau, Wayne State
	5	//
	6	// system: 0: PbPb 1: pPb
	7	// singlesOnly: 0: full correlations 1: singles only
	8	// useWeights: 0: no 1: yes
	9	// centralityMethod: 3: track count 4: V0 centrality 7: V0A centrality for pPb
	10	// chargeSet: 0: ++ 1: +- 2: -+ 3: --
	11	/////////////////////////////////////////////////////////////////////////////////
	12
	13	AliDptDptInMC *AddTaskCorrMC
947cb887	14	(int singlesOnly = 1,
947cb887	15	int useWeights = 0,
1486c432 MW	16	int centralityMethod = 4,
	17	int chargeSet = 1,
	18	int trackFilterBit = 128,
	19	double etaMin = -0.8,
	20	double etaMax = 0.8,
	21	double dcaZMin = -3.0,
	22	double dcaZMax = 3.0,
	23	double dcaXYMin = -2.4,
	24	double dcaXYMax = 2.4,
947cb887	25	int nCentrality = 5,
947cb887	26	TString anadata = "MCAODreco",
1486c432 MW	27	Bool_t NoResonances = kTRUE,
1486c432 MW	28	Bool_t NoElectron = kTRUE,
9aba9023	29	const char* taskname = "MCTruth",
1486c432 MW	30	char *inputHistogramFileName = "alien:///alice/cern.ch/user/p/prabhat/MCCalib/MCTruthCaalib.root")
	31
	32	{
	33
	34	// Set Default Configuration of this analysis
	35	// ==========================================
	36	int debugLevel = 0;
	37	int rejectPileup = 1;
	38	int rejectPairConversion = 1;
	39	int sameFilter = 1;
	40
	41
	42	//int nCentrality;
	43	double minCentrality[10];
	44	double maxCentrality[10];
	45
	46
	47	if (centralityMethod == 4)
	48	{
	49
9aba9023	50	minCentrality[0] = 0.0; maxCentrality[0] = 10.0;
947cb887	51	minCentrality[1] = 10.0; maxCentrality[1] = 20.0;
	52	minCentrality[2] = 20.0; maxCentrality[2] = 30.0;
	53	minCentrality[3] = 30.0; maxCentrality[3] = 40.0;
	54	minCentrality[4] = 60.0; maxCentrality[4] = 70.0;
1486c432 MW	55
	56	}
	57	else
	58	{
	59	cout << "-F- AddTaskDptDptCorrelations() system:" << system << ". centralityMethod:" << centralityMethod << " Option NOT AVAILABLE. ABORT."
	60	return 0;
	61	}
	62
	63
	64	double zMin = -10.;
	65	double zMax = 10.;
	66	double ptMin = 0.2;
	67	double ptMax = 2.0;
	68	double dedxMin = 0.0;
	69	double dedxMax = 20000.0;
	70	int requestedCharge1 = 1; //default
	71	int requestedCharge2 = -1; //default
	72
	73
	74	// Get the pointer to the existing analysis manager via the static access method.
	75	// ==============================================================================
	76	AliAnalysisManager *analysisManager = AliAnalysisManager::GetAnalysisManager();
	77
	78	if (!analysisManager)
	79	{
	80	::Error("AddTaskDptDptCorrelations", "No analysis manager to connect to.");
	81	return NULL;
	82	}
	83
	84	TString part1Name;
	85	TString part2Name;
	86	TString eventName;
	87	TString prefixName = "Corr_";
	88	TString pileupRejecSuffix = "_PileupRejec";
	89	TString pairRejecSuffix = "_PairRejec";
	90	TString calibSuffix = "_calib";
	91	TString singlesOnlySuffix = "_SO";
	92	TString suffix;
	93
	94	TString inputPath = ".";
	95	TString outputPath = ".";
	96	TString baseName;
	97	TString listName;
	98	TString taskName;
	99	//TString inputHistogramFileName;
	100	TString outputHistogramFileName;
	101
	102	// Create the task and add subtask.
	103	// ===========================================================================
	104	int iTask = 0; // task counter
	105	AliAnalysisDataContainer *taskInputContainer;
	106	AliAnalysisDataContainer *taskOutputContainer;
	107	AliDptDptInMC* task;
	108
	109	for (int iCentrality=0; iCentrality < nCentrality; ++iCentrality)
	110	{
	111	switch (chargeSet)
	112	{
	113	case 0: part1Name = "P_"; part2Name = "P_"; requestedCharge1 = 1; requestedCharge2 = 1; sameFilter = 1; break;
	114	case 1: part1Name = "P_"; part2Name = "M_"; requestedCharge1 = 1; requestedCharge2 = -1; sameFilter = 0; break;
	115	case 2: part1Name = "M_"; part2Name = "P_"; requestedCharge1 = -1; requestedCharge2 = 1; sameFilter = 0; break;
	116	case 3: part1Name = "M_"; part2Name = "M_"; requestedCharge1 = -1; requestedCharge2 = -1; sameFilter = 1; break;
	117	}
	118	//part1Name += int(1000*etaMin);
119	part1Name += "eta";
120	part1Name += int(1000*etaMax);
121	part1Name += "_";
122	part1Name += int(1000*ptMin);
123	part1Name += "pt";
124	part1Name += int(1000*ptMax);
125	part1Name += "_";
126	part1Name += int(1000*dcaZMin);
127	part1Name += "DCA";
128	part1Name += int(1000*dcaZMax);
129	part1Name += "_";
130
131
132	//part2Name += int(1000*etaMin);
133	part2Name += "eta";
134	part2Name += int(1000*etaMax);
135	part2Name += "_";
136	part2Name += int(1000*ptMin);
137	part2Name += "pt";
138	part2Name += int(1000*ptMax);
139	part2Name += "_";
140	part2Name += int(1000*dcaZMin);
141	part2Name += "DCA";
142	part2Name += int(1000*dcaZMax);
143	part2Name += "_";
144
145	eventName = "";
146	eventName += int(10.*minCentrality[iCentrality] );
147	eventName += "Vo";
148	eventName += int(10.*maxCentrality[iCentrality] );
149
150
151	baseName = prefixName;
152	baseName += part1Name;
153	baseName += part2Name;
154	baseName += eventName;
155	listName = baseName;
156	taskName = baseName;
157
158
159	outputHistogramFileName = baseName;
160	if (singlesOnly) outputHistogramFileName += singlesOnlySuffix;
161	outputHistogramFileName += ".root";
162
163
164	TFile * inputFile = 0;
165	TList * histoList = 0;
166	TH3F * weight_1 = 0;
167	TH3F * weight_2 = 0;
168	if (useWeights)
169	{
170	TGrid::Connect("alien:");
171	inputFile = TFile::Open(inputHistogramFileName,"OLD");
172	if (!inputFile)
173	{
174	cout << "Requested file:" << inputHistogramFileName << " was not opened. ABORT." << endl;
175	return;
176	}
177	TString nameHistoBase = "correction_";
178	TString nameHisto;
179	nameHistoBase += eventName;
180	if (requestedCharge1 == 1)
181	{
182	nameHisto = nameHistoBase + "_p";
183	cout << "Input Histogram named: " << nameHisto << endl;
184	weight_1 = (TH3F *) inputFile->Get(nameHisto);
185	}
186	else
187	{
188	nameHisto = nameHistoBase + "_m";
189	cout << "Input Histogram named: " << nameHisto << endl;
190	weight_1 = (TH3F *) inputFile->Get(nameHisto);
191	}
192	if (!weight_1)
193	{
194	cout << "Requested histogram 'correction_p/m' was not found. ABORT." << endl;
195	return 0;
196	}
197
198	if (!sameFilter)
199	{
200	weight_2 = 0;
201	if (requestedCharge2 == 1)
202	{
203	nameHisto = nameHistoBase + "_p";
204	cout << "Input Histogram named: " << nameHisto << endl;
205	weight_2 = (TH3F *) inputFile->Get(nameHisto);
206	}
207	else
208	{
209	nameHisto = nameHistoBase + "_m";
210	cout << "Input Histogram named: " << nameHisto << endl;
211	weight_2 = (TH3F *) inputFile->Get(nameHisto);
212	}
213	if (!weight_2)
214	{
215	cout << "Requested histogram 'correction_p/m' was not found. ABORT." << endl;
216	return 0;
217	}
218	}
219	}
220
221	task = new AliDptDptInMC(taskName);
222	//configure my task
223	task->SetDebugLevel( debugLevel );
224	task->SetSameFilter( sameFilter );
225	task->SetSinglesOnly( singlesOnly );
226	task->SetUseWeights( useWeights );
227	task->SetRejectPileup( rejectPileup );
228	task->SetRejectPairConversion(rejectPairConversion);
229	task->SetVertexZMin( zMin );
230	task->SetVertexZMax( zMax );
231	task->SetVertexXYMin( -1. );
232	task->SetVertexXYMax( 1. );
233	task->SetCentralityMethod( centralityMethod);
234	task->SetCentrality( minCentrality[iCentrality], maxCentrality[iCentrality]);
235	task->SetPtMin1( ptMin );
236	task->SetPtMax1( ptMax );
237	task->SetEtaMin1( etaMin );
238	task->SetEtaMax1( etaMax );
239	task->SetPtMin2( ptMin );
240	task->SetPtMax2( ptMax );
241	task->SetEtaMin2( etaMin );
242	task->SetEtaMax2( etaMax );
243	task->SetDcaZMin( dcaZMin );
244	task->SetDcaZMax( dcaZMax );
245	task->SetDcaXYMin( dcaXYMin );
246	task->SetDcaXYMax( dcaXYMax );
247	task->SetDedxMin( dedxMin );
248	task->SetDedxMax( dedxMax );
249	//task->SetNClusterMin( nClusterMin );
250	task->SetTrackFilterBit( trackFilterBit );
251	task->SetRequestedCharge_1( requestedCharge1);
252	task->SetRequestedCharge_2( requestedCharge2);
253	task->SetWeigth_1( weight_1 );
254	task->SetWeigth_2( weight_2 );
255	task->SetAnalysisType(anadata);
256
257	task->SelectCollisionCandidates(AliVEvent::kMB);
258
259	task->SetResonancesCut(NoResonances);
260	task->SetElectronCut(NoElectron);
261
262	cout << "Creating task output container" << endl;
263
264	taskOutputContainer = analysisManager->CreateContainer(listName,
265	TList::Class(),
266	AliAnalysisManager::kOutputContainer,
267	Form("%s:%s", AliAnalysisManager::GetCommonFileName(),taskname));
268	cout << "Add task to analysis manager and connect it to input and output containers" << endl;
269
270
271	analysisManager->AddTask(task);
272	analysisManager->ConnectInput( task, 0, analysisManager->GetCommonInputContainer());
273	analysisManager->ConnectOutput(task, 0, taskOutputContainer );
274	cout << "Task added ...." << endl;
275
276	iTask++;
277
278	}
279
280	return task;
281	}