updated test classes (alberto)
[u/mrichter/AliRoot.git] / SHUTTLE / test / AliTPCDataDCS.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
17 //
18 // This is an example of a container class
19 // of the data retreieved from the DCS archive DB.
20 // It is called by the detector's Preprocessor and
21 // it is stored in the CDB.
22 //
23
24 #include "AliTPCDataDCS.h"
25
26 #include "AliCDBMetaData.h"
27 #include "AliDCSValue.h"
28 #include "AliLog.h"
29
30 #include <TTimeStamp.h>
31 #include <TObjString.h>
32 #include <TH2F.h>
33 #include <TProfile.h>
34 #include <TGraph.h>
35 #include <TDatime.h>
36 #include <TStyle.h>
37 #include <TCanvas.h>
38
39 ClassImp(AliTPCDataDCS)
40
41 //---------------------------------------------------------------
42 AliTPCDataDCS::AliTPCDataDCS():
43         TObject(),
44         fRun(0),
45         fStartTime(0),
46         fEndTime(0),
47         fGraphs("TGraph",kNGraphs),
48         fIsProcessed(kFALSE)
49 {
50 // default constructor
51
52         for(int i=0;i<kNHistos;i++) fHv[i]=0x0;
53         fFunc = 0;
54 }
55
56 //---------------------------------------------------------------
57 AliTPCDataDCS::AliTPCDataDCS(Int_t nRun, UInt_t startTime, UInt_t endTime):
58         TObject(),
59         fRun(nRun),
60         fStartTime(startTime),
61         fEndTime(endTime),
62         fGraphs("TGraph",kNGraphs),
63         fIsProcessed(kFALSE)
64 {
65 // constructor
66
67         AliInfo(Form("\n\tRun %d \n\tStartTime %s \n\tEndTime %s", nRun,
68         TTimeStamp(startTime).AsString(),
69         TTimeStamp(endTime).AsString()));
70
71         fFunc = 0;
72         Init();
73
74 }
75
76 //---------------------------------------------------------------
77 AliTPCDataDCS::~AliTPCDataDCS()
78 {
79 // destructor
80
81         for(int i=0;i<kNHistos;i++) {delete fHv[i]; fHv[i]=0;}
82         fGraphs.Clear("C");
83         fFunc=0;
84 }
85
86 //---------------------------------------------------------------
87 void AliTPCDataDCS::ProcessData(TMap& aliasMap)
88 {
89 // process TMap of alias values retrieved from the DCS archive DB
90
91         if(!(fHv[0])) Init();
92
93         TObjArray *aliasArr;
94         AliDCSValue* aValue;
95         for(int j=0; j<kNAliases; j++){
96                 aliasArr = (TObjArray*) aliasMap.GetValue(fAliasNames[j].Data());
97                 if(!aliasArr){
98                         AliError(Form("Alias %s not found!", fAliasNames[j].Data()));
99                         continue;
100                 }
101                 Introduce(j, aliasArr);
102
103                 if(aliasArr->GetEntries()<2){
104                         AliError(Form("Alias %s has just %d entries!",
105                                         fAliasNames[j].Data(),aliasArr->GetEntries()));
106                         continue;
107                 }
108
109                 TIter iterarray(aliasArr);
110
111                 Double_t *time = new Double_t[aliasArr->GetEntries()];
112                 Double_t *val = new Double_t[aliasArr->GetEntries()];
113
114                 UInt_t ne=0;
115                 while ((aValue = (AliDCSValue*) iterarray.Next())) {
116                 val[ne] = aValue->GetSimpleValue().GetFloat();
117                 time[ne] = (Double_t) (aValue->GetTimeStamp());
118                 // fill histos (alias 0-2)
119                 if(j < 3) fHv[j]->Fill(val[ne]);
120                 ne++;
121                 }
122                 // fill graphs (alias 3-5)
123                 if(j >= 3) CreateGraph(j, aliasArr->GetEntries(), time, val);
124                 delete[] val;
125                 delete[] time;
126         }
127
128         // calculate mean and rms of the first two histos
129         for(int i=0;i<kNHistos;i++){
130                 fMean[i] = fHv[i]->GetMean();
131                 fWidth[i] = fHv[i]->GetRMS();
132         }
133
134         // pol1 fit of the first graph
135         if(fGraphs.GetEntries() > 0){
136                 ((TGraph*) fGraphs.UncheckedAt(0))->Fit("pol1");
137                 fFunc = ((TGraph*) fGraphs.UncheckedAt(0))->GetFunction("pol1");
138         }
139
140         fIsProcessed=kTRUE;
141
142
143 }
144
145 //---------------------------------------------------------------
146 void AliTPCDataDCS::Init()
147 {
148 // Init alias names and histos
149
150         TH1::AddDirectory(kFALSE);
151
152         fGraphs.SetOwner(1);
153
154         for(int i=0;i<kNAliases;i++){
155                 fAliasNames[i] = "TpcHvSect0";
156                 fAliasNames[i] += i;
157                 fAliasNames[i] += ".FloatValue";
158         }
159
160         for(int i=0;i<kNHistos;i++){
161                 fHv[i] = new TH1F(fAliasNames[i].Data(),fAliasNames[i].Data(), 20, kHvMin, kHvMax);
162                 fHv[i]->GetXaxis()->SetTitle("Hv");
163         }
164 }
165
166 //---------------------------------------------------------------
167 void AliTPCDataDCS::Introduce(UInt_t numAlias, const TObjArray* aliasArr)
168 {
169 // Say something about alias value array being processed
170
171         int entries=aliasArr->GetEntries();
172         AliInfo(Form("************ Alias: %s **********",fAliasNames[numAlias].Data()));
173         AliInfo(Form("          %d DP values collected",entries));
174
175 }
176
177 //---------------------------------------------------------------
178 void AliTPCDataDCS::CreateGraph(int i, int dim, const Double_t *x, const Double_t *y)
179 {
180 // create TGraph to be filled with DP values
181
182         TGraph *gr = new(fGraphs[fGraphs.GetEntriesFast()]) TGraph(dim, x, y);
183
184         gr->GetXaxis()->SetTimeDisplay(1);
185         gr->SetTitle(fAliasNames[i].Data());
186
187         AliInfo(Form("Array entries: %d",fGraphs.GetEntriesFast()));
188
189
190 }
191
192 //---------------------------------------------------------------
193 void AliTPCDataDCS::Draw(const Option_t* /*option*/)
194 {
195 // Draw all histos and graphs
196
197   if(!fIsProcessed) return;
198
199   TCanvas *ch;
200   TString canvasHistoName="Histos";
201   ch=new TCanvas(canvasHistoName,canvasHistoName,20,20,600,600);
202   ch->Divide(2,2);
203   ch->cd(1);
204   fHv[0]->Draw();
205   ch->cd(2);
206   fHv[1]->Draw();
207   ch->cd(3);
208   fHv[2]->Draw();
209
210
211   if(fGraphs.GetEntries() == 0) return;
212
213   TCanvas *cg;
214   TString canvasGraphName="Graphs";
215   cg=new TCanvas(canvasGraphName,canvasGraphName,40,40,600,600);
216   cg->Divide(2,2);
217   cg->cd(1);
218   ((TGraph*) fGraphs.UncheckedAt(0))->Draw("alp");
219
220   cg->cd(2);
221   ((TGraph*) fGraphs.UncheckedAt(1))->Draw("alp");
222   cg->cd(3);
223   ((TGraph*) fGraphs.UncheckedAt(2))->Draw("alp");
224
225   if(fFunc){
226         cg->cd(4);
227         fFunc->Draw("l");
228   }
229
230 }
231