326 std::vector<DetId> boundaryDetIds;
327 std::vector<int> stati;
329 double boundaryEnergy = 0;
330 double boundaryET = 0;
331 int beCellCounter = 0;
332 bool nextToBorder =
false;
336 boundaryEnergy +=
hit->energy();
337 EcalDetId hitdetid = (EcalDetId)
hit->
id();
338 boundaryDetIds.push_back(hitdetid);
341 double eta = cellGeom->getPosition().eta();
342 boundaryET +=
hit->energy() / cosh(
eta);
345 edm::LogInfo(
"EcalBoundaryInfoCalculator") <<
"Find Boundary RecHits...";
349 <<
"Starting at : (" << ((
EBDetId)hitdetid).ieta() <<
"," << ((
EBDetId)hitdetid).iphi() <<
")";
352 <<
"Starting at : (" << ((
EEDetId)hitdetid).ix() <<
"," << ((
EEDetId)hitdetid).iy() <<
")";
359 bool reverseOrientation =
false;
362 EcalDetId
start = hitdetid;
363 EcalDetId current =
start;
364 int current_status = 0;
367 bool startAlgo =
false;
371 theEcalNav->setHome(current);
374 int status = (chit != ecalStatus->
end()) ? chit->getStatusCode() & 0x1F : -1;
380 currDirection =
turnLeft(currDirection, reverseOrientation);
383 throw cms::Exception(
"NoStartingDirection") <<
"EcalBoundaryInfoCalculator: No starting direction can be found: "
384 "This should never happen if RecHit has a dead neighbour!";
389 currDirection =
turnRight(currDirection, reverseOrientation);
392 bool nextIsStart =
false;
393 bool atBorder =
false;
395 while (!nextIsStart) {
396 bool nextStepFound =
false;
399 while (!nextStepFound) {
401 theEcalNav->setHome(current);
404 status = (chit != ecalStatus->
end()) ? chit->getStatusCode() & 0x1F : -1;
407 bool present =
false;
408 for (std::vector<int>::const_iterator
s = stati.begin();
s != stati.end(); ++
s) {
418 nextStepFound =
true;
420 currDirection =
turnRight(currDirection, reverseOrientation);
422 }
else if (
next == EcalDetId(0)) {
424 currDirection =
turnLeft(currDirection, reverseOrientation);
427 nextStepFound =
true;
432 <<
"EcalBoundaryInfoCalculator: No valid next direction can be found: This should never happen!";
442 edm::LogInfo(
"EcalBoundaryInfoCalculator") <<
"Boundary path reached starting position!";
447 <<
"Next step: " << (EcalDetId)
next <<
" Status: " <<
status <<
" Start: " << (EcalDetId)
start;
450 if ((!atBorder ||
status == 0) && !nextIsStart) {
451 boundaryDetIds.push_back(
next);
456 boundaryEnergy += nexthit.
energy();
458 eta = cellGeom->getPosition().eta();
459 boundaryET += nexthit.
energy() / cosh(
eta);
463 if (current_status == 0 &&
status == 0 && atBorder) {
465 currDirection =
turnRight(currDirection, reverseOrientation);
468 if (
status == 0 && atBorder) {
470 currDirection =
turnLeft(currDirection, reverseOrientation);
474 currDirection =
turnLeft(currDirection, reverseOrientation);
477 currDirection =
turnRight(currDirection, reverseOrientation);
487 edm::LogInfo(
"EcalBoundaryInfoCalculator") <<
"<<<<<<<<<<<<<<< Final Boundary object <<<<<<<<<<<<<<<";
489 edm::LogInfo(
"EcalBoundaryInfoCalculator") <<
"no of neighbouring DetIds: " << boundaryDetIds.size();
490 edm::LogInfo(
"EcalBoundaryInfoCalculator") <<
"boundary energy: " << boundaryEnergy;
491 edm::LogInfo(
"EcalBoundaryInfoCalculator") <<
"boundary ET: " << boundaryET;
492 edm::LogInfo(
"EcalBoundaryInfoCalculator") <<
"no of cells contributing to boundary energy: " << beCellCounter;
493 edm::LogInfo(
"EcalBoundaryInfoCalculator") <<
"Channel stati: ";
494 for (std::vector<int>::iterator it = stati.begin(); it != stati.end(); ++it) {
495 edm::LogInfo(
"EcalBoundaryInfoCalculator") << *it <<
" ";
501 boundInfo.
subdet = hitdetid.subdet();
502 boundInfo.
detIds = boundaryDetIds;