]>
Commit | Line | Data |
---|---|---|
4b4eb769 | 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 |