]>
Commit | Line | Data |
---|---|---|
80fb7693 | 1 | // $Id$ |
2 | ||
3 | //************************************************************************** | |
4 | //* This file is property of and copyright by the ALICE HLT Project * | |
5 | //* ALICE Experiment at CERN, All rights reserved. * | |
6 | //* * | |
7 | //* Primary Authors: Matthias Richter <Matthias.Richter@ift.uib.no> * | |
8 | //* for The ALICE HLT Project. * | |
9 | //* * | |
10 | //* Permission to use, copy, modify and distribute this software and its * | |
11 | //* documentation strictly for non-commercial purposes is hereby granted * | |
12 | //* without fee, provided that the above copyright notice appears in all * | |
13 | //* copies and that both the copyright notice and this permission notice * | |
14 | //* appear in the supporting documentation. The authors make no claims * | |
15 | //* about the suitability of this software for any purpose. It is * | |
16 | //* provided "as is" without express or implied warranty. * | |
17 | //************************************************************************** | |
18 | ||
19 | /// @file AliHLTDataDeflaterSimple.cxx | |
20 | /// @author Matthias Richter | |
21 | /// @date 2011-08-10 | |
22 | /// @brief Simple deflater implementation storing frequent values below a | |
23 | /// maximum value with a reduced bit number and others with the full | |
24 | /// number of bits. | |
25 | ||
26 | #include "AliHLTDataDeflaterSimple.h" | |
a0d59d54 | 27 | #include "TFile.h" |
28 | #include "TObjArray.h" | |
29 | #include "TH1I.h" | |
30 | #include "TH2F.h" | |
80fb7693 | 31 | #include <memory> |
32 | #include <algorithm> | |
33 | #include <iostream> | |
34 | ||
35 | /** ROOT macro for the implementation of ROOT specific class methods */ | |
36 | ClassImp(AliHLTDataDeflaterSimple) | |
37 | ||
38 | AliHLTDataDeflaterSimple::AliHLTDataDeflaterSimple() | |
39 | : AliHLTDataDeflater() | |
40 | , fParameterDefinitions() | |
a0d59d54 | 41 | , fHistograms(NULL) |
80fb7693 | 42 | { |
43 | // see header file for class documentation | |
44 | // or | |
45 | // refer to README to build package | |
46 | // or | |
47 | // visit http://web.ift.uib.no/~kjeks/doc/alice-hlt | |
a0d59d54 | 48 | if (fHistograms) fHistograms->SetOwner(kTRUE); |
80fb7693 | 49 | } |
50 | ||
51 | AliHLTDataDeflaterSimple::~AliHLTDataDeflaterSimple() | |
52 | { | |
53 | // destructor | |
54 | Clear(); | |
a0d59d54 | 55 | |
56 | if (fHistograms) { | |
57 | delete fHistograms; | |
58 | } | |
59 | fHistograms=NULL; | |
60 | ||
61 | } | |
62 | ||
63 | int AliHLTDataDeflaterSimple::AddParameterDefinition(const char* name, int bitLength, int reducedBitLength) | |
64 | { | |
65 | /// add a parameter definition to the configuration, return reference id | |
66 | fParameterDefinitions.push_back(AliHLTDataDeflaterParameter(name, bitLength, reducedBitLength)); | |
67 | if (fHistograms) { | |
68 | if (!fHistograms->FindObject(name)) { | |
69 | fHistograms->Add(new TH1I(name, name, 100, 0, 0x1<<bitLength)); | |
70 | } | |
71 | } | |
72 | return fParameterDefinitions.size()-1; | |
73 | } | |
74 | ||
75 | int AliHLTDataDeflaterSimple::AddHistogram(TH1* h) | |
76 | { | |
77 | /// add a histogram for deflater statistic of the corresponding parameter | |
78 | if (!fHistograms) fHistograms=new TObjArray; | |
79 | if (!fHistograms) return -ENOMEM; | |
80 | if (h!=NULL && fHistograms->FindObject(h->GetName())) { | |
81 | HLTWarning("parameter with name %s already existing, skipping histogram", h->GetName()); | |
82 | return -EEXIST; | |
83 | } | |
84 | if (h) fHistograms->Add(h); | |
85 | return 0; | |
80fb7693 | 86 | } |
87 | ||
88 | bool AliHLTDataDeflaterSimple::OutputParameterBits( int memberId, AliHLTUInt64_t const & value ) | |
89 | { | |
90 | // write bit pattern of a member to the current byte and position | |
91 | if (memberId>=(int)fParameterDefinitions.size()) return false; | |
92 | ||
93 | AliHLTUInt32_t switchBit=fParameterDefinitions[memberId].SwitchBit(value); // 0 -> reduced, 1 -> full | |
94 | AliHLTUInt64_t v=fParameterDefinitions[memberId].Value(value); | |
95 | AliHLTUInt32_t length=fParameterDefinitions[memberId].ValueLength(value); | |
96 | fParameterDefinitions[memberId].IncrementBitCount(value); | |
97 | ||
a0d59d54 | 98 | if (fHistograms) { |
99 | TObject* o=fHistograms->FindObject(fParameterDefinitions[memberId].GetName()); | |
100 | if (o) { | |
101 | TH1* h=dynamic_cast<TH1*>(o); | |
102 | if (h) { | |
103 | h->Fill(v); | |
104 | } | |
105 | } | |
106 | } | |
107 | ||
80fb7693 | 108 | if (!OutputBit(switchBit)) return false; |
109 | return OutputBits(v, length); | |
110 | } | |
111 | ||
f95bc7cd | 112 | void AliHLTDataDeflaterSimple::Clear(Option_t * option) |
80fb7693 | 113 | { |
114 | // internal cleanup | |
a0d59d54 | 115 | TH2F* hParameterCompression=NULL; |
116 | TH2F* hParameterByteSaving=NULL; | |
117 | if (fHistograms) { | |
118 | int bins=fParameterDefinitions.size(); | |
119 | TObject* o=NULL; | |
120 | o=fHistograms->FindObject("ParameterCompression"); | |
121 | if (o) { | |
122 | hParameterCompression=dynamic_cast<TH2F*>(o); | |
123 | } else { | |
124 | hParameterCompression=new TH2F("ParameterCompression", "ParameterCompression", bins, 0, bins, 100, 0., 1.1); | |
125 | if (hParameterCompression) fHistograms->Add(hParameterCompression); | |
126 | } | |
127 | /* | |
128 | o=fHistograms->FindObject("ParameterByteSaving"); | |
129 | if (o) { | |
130 | hParameterByteSaving=dynamic_cast<TH2F*>(o); | |
131 | } else { | |
132 | hParameterByteSaving=new TH2F("ParameterByteSaving", "ParameterByteSaving", bins, 0, bins, 10, 0., 1.1); | |
133 | if (hParameterByteSaving) fHistograms->Add(hParameterByteSaving); | |
134 | } | |
135 | */ | |
136 | } | |
137 | unsigned i=0; | |
f95bc7cd | 138 | for (vector<AliHLTDataDeflaterParameter>::iterator m=fParameterDefinitions.begin(); |
a0d59d54 | 139 | m!=fParameterDefinitions.end(); m++, i++) { |
140 | int bitLength=m->GetBitLength(); | |
141 | int valueCount=m->GetValueCount(); | |
142 | if (bitLength==0 || valueCount==0) continue; | |
143 | float ratio=(float)m->GetBitCount(); | |
144 | ratio/=bitLength*valueCount; | |
145 | if (hParameterCompression) | |
146 | hParameterCompression->Fill(i, ratio); | |
147 | ratio=(1-ratio)*valueCount*bitLength/8; | |
148 | if (hParameterByteSaving) | |
149 | hParameterByteSaving->Fill(i, ratio); | |
150 | ||
f95bc7cd | 151 | m->ResetBitCount(); |
152 | } | |
153 | AliHLTDataDeflater::Clear(option); | |
80fb7693 | 154 | } |
155 | ||
a0d59d54 | 156 | void AliHLTDataDeflaterSimple::SaveAs(const char *filename,Option_t */*option*/) const |
157 | { | |
158 | // safe histograms to file | |
159 | std::auto_ptr<TFile> file(TFile::Open(filename, "RECREATE")); | |
160 | if (!file.get() || file->IsZombie()) { | |
161 | HLTError("can not open file %s", filename);; | |
162 | return; | |
163 | } | |
164 | file->cd(); | |
165 | if (fHistograms) { | |
166 | fHistograms->Write(); | |
167 | } | |
168 | ||
169 | file->Close(); | |
170 | } | |
171 | ||
80fb7693 | 172 | void AliHLTDataDeflaterSimple::Print(Option_t *option) const |
173 | { | |
174 | // print info | |
175 | Print(cout, option); | |
176 | } | |
177 | ||
178 | void AliHLTDataDeflaterSimple::Print(ostream& out, Option_t *option) const | |
179 | { | |
180 | // print to stream | |
181 | out << "AliHLTDataDeflaterSimple:" << endl; | |
182 | AliHLTUInt64_t bitCount=0; | |
183 | AliHLTUInt64_t fullSize=0; | |
184 | for (vector<AliHLTDataDeflaterParameter>::const_iterator m=fParameterDefinitions.begin(); | |
185 | m!=fParameterDefinitions.end(); m++) { | |
186 | cout << " "; m->Print(option); | |
187 | bitCount+=m->GetBitCount(); | |
188 | fullSize+=m->GetValueCount()*m->GetBitLength(); | |
189 | } | |
190 | out << " total: " << bitCount << "/" << fullSize << " " << (fullSize>0?float(bitCount)/fullSize:0.0) << endl; | |
191 | } | |
192 | ||
193 | void AliHLTDataDeflaterSimple::AliHLTDataDeflaterParameter::Print(const char* /*option*/) const | |
194 | { | |
195 | // print info | |
196 | cout << fName << " (" << fFullBitLength << "," << fReducedBitLength << "): " | |
197 | << fValueCount << " entries " | |
198 | << fBitCount << "/" << fFullBitLength*fValueCount; | |
199 | if (fFullBitLength && fValueCount) { | |
200 | cout << " " << float(fBitCount)/(fValueCount*fFullBitLength); | |
201 | } | |
202 | cout << endl; | |
203 | } | |
204 | ||
205 | ostream& operator<<(ostream &out, const AliHLTDataDeflaterSimple& me) | |
206 | { | |
207 | me.Print(out); | |
208 | return out; | |
209 | } |