84 input_(iConfig.getUntrackedParameter<std::string>(
"inputprocess",
"castorreco")),
85 towercut_(iConfig.getUntrackedParameter<double>(
"towercut",1.)),
86 mintime_(iConfig.getUntrackedParameter<double>(
"mintime",-999)),
87 maxtime_(iConfig.getUntrackedParameter<double>(
"maxtime",999))
90 produces<CastorTowerCollection>();
110 using namespace reco;
111 using namespace TMath;
121 int nRecHits = InputRecHits->size();
124 <<
"2. entering CastorTowerProducer"<<std::endl;
127 LogDebug(
"CastorTowerProducer") <<
"Warning: You are trying to run the Tower algorithm with 0 input rechits.";
132 double poscastortowerarray[4][16];
133 double negcastortowerarray[4][16];
139 for (
int j = 0;
j < 16;
j++) {
140 poscastortowerarray[3][
j] = -2.94524 +
j*0.3927;
141 poscastortowerarray[0][
j] = 0.;
142 poscastortowerarray[1][
j] = 0.;
143 poscastortowerarray[2][
j] = 0.;
145 negcastortowerarray[3][
j] = -2.94524 +
j*0.3927;
146 negcastortowerarray[0][
j] = 0.;
147 negcastortowerarray[1][
j] = 0.;
148 negcastortowerarray[2][
j] = 0.;
152 for (
unsigned int i = 0;
i < InputRecHits->size();
i++) {
156 double Erechit = rechit_p->energy();
159 int sector =
id.sector();
161 if (module < 3) zrechit = -14390 - 24.75 - 49.5*(module-1);
162 if (module > 2) zrechit = -14390 - 99 - 49.5 - 99*(module-3);
163 double phirechit = -100;
164 if (sector < 9) phirechit = 0.19635 + (sector-1)*0.3927;
165 if (sector > 8) phirechit = -2.94524 + (sector - 9)*0.3927;
171 for (
int j=0;
j<16;
j++) {
174 if (TMath::Abs(phirechit - poscastortowerarray[3][
j]) < 0.0001) {
178 poscastortowerarray[0][
j]+=Erechit;
179 if (module < 3) {poscastortowerarray[1][
j]+=Erechit;}
else {poscastortowerarray[2][
j]+=Erechit;}
180 poscastorusedrechits[
j].push_back(rechit_p);
182 negcastortowerarray[0][
j]+=Erechit;
183 if (module < 3) {negcastortowerarray[1][
j]+=Erechit;}
else {negcastortowerarray[2][
j]+=Erechit;}
184 negcastorusedrechits[
j].push_back(rechit_p);
194 double fem, Ehot, depth;
195 double rhoTower = 88.5;
198 for (
int k=0;
k<16;
k++) {
207 fem = poscastortowerarray[1][
k]/poscastortowerarray[0][
k];
212 <<
"tower "<<
k+1<<
": fem = "<<fem<<
" ,depth = "<<depth<<
" ,Ehot = "<<Ehot<<std::endl;
214 TowerPoint temptowerposition(rhoTower,5.9,poscastortowerarray[3][
k]);
215 Point towerposition(temptowerposition);
217 CastorTower newtower(poscastortowerarray[0][k],towerposition,poscastortowerarray[1][k],poscastortowerarray[2][k],fem,depth,Ehot,
218 poscastorusedrechits[k]);
219 OutputTowers->push_back(newtower);
225 fem = negcastortowerarray[1][
k]/negcastortowerarray[0][
k];
230 <<
"tower "<<
k+1 <<
" energy = " << negcastortowerarray[0][
k] <<
"EM = " << negcastortowerarray[1][
k] <<
"HAD = " << negcastortowerarray[2][
k] <<
"phi = " << negcastortowerarray[3][
k] <<
": fem = "<<fem<<
" ,depth = "<<depth<<
" ,Ehot = "<<Ehot<<std::endl;
232 TowerPoint temptowerposition(rhoTower,-5.9,negcastortowerarray[3][
k]);
233 Point towerposition(temptowerposition);
235 CastorTower newtower(negcastortowerarray[0][k],towerposition,negcastortowerarray[1][k],negcastortowerarray[2][k],fem,depth,Ehot,
236 negcastorusedrechits[k]);
237 OutputTowers->push_back(newtower);
242 iEvent.
put(OutputTowers);
249 <<
"Starting CastorTowerProducer";
255 <<
"Ending CastorTowerProducer";
260 using namespace reco;
268 double Erechit = rechit_p->energy();
272 if (module < 3) zrechit = -14390 - 24.75 - 49.5*(module-1);
273 if (module > 2) zrechit = -14390 - 99 - 49.5 - 99*(module-3);
275 if(Erechit > Ehot) Ehot = Erechit;
276 depth+=Erechit*zrechit;
edm::SortedCollection< CastorRecHit > CastorRecHitCollection
ROOT::Math::RhoZPhiPoint CellPoint
#define DEFINE_FWK_MODULE(type)
std::vector< reco::CastorTower > CastorTowerCollection
ROOT::Math::RhoEtaPhiPoint TowerPoint
OrphanHandle< PROD > put(std::auto_ptr< PROD > product)
Put a new product.
virtual void produce(edm::Event &, const edm::EventSetup &)
edm::RefVector< CastorRecHitCollection > CastorRecHitRefVector
bool getByLabel(InputTag const &tag, Handle< PROD > &result) const
ROOT::Math::PositionVector3D< ROOT::Math::Cartesian3D< double > > XYZPointD
point in space with cartesian internal representation
CastorTowerProducer(const edm::ParameterSet &)
virtual void ComputeTowerVariable(const edm::RefVector< edm::SortedCollection< CastorRecHit > > &usedRecHits, double &Ehot, double &depth)
ProductID id() const
Accessor for product ID.