////////////////////////////////////////////////////////////////////////////////////////////////////////
//
// Component: CalDataMaker.cc
// Purpose: Implementation of CalDataMaker class. 
//
// Created: 15/04/99 Pierre Savard (adapted from J. Lamoureux's 
//                                  CalBankMaker)
// History: 
//          10/25/01 Beate Heinmann: abort Job if no DB access
//          02/18/01 Beate Heinmann: bug fix: if the HTDC time was non-zero for two 
//                                   hadronic calorimeters (e.g. CHA and WHA) the information
//                                   was only stored correctly for one of them but 
//                                   not the other. This is corrected now
//          02/19/01 Beate Heinmann: slewing corrections are applied as function of the 
//                                   geometric mean rather than arithmetic mean
//          06/14/01 Beate Heinmann: add MiniPlug data
//          27/10/02 Angela Wyatt: set validCode for bad TDC efficiencies
// 	    06/30/03 Willis Sakumoto: standarize/validate LER corrections from DB
//          17/09/04 Angela Wyatt: Add new MiniPlug inst. lumi. dependent calibrations
//          09/01/06 Emily Nurse: Get rid of severe error message when inst. lumi is not found for MC in MiniPlug energy correction
////////////////////////////////////////////////////////////////////////////////////////////////////////
#include <iostream>
#include <iomanip>
#include <cstdlib>
#include <ctime>
#include <math.h>
#include "inc/bcs.h"
#include "BaBar/Cdf.hh"
#include "ErrorLogger_i/gERRLOG.hh"

// Calorimeter D bank headers
#include "RawDataBanks/RawQieStorableBank.hh"
#include "RawDataBanks/CEMD_StorableBank.hh"
#include "RawDataBanks/PEMD_StorableBank.hh"
#include "RawDataBanks/CHAD_StorableBank.hh"
#include "RawDataBanks/WHAD_StorableBank.hh"
#include "RawDataBanks/PHAD_StorableBank.hh"
#include "RawDataBanks/PPRD_StorableBank.hh"
#include "RawDataBanks/HATD_StorableBank.hh"
#include "RawDataBanks/EMTD_StorableBank.hh"
#include "RawDataBanks/MPAD_StorableBank.hh"
#include "RawDataBanks/CLLD_StorableBank.hh"
#include "HeaderObjects/LRID_StorableBank_utilities.hh"
#include "CalorGeometry/CalorKey.hh"
#include "Edm/Handle.hh"
#include "AbsEnv/AbsEnv.hh"
// Geometry headers
#include "GeometryBase/CdfDetector.hh"
#include "CalorGeometry/AbsCalorDetectorNode.hh"
#include "CalorGeometry/HWCalorDetectorNode.hh"

// Calorimetry class headers
#include "CalorObjects/CalTower.hh"
#include "CalorObjects/TowerType.hh"
#include "Calor/CalDataMaker.hh"
#include "CalorObjects/TdcKey.hh"
#include "CalorObjects/CalData.hh"
using namespace std;

// namespace calor {


  //
  // Memory management
  //

CalDataMaker::CalDataMaker(){
  
  _calData = 0;  
  _calib = Calib();
  _dbInitialized = false;
  initJob(); 
}

//-------------------------------------------------------------------
//
//  Destructor
//
//-------------------------------------------------------------------
CalDataMaker::~CalDataMaker() {
}

void CalDataMaker::initJob() {
  for (int ie=0;ie<TENETA;ie++) {
    for (int ip=0;ip<TENPHI;ip++) {
      cemCal[ie][ip]=1.0;
      pemCal[ie][ip]=1.0;
      pprCal[ie][ip]=1.0;
      chaCal[ie][ip]=1.0;
      phaCal[ie][ip]=1.0;
      whaCal[ie][ip]=1.0;
    }
  }
}

void CalDataMaker::initRun(bool verbose) {
  initRun(verbose,true);
}
void CalDataMaker::initRun(bool verbose, bool accDB) {
  initRun(verbose,true,false);
}

void CalDataMaker::initRun(bool verbose, bool accDB, bool offScl) {

// Pick up the calorimetry calibrations from the database. It is not
// combined with other calorimetry classes that use these calibrations
// so that CalDataMaker can be used by itself.         WKS:24Jun03 

  _dbInitialized = false;

  if (accDB) {
    bool result = _calib.initializeDB();
    if(result) _dbInitialized = true;
    if(!result &&  !gblEnv->monteFlag()   )  ERRLOG(ELabort,"CalDataMaker: ") <<
		   "Error initialising DB!" << endmsg;
    if(_dbInitialized &&  !gblEnv->monteFlag()){
	bool result = _calib.loadCalibConsts();
	if(!result)  ERRLOG(ELabort,"CalDataMaker: ") <<
	             "Error fetching calibration constants from DB!" << endmsg;
    }
  } else {
    // load fake constants if only running for test purposes
    _calib.loadFakeConsts();
    ERRLOG(ELwarning,"CalDataMaker: ") <<
	   "Using fake calibration constants !" << endmsg;
  }

  if(gblEnv->monteFlag() ||
     AbsEnv::instance()->runType()==LRID_run_type::monte_carlo_embedded)
  {
    _calib.loadMCConsts();
    // set tower-by-tower calibrations to unity for simulation
    for (int ie=0;ie<TENETA;ie++){
      for (int ip=0;ip<TENPHI;ip++){
	cemCal[ie][ip]=1.0;
	pemCal[ie][ip]=1.0;
	pprCal[ie][ip]=1.0;
	chaCal[ie][ip]=1.0;
	whaCal[ie][ip]=1.0;
	phaCal[ie][ip]=1.0;
      }
    }
  } else {

    // For data requiring pmt-by-pmt LER calibrations from the
    // database. Otherwise, initJob() default: LER=1
 
    if ( _dbInitialized ) {
      for ( int ie = 0; ie < TENETA; ++ie ) {
      for ( int ip = 0; ip < TENPHI; ++ip ) {
        cemCal[ie][ip] = _calib.cemLer( ie, ip );
        pemCal[ie][ip] = _calib.pemLer( ie, ip );
        pprCal[ie][ip] = _calib.pprLer( ie, ip );
        chaCal[ie][ip] = _calib.chaLer( ie, ip );
        whaCal[ie][ip] = _calib.whaLer( ie, ip );
        phaCal[ie][ip] = _calib.phaLer( ie, ip );
      }
      }
    }
  }
  
  if(verbose) {
    _calib.print(std::cout);
    //_calib.printEffic(std::cout);}
    //_calib.printBadTowers(std::cout);
    //_calib.printLers(std::cout);
  }
}



void CalDataMaker::initEvent(AbsEvent* anEvent) {
  
  
  // make a CalData object
  // CalData(anEvent);
  _calData    =  new CalData();
}

//_____________________________________________________________________________
void CalDataMaker::initEvent(CalData* caldata) {
  // make a copy of the old CalData object for the simulation 

  CalDataCalib* clb;

  _calData       = new CalData(*caldata);
//-----------------------------------------------------------------------------
// set scales/pedestals to be used by the maker
//-----------------------------------------------------------------------------
  clb   = caldata->getCalibrations();
  _calib.loadConsts(clb);

  _dbInitialized = true;
}

void CalDataMaker::endEvent(AbsEvent* anEvent, bool verbose){

  this->endEvent(anEvent, 0, 0, verbose);
}

void CalDataMaker::endEvent(AbsEvent* anEvent, int skcut, float etmin, bool verbose){
  // copy calibration constants inside CalData
  _calData->storeCalibrations(_calib);

  // filter CalData - and make emEnergy and hadEnergy.
  _calData->cleanup(skcut,etmin);
  if (verbose) { 
    _calData->printSummary();    
    _calData->print(std::cout);
  }
}

int CalDataMaker::writeCalData(AbsEvent* anEvent, bool verbose){

  // EDM stuff goes here:
  // Caldata is now made - hook the thing to the event
  GenericConstHandle writeHandle(anEvent->append(_calData));
  _calData = 0;
  return 1;
}

//-------------------------------------------------------------------
//   Make:
//   read in calorimeter D banks and make the CalData object
//
//-------------------------------------------------------------------

int CalDataMaker::makeCalData(AbsEvent* anEvent,bool verbose, bool thresh) {
 
  int result = makeCalData(anEvent,verbose,thresh,false);
  return result;
}

int CalDataMaker::makeCalData(AbsEvent* anEvent,bool verbose, bool thresh, bool plugGain) {

  int result = makeCalData(anEvent,verbose,thresh,plugGain,false);
  return result;
}

int CalDataMaker::makeCalData(AbsEvent* anEvent,bool verbose, bool thresh, bool plugGain, bool useMPALumiCalib) {
  
  CalorKey calKey;
  CellKey    cKey;

  //
  // CEM
  //

  EventRecord::ConstIterator p_cemdBank( anEvent,  "CEMD_StorableBank" );
  if(p_cemdBank.is_valid()){ 

    ConstHandle<RawQieStorableBank> bank_handle(p_cemdBank) ;
    
    //  Loop over the bank we just found. 
    
    for(ConstGrandBankIterI2s1<RawQieStorableBank> data_iter(bank_handle); 
	data_iter.is_valid(); ++data_iter) {
	      
      //  convert address 
      
      calKey  = bank_handle->get_calKey( data_iter );
      if (calKey.validKey(calKey)) {
	cKey = CellKey(calKey);
	
	if ( !cKey.isValid() ) {
	  ERRLOG( ELsevere, "CalDataMaker" ) << " unpack CEMD bank : CellKey is inValid => Data probably corrupted" << endmsg;
	}
	else {
	  
	  //  Apply SCL and LER calibrations to convert from digitisation units to 
	  //  energy.
	  
	  float energy = (bank_handle->get_energy( data_iter ) -_calib.pedCem())* _calib.detCem();
	  
	  uint2 iData = bank_handle->get_data( data_iter ); 
	  int   ipmt  = bank_handle->get_pmt(  data_iter );
          int   wedge = bank_handle->get_phi(  data_iter );

	  if (TEPSEG[cKey.towerType()]==48) ipmt = 0;
	  // Apply threshold and remove chimney towers (is this okay?)
	  if(thresh == false || (energy > CEMMIN && !(cKey.iPhi()==5 && cKey.iEta()>33 &&
						      cKey.iEta()<36 ))){
	    
	    // Put the energy into Caldata:
	    CalTower* t = _calData->towerMake(cKey);
	    // apply offline LERs here only to energy not to iData, i.e. leave raw data untouched
	    energy *= cemCal[cKey.iEta()][2*wedge + ipmt];			// WKS:21Jun03
	    t->addPmtEnergy(calKey.detector(),ipmt,energy,iData);
	  //		  calKey.print();
	  }
	}
      } else {
	ERRLOG( ELsevere, "CalDataMaker" ) << " unpack CEMD bank : invalid tower Key  => Data probably corrupted" << endmsg;
      }
	      
    }                         // end loop over data
  }
  
  //
  // PEM
  //

  EventRecord::ConstIterator p_pemdBank( anEvent,  "PEMD_StorableBank" );	
  if(p_pemdBank.is_valid()){ 
    ConstHandle<RawQieStorableBank> bank_handle(p_pemdBank) ;
	    
    //  Loop over the bank we just found. 
	    
    for(ConstGrandBankIterI2s1<RawQieStorableBank> data_iter(bank_handle);
	data_iter.is_valid(); ++data_iter) {
	      
      //  convert address 
	      
      calKey  = bank_handle->get_calKey( data_iter );
      if (calKey.validKey(calKey)) {
	cKey = CellKey(calKey);
	if ( !cKey.isValid() ) {
	  ERRLOG( ELsevere, "CalDataMaker" ) << " unpack PEMD bank : CellKey is inValid => Data probably corrupted" << endmsg;
	}
	else {
	  float energy = (bank_handle->get_energy( data_iter )-_calib.pedPem())* _calib.detPem();

	  uint2 iData = bank_handle->get_data( data_iter );  
	  int   ipmt  = bank_handle->get_pmt(  data_iter );
          int   wedge = bank_handle->get_phi(  data_iter );
	  
	  // Put the energy in CalData. The LERs can be traced by setting energy = LER
	  
	  if(thresh == false || (energy >PEMMIN)){
	    CalTower* t = _calData->towerMake(cKey);
	    // apply offline LERs here only to energy not to iData, i.e. leave raw data untouched
	    energy *= pemCal[cKey.iEta()][2*wedge + ipmt];			// WKS:21Jun03
	    t->addPmtEnergy(calKey.detector(),ipmt,energy,iData);
	  }
	}
      } else {
	ERRLOG( ELsevere, "CalDataMaker" ) << " unpack PEMD bank : invalid tower Key => Data probably corrupted" << endmsg;
      }
    }
  }

  //
  // PPR
  //
  EventRecord::ConstIterator p_pprdBank( anEvent,  "PPRD_StorableBank" );
  if(p_pprdBank.is_valid()){ 

    ConstHandle<RawQieStorableBank> bank_handle(p_pprdBank) ;
	    
    //  Loop over the bank we just found. 
    for(ConstGrandBankIterI2s1<RawQieStorableBank> data_iter(bank_handle);
	data_iter.is_valid(); ++data_iter) {
      
      //  convert address 
      calKey  = bank_handle->get_calKey( data_iter );

      if (calKey.validKey(calKey)) {
	cKey = CellKey(calKey);
	if ( !cKey.isValid() ) {
	  ERRLOG( ELsevere, "CalDataMaker" ) << " unpack PPRD bank : CellKey is inValid => Data probably corrupted" << endmsg;
	}
	else {

	  // For the PPR, no online LERs are applied. All corrections are done offline

	  float energy = (bank_handle->get_energy( data_iter )-_calib.pedPpr())* _calib.detPpr();

	  uint2 iData = bank_handle->get_data( data_iter );  
	  int   ipmt  = bank_handle->get_pmt(  data_iter );
          int   wedge = bank_handle->get_phi(  data_iter );

	  //
	  // WARNING hardwired 10 counts threshold (after talking to Willis)
	  // We need a better number and we need to put this in CalThreshold
	  // once we know what the number is
	  
	  if(thresh == false || (energy >0.012)){
	    CalTower* t = _calData->towerMake(cKey);
	    // apply offline LERs here only to energy not to iData, i.e. leave raw data untouched
	    energy *= pprCal[cKey.iEta()][2*wedge + ipmt];			// WKS:21Jun03
	    t->addPmtEnergy(calKey.detector(),ipmt,energy,iData);
	  }
	}
      } else { 
	ERRLOG( ELsevere, "CalDataMaker" ) << " unpack PPRD bank : invalid tower Key => Data probably corrupted" << endmsg;
      }      
    }  
  }
	  
	  
  //
  // CHA
  //
  EventRecord::ConstIterator p_chadBank( anEvent,  "CHAD_StorableBank" );	  
  if(p_chadBank.is_valid()){ 
    
    ConstHandle<RawQieStorableBank> bank_handle(p_chadBank) ;
    
    //  Loop over the bank we just found. 
	    
    for(ConstGrandBankIterI2s1<RawQieStorableBank> data_iter(bank_handle);
	data_iter.is_valid(); ++data_iter) {
	      
      //  convert address 
      
      calKey  = bank_handle->get_calKey( data_iter );

      if (calKey.validKey(calKey)) {
	cKey = CellKey(calKey);
	if ( !cKey.isValid() ) {
	  ERRLOG( ELsevere, "CalDataMaker" ) << " unpack CHAD bank : CellKey is inValid => Data probably corrupted" << endmsg;
	}
	else {
	  float energy = (bank_handle->get_energy( data_iter ) -_calib.pedCha())* _calib.detCha();

	  uint2 iData = bank_handle->get_data( data_iter ); 
	  int   ipmt  = bank_handle->get_pmt(  data_iter );
          int   wedge = bank_handle->get_phi(  data_iter );

	  if (TEPSEG[cKey.towerType()]==48) ipmt = 0;
	  
	  // Put the energy in CalData:
	  
	  if( thresh == false || 
             (energy > CHAMIN && !( (wedge==5 && ipmt==1) &&
                                    (cKey.iEta()>=30 && cKey.iEta()<=33) ) ) ) {
	    energy *= chaCal[cKey.iEta()][2*wedge + ipmt];		// WKS:21Jun03
	    CalTower* t = _calData->towerMake(cKey);
	    t->addPmtEnergy(calKey.detector(),ipmt,energy,iData); 
          }
	}
      } else {
	ERRLOG( ELsevere, "CalDataMaker" ) << " unpack CHAD bank : invalid tower Key => Data probably corrupted" << endmsg;
      }      
    }  
  }    
	  
	  
  //
  // WHA
  //

  EventRecord::ConstIterator p_whadBank( anEvent,  "WHAD_StorableBank" );
  if(p_whadBank.is_valid()){ 

    int  runNO = AbsEnv::instance()->runNumber();		// WKS: 30Jun03
    bool mcRun = AbsEnv::instance()->monteFlag();
    
    ConstHandle<RawQieStorableBank> bank_handle(p_whadBank) ;
	    
    //  Loop over the bank we just found. 
	    
    for(ConstGrandBankIterI2s1<RawQieStorableBank> data_iter(bank_handle);
	data_iter.is_valid(); ++data_iter) {
      
      //  convert address 
      
      calKey  = bank_handle->get_calKey( data_iter );
      
      if (calKey.validKey(calKey)) {
	cKey = CellKey(calKey);
	if ( !cKey.isValid() ) {
	  ERRLOG( ELsevere, "CalDataMaker" ) << " unpack WHAD bank : CellKey is inValid => Data probably corrupted" << endmsg;
	}
	else {
	  
	  float energy = (bank_handle->get_energy( data_iter ) -_calib.pedWha())* _calib.detWha();

	  uint2 iData = bank_handle->get_data( data_iter );  
	  int   ipmt  = bank_handle->get_pmt(  data_iter );
          int   wedge = bank_handle->get_phi(  data_iter );

	  if (TEPSEG[cKey.towerType()]==48) ipmt = 0;
	  
	  // Put the energy in CalData:
	  
	  if(thresh == false || (energy > WHAMIN )){
	
	    // Correct the west W11 6R <--> 7L swap for runs before the
	    // hardware fix.
	    //                                 Nikolay, RT     21Feb02

	    if ( !mcRun && runNO < 139608 ) {
	      if ( wedge == 11 && calKey.iew() == 0 ) {
		if      ( calKey.itower() == 6 && ipmt == 0 ) {
		  calKey.setKey(calKey.iew(),wedge,7,1,calKey.detectorName());
		  ipmt = 1;
		}
		else if ( calKey.itower() == 7 && ipmt == 1 ) {
		  calKey.setKey(calKey.iew(),wedge,6,0,calKey.detectorName());
		  ipmt = 0;
		}
	      }
	      cKey = CellKey(calKey);
	    }

	    // Correct the west W4  6R <--> 7R swap for runs before the 
	    // hardware fix.
	    //                                 Fotis, Nikolay  23Oct03

	    if ( !mcRun && runNO < 171885 ) {
	      if (ipmt==0 && wedge==4 && calKey.iew()==0){
		if (calKey.itower()==6) {
		  calKey.setKey(calKey.iew(),wedge ,7 , ipmt , calKey.detectorName());
		} else if (calKey.itower()==7) {
		  calKey.setKey(calKey.iew(),wedge ,6 , ipmt , calKey.detectorName());
		} 
	      } 
	      cKey = CellKey(calKey);
	    }
	    CalTower* t = _calData->towerMake(cKey);
	    energy *= whaCal[cKey.iEta()][2*wedge + ipmt];		// WKS:21Jun03
	    t->addPmtEnergy(calKey.detector(),ipmt,energy,iData); 
	  }
	}
      } else {
	ERRLOG( ELsevere, "CalDataMaker" ) << " unpack WHAD bank : invalid tower Key => Data probably corrupted" << endmsg;
      }
    }
  }  
	      
	  
  //
  // PHA
  //

  EventRecord::ConstIterator p_phadBank( anEvent,  "PHAD_StorableBank" );
  if(p_phadBank.is_valid()){ 
	    
    ConstHandle<RawQieStorableBank> bank_handle(p_phadBank) ;
    
    //  Loop over the bank we just found. 
    
    for(ConstGrandBankIterI2s1<RawQieStorableBank> data_iter(bank_handle);
	data_iter.is_valid(); ++data_iter) {
	      
      //  convert address 
      
      calKey  = bank_handle->get_calKey( data_iter );

      if (calKey.validKey(calKey)) {
	cKey = CellKey(calKey);
	if ( !cKey.isValid() ) {
	  ERRLOG( ELsevere, "CalDataMaker" ) << " unpack PHAD bank : CellKey is inValid => Data probably corrupted" << endmsg;
	}
	else {	
	  float energy = (bank_handle->get_energy( data_iter )-_calib.pedPha())* _calib.detPha();

	  uint2 iData = uint2(bank_handle->get_data( data_iter )); 
	  int   ipmt  = bank_handle->get_pmt(  data_iter );
          int   wedge = bank_handle->get_phi(  data_iter );
	  
	  // Put the energy in CalData:
	  
	  if(thresh == false || (energy > PHAMIN )){
	    CalTower* t = _calData->towerMake(cKey);
	    energy *= phaCal[cKey.iEta()][2*wedge + ipmt];		// WKS:21Jun03
	    t->addPmtEnergy(calKey.detector(),ipmt,energy,iData); 
	  }
	}
      } else {
	ERRLOG( ELsevere, "CalDataMaker" ) << " unpack PHAD bank : invalid tower Key => Data probably corrupted" << endmsg;
      }      
    }   
  }

  // MiniPlug
  float ilum = 2.;
  if (useMPALumiCalib  && _calib.mpaOffset(0,0) != -9999) {
    EventRecord::ConstIterator clld_iter(anEvent, "CLLD_StorableBank");
    if (clld_iter.is_valid()) {
    ConstHandle<CLLD_StorableBank> clld(clld_iter);
    ilum = (clld->get_l_inst())/1.0e30;
    }  
    // don't print the error for MC
    else if (!gblEnv->monteFlag()) ERRLOG( ELsevere, "CalDataMaker" ) << " No inst lumi "  << endmsg;
    _calib.setDetMpa(ilum);
  }

  const Float_t TOW2DEP_CORR[14] = {
    0.2672,0.2548,
    0.2137,0.2179,0.2198,
    0.2309,0.2343,0.2277,0.2399,
    0.2015,0.2496,0.2472,0.2492,0.2495
  };
  const Double_t MC_EXP1_NEW[4] = {-0.099, -0.156, -0.171, -0.199};

  EventRecord::ConstIterator p_mpadBank( anEvent,  "MPAD_StorableBank" );
  if(p_mpadBank.is_valid()){ 
    ConstHandle<MPAD_StorableBank> bank_handle(p_mpadBank) ;
            
    //  Loop over the bank we just found. 
            
    for(ConstGrandBankIterI2s1<MPAD_StorableBank> data_iter(bank_handle);
        data_iter.is_valid(); ++data_iter) {
      // ignore sum towers
      int itype = bank_handle->get_type( data_iter );
      if (itype == 1) {
	//  convert address  
	calKey  = bank_handle->get_calKey( data_iter );
	if (calKey.validKey(calKey)) {
	  cKey = CellKey(calKey);
	  if ( !cKey.isValid() ) {
	    ERRLOG( ELsevere, "CalDataMaker" ) << " unpack MPAD bank : CellKey is inValid => Data probably corrupted" << endmsg;
	  }
	  else {
	    float energy = 0.;
	    if (useMPALumiCalib && _calib.mpaOffset(0,0) != -9999) {
	      int is =  bank_handle->get_we( data_iter );
	      int itow = bank_handle->get_sum_or_tower( data_iter );
	      int layer = 25 - calKey.itower();
	      
	      float mpCorr = (_calib.mpaOffset(is,itow) + _calib.mpaSlope(is,itow)*ilum) /(MC_EXP1_NEW[layer]); 
	      mpCorr = mpCorr*TOW2DEP_CORR[itow%14];
	      if (mpCorr < 0.) mpCorr = 0.;
	      int data = bank_handle->get_energy( data_iter );
	      energy = (data - _calib.pedMpa())*mpCorr;
	      if (energy < 0.) energy = 0.;                     // to correct for pedestal width
	    }
	    else {
	      energy = (bank_handle->get_energy( data_iter ) -_calib.pedMpa())* _calib.detMpa();
	      if (energy < 0.) energy = 0.; 
	      //          float energy = (bank_handle->get_energy( data_iter ) - 200)* 0.001;
	    }
	    uint2 iData = bank_handle->get_data( data_iter ); 
	    int ipmt = calKey.ipmt();
	    if (TEPSEG[cKey.towerType()]==48) ipmt = 0;
	    
	    // Put the energy in CalData:
	    if(thresh == false || (energy >MPAMIN)){
	      CalTower* t = _calData->towerMake(cKey);
	      t->addPmtEnergy(calKey.detector(),ipmt,energy,iData);
	    }
	  }
	} else {
	  ERRLOG( ELsevere, "CalDataMaker" ) << " unpack MPAD bank : invalid tower Key => Data probably corrupted" << endmsg;
	}
      }
    }
  }


  //
  // HAT (hadron TDC banks) and EMTD (EM TDC bank)

  int1 nhit_tdc_cha[TENETA][TENPHI];
  int1 nhit_tdc_wha[TENETA][TENPHI];
  int1 nhit_tdc_pha[TENETA][TENPHI];
  int1 nhit_tdc_cem[TENETA][TENPHI];
  int1 nhit_tdc_pem[TENETA][TENPHI];

  Int_t ihit_tdc_cha[TENETA][TENPHI][CAL_TDC_NHITS];
  Int_t ihit_tdc_wha[TENETA][TENPHI][CAL_TDC_NHITS];
  Int_t ihit_tdc_pha[TENETA][TENPHI][CAL_TDC_NHITS];
  Int_t ihit_tdc_cem[TENETA][TENPHI][CAL_TDC_NHITS];
  Int_t ihit_tdc_pem[TENETA][TENPHI][CAL_TDC_NHITS];

  int det_tdc[TENETA][TENPHI];
  memset( det_tdc,0,TENETA*TENPHI*sizeof(int));

  memset(nhit_tdc_cha,0,TENETA*TENPHI*sizeof(int1));
  memset(nhit_tdc_wha,0,TENETA*TENPHI*sizeof(int1));
  memset(nhit_tdc_pha,0,TENETA*TENPHI*sizeof(int1));
  memset(nhit_tdc_cem,0,TENETA*TENPHI*sizeof(int1));
  memset(nhit_tdc_pem,0,TENETA*TENPHI*sizeof(int1));

  memset(ihit_tdc_cha,0,TENETA*TENPHI*CAL_TDC_NHITS*sizeof(Int_t));
  memset(ihit_tdc_wha,0,TENETA*TENPHI*CAL_TDC_NHITS*sizeof(Int_t));
  memset(ihit_tdc_pha,0,TENETA*TENPHI*CAL_TDC_NHITS*sizeof(Int_t));
  memset(ihit_tdc_cem,0,TENETA*TENPHI*CAL_TDC_NHITS*sizeof(Int_t));
  memset(ihit_tdc_pem,0,TENETA*TENPHI*CAL_TDC_NHITS*sizeof(Int_t)); 

  // Locate a EMTD bank in the event record
  EventRecord::ConstIterator iEMTD( anEvent, "EMTD_StorableBank");
  // Construct a handle to the EMTD bank from the iterator.  
  ConstHandle<EMTD_StorableBank> EMTD;
  bool haveEMTD = false;
  if (! iEMTD.is_valid()){
    //try a plain bank
    for (EventRecord::ConstIterator jEMTD(anEvent,"StorableBank") ;
         jEMTD.is_valid() ; ++jEMTD){
      ConstHandle<StorableBank> handle(jEMTD);
      if (handle->bank_name() == std::string("EMTD")){
        EMTD.assign(new EMTD_StorableBank((const EMTD_StorableBank&)(handle.operator*())));    
	haveEMTD = true;
      }
    }
  } else {
    haveEMTD = true;
    EMTD.assign(iEMTD);
  }
  if ( haveEMTD ){
      
    // Loop over the blocks in the D bank.
    for (ConstBankIterTDC block(EMTD); block.is_valid(); ++block) {
      
      // Loop over the TDC modules in Block
      for (ConstBlockIterTDC module(block) ; module.is_valid() ; ++module) {
	
	// Loop over hits in each TDC module
	for (ConstModuleIterTDC hit(module) ;hit.is_valid() ; ++hit) {
	  
	  // Unpack the hit information (extract the bitfields)
	  // Same as for the header word.  Note the argument is "hit". 
	  int channel = EMTD->get_channel_ID(hit);  // all 9 bits of it
	  int leadingEdge = EMTD->get_start(hit);
	  int firstFlag = EMTD->get_first_flag(hit);
	  int side, wedge, rapidity, ierror;
	  
	  EMTD->get_geometric_indices(hit,side,wedge,rapidity,ierror);
	  Detector subcal = EMTD->get_detector(hit);     
	  TowerKey twKey = EMTD->getTowerKey(hit,ierror);
	  if (ierror != 0){
	    ERRLOG(ELwarning, "CalDataMaker")
	      << "failed to get TowerKey for EMTD bank hit" << endmsg;
	    continue;
	  }
	  size_t ieta=twKey.iEta();
	  size_t iphi=twKey.iPhi();
	  
	  int detindex = det_tdc[ieta][iphi];
	  det_tdc[ieta][iphi] = detindex | (1 << subcal);
	  if (twKey.towerType()==5 && iphi%2 == 0 && iphi+1 < TENPHI)  
	    det_tdc[ieta][iphi+1] = detindex | (1 << subcal);
	  
	  int nhitcem=nhit_tdc_cem[ieta][iphi];
	  int nhitpem=nhit_tdc_pem[ieta][iphi];
	  if (subcal==CEM) {
	    if (nhitcem<CAL_TDC_NHITS) {
	      ihit_tdc_cem[ieta][iphi][nhitcem]=leadingEdge;
	    } else {
	      ERRLOG( ELsevere, "CalDataMaker" ) << 
		" unpack EMTD bank : more than 8 hits in CEM " << endmsg;
	    }
	    nhit_tdc_cem[ieta][iphi]++;
            if (nhit_tdc_cem[ieta][iphi]>CAL_TDC_NHITS) nhit_tdc_cem[ieta][iphi]=CAL_TDC_NHITS;
	  } else if (subcal==PEM) {
	    if (nhitpem<CAL_TDC_NHITS) {
	      ihit_tdc_pem[ieta][iphi][nhitpem]=leadingEdge;
	      //Fill a copy to the odd phi tower in type5; with some sanity check
	      if (twKey.towerType()==5 && iphi%2 == 0 && iphi+1 < TENPHI)  
		ihit_tdc_pem[ieta][iphi+1][nhitpem]=leadingEdge;
	    } else {
	      ERRLOG( ELsevere, "CalDataMaker" ) << 
		" unpack EMTD bank : more than 8 hits in PEM " << endmsg;
	    }
	    nhit_tdc_pem[ieta][iphi]++;
            if (nhit_tdc_pem[ieta][iphi]>CAL_TDC_NHITS) nhit_tdc_pem[ieta][iphi]=CAL_TDC_NHITS;

	    //Fill a copy to the odd phi tower in type5; with some sanity check
	    if (twKey.towerType()==5 && iphi%2 == 0 && iphi+1 < TENPHI)  
	      nhit_tdc_pem[ieta][iphi+1]=nhit_tdc_pem[ieta][iphi];
	    
	  }
        }  // end loop over hits in TDC
      }  // end loop over TDCs in a block
    }  // end loop over blocks in bank
  }
  
  // Locate a HATD bank in the event record
  EventRecord::ConstIterator iHATD( anEvent, "HATD_StorableBank");
  if ( iHATD.is_valid( ) ){
    
     // Construct a handle to the HATD bank from the iterator.  

     ConstHandle<HATD_StorableBank> HATD( iHATD );
     
     // Loop over the blocks in the D bank.

     for (ConstBankIterTDC block(HATD); block.is_valid(); ++block) {
       
       
       // Loop over the TDC modules in Block
       
       for (ConstBlockIterTDC module(block) ; module.is_valid() ; ++module) {
	 
	 // Loop over hits in each TDC module
	 
	 for (ConstModuleIterTDC hit(module) ;hit.is_valid() ; ++hit) {
	  
	  // Unpack the hit information (extract the bitfields)
	  // Same as for the header word.  Note the argument is "hit". 
	  int channel = HATD->get_channel_ID(hit);  // all 9 bits of it
	  int leadingEdge = HATD->get_start(hit);
	  int firstFlag = HATD->get_first_flag(hit);
	  int side, wedge, rapidity, depth_or_tile, ierror;

	  HATD->get_geometric_indices(hit,
				      side,
				      wedge,
				      rapidity,
				      depth_or_tile,
				      ierror);

	  TdcKey tdccha(side, wedge, rapidity, depth_or_tile);
	  Detector subcal = tdccha.detector();     

	  size_t ieta=tdccha.ixEta();
	  size_t iphi=tdccha.ixPhi();

	  int detindex = det_tdc[ieta][iphi];
	  det_tdc[ieta][iphi] = detindex | (1 << subcal);

	  int nhitcha=nhit_tdc_cha[ieta][iphi];
	  int nhitwha=nhit_tdc_wha[ieta][iphi];
	  int nhitpha=nhit_tdc_pha[ieta][iphi];
	  if (subcal==CHA) {
	    if (nhitcha<CAL_TDC_NHITS) {
	      ihit_tdc_cha[ieta][iphi][nhitcha]=leadingEdge;
	    } else {
	      ERRLOG( ELsevere, "CalDataMaker" ) << 
		" unpack HATD bank : more than "<<CAL_TDC_NHITS<<" hits in CHA " << endmsg;
	    }
	    nhit_tdc_cha[ieta][iphi]++;
            if (nhit_tdc_cha[ieta][iphi]>CAL_TDC_NHITS) nhit_tdc_cha[ieta][iphi]=CAL_TDC_NHITS;
	  } else if (subcal==WHA) {
	    if (nhitwha<CAL_TDC_NHITS) {
	      ihit_tdc_wha[ieta][iphi][nhitwha]=leadingEdge;
	    } else {
	      ERRLOG( ELsevere, "CalDataMaker" ) << 
		" unpack HATD bank : more than "<<CAL_TDC_NHITS<<" hits in WHA " << endmsg;
	    }
	    nhit_tdc_wha[ieta][iphi]++;
            if (nhit_tdc_wha[ieta][iphi]>CAL_TDC_NHITS) nhit_tdc_wha[ieta][iphi]=CAL_TDC_NHITS;
	  } else if (subcal==PHA) {
	    if (nhitpha<CAL_TDC_NHITS) {
	      ihit_tdc_pha[ieta][iphi][nhitpha]=leadingEdge;
	    } else {
	      ERRLOG( ELsevere, "CalDataMaker" ) 
		<< " unpack HATD bank : more than "<<CAL_TDC_NHITS<<" hits in PHA " << endmsg;
	    }
	    nhit_tdc_pha[ieta][iphi]++;
            if (nhit_tdc_pha[ieta][iphi]>CAL_TDC_NHITS) nhit_tdc_pha[ieta][iphi]=CAL_TDC_NHITS;
	  }
        }  // end loop over hits in TDC
      }  // end loop over TDCs in a block
    }  // end loop over blocks in bank
  }


  // apply the slewing corrections, and T0 calibration and store them in CalData
  for (int ieta=0; ieta<TENETA; ieta++) {
    for (int iphi=0; iphi<TENPHI; iphi++) {
      // get detector
      for(int detindex=0; detindex < 5; detindex++){
	
	Detector Det=Detector(detindex);
	cKey = CellKey(ieta,iphi,Det);
      	// check whether tower exists already. Else we don't store the timing information
	CalTower* t = _calData->tower(ieta,iphi);
	if (t>0) {
	  
	  // set ValidCode

	  //cem
	  //if (cKey.detector()==CEM){
	    //	    float effic = _calib.efficTdcCEM(ieta,iphi);
	    //      if (effic > 2000.)  t->setBad(CEMTDC);
	  //}
	  
	  //pem
	  //if(cKey.detector() == PEM){
	    //	    float effic = _calib.efficTdcPEM(ieta,iphi);
	    //	    if (effic > 2000.)  t->setBad(PEMTDC);
	  //}
	  
	  
	  // hadronic
	  for (int ihad = 0; ihad != cKey.numHadDepth(); ++ihad) {
	  	    
	    if (cKey.hadDetector(ihad)==CHA && cKey.detector()==CHA){
	
	      float effic = _calib.efficTdcCHA(ieta,iphi);
	      if (effic > 2000.) t->setBad(CHATDC);
	 
	    } else if (cKey.hadDetector(ihad)==WHA && cKey.detector()==WHA) {

	      float effic = _calib.efficTdcWHA(ieta,iphi);
	      if (effic > 2000.) t->setBad(WHATDC);
	     	      
	    } else if (cKey.hadDetector(ihad)==PHA && cKey.detector()==PHA){
	    
	      float effic = _calib.efficTdcPHA(ieta,iphi);
	      if (effic > 2000.)  t->setBad(PHATDC);
	
	    }
	  }
	  
	  int result =(1 << detindex) & det_tdc[ieta][iphi];
	  if(! (result > 0) ) continue; 

	  if (nhit_tdc_cha[ieta][iphi]+nhit_tdc_wha[ieta][iphi]+nhit_tdc_pha[ieta][iphi]+
	      nhit_tdc_cem[ieta][iphi]+nhit_tdc_pem[ieta][iphi]>0) {
	    float pmt0=t->pmtEnergy(cKey,0);
	    float pmt1=t->pmtEnergy(cKey,1);
	    int iht=0;
	    if (cKey.detector()==CHA){
	      t->setTDCHits(cKey,nhit_tdc_cha[ieta][iphi]);
	      // store them in reverse order which is then the real time order
	      for (int ih=nhit_tdc_cha[ieta][iphi]-1;ih>-1; ih--) {
		float time   = 0;
		float slew = _calib.slewTdcCHA(ieta,iphi);
		float t0   = _calib.t0TdcCHA(ieta,iphi);
		bool bpmt0(_calib.badChannelsCHA(ieta, 2*iphi));
		bool bpmt1(_calib.badChannelsCHA(ieta, 2*iphi+1));
		float energy( 0.0 );
		if ( !bpmt0 && !bpmt1 ) {
		  energy = sqrt(pmt0*pmt1);
		  if (energy == 0) energy = pmt0+pmt1; // 14-Jan-04 If one PMT is zero, use other.
		}
		else if ( !bpmt0 && bpmt1 ) energy = pmt0;
		else if ( !bpmt1 && bpmt0 ) energy = pmt1;
	
		if (energy>0) {
		  time=float(ihit_tdc_cha[ieta][iphi][ih])-t0-slew/sqrt(energy)+float(HTDC_Offset);
		} else {
		  time=ihit_tdc_cha[ieta][iphi][ih]-t0+HTDC_Offset;
		}
		t->setTime(cKey,int(time*10),iht);
		iht++;
	      }
	    } else if (cKey.detector()==WHA){
	      int iht=0;
	      t->setTDCHits(cKey,nhit_tdc_wha[ieta][iphi]);
	      // store them in reverse order which is then the real time order
	      for (int ih=nhit_tdc_wha[ieta][iphi]-1;ih>-1; ih--) {
		float time   = 0;
		float slew = _calib.slewTdcWHA(ieta,iphi);
		float t0     = _calib.t0TdcWHA(ieta,iphi);
		bool bpmt0(_calib.badChannelsWHA(ieta, 2*iphi));
		bool bpmt1(_calib.badChannelsWHA(ieta, 2*iphi+1));
		float energy( 0.0 );
		if ( !bpmt0 && !bpmt1 ) {
		  energy = sqrt(pmt0*pmt1);
		  if (energy == 0) energy = pmt0+pmt1; // 14-Jan-04 If one PMT is zero, use other.
		}
		else if ( !bpmt0 && bpmt1 ) energy = pmt0;
		else if ( !bpmt1 && bpmt0 ) energy = pmt1;

		if (energy>0){
		  time=float(ihit_tdc_wha[ieta][iphi][ih])-t0-slew/sqrt(energy)+float(HTDC_Offset);
		} else {
		  time=float(ihit_tdc_wha[ieta][iphi][ih])-t0+float(HTDC_Offset);
		}
		t->setTime(cKey,int(time*10),iht);
		iht++;
	      }
	    } else if (cKey.detector()==PHA){
	      int iht=0;
	      t->setTDCHits(cKey,nhit_tdc_pha[ieta][iphi]);
	      // store them in reverse order which is then the real time order
	      for (int ih=nhit_tdc_pha[ieta][iphi]-1;ih>-1; ih--) {
		float time   = 0;
		float slew = _calib.slewTdcPHA(ieta,iphi);
		float t0     = _calib.t0TdcPHA(ieta,iphi);
		int TType( TOWER_TYPE[ieta] );
		bool bpmt0( false );
		if ( TType == 5 )
		  bpmt0 = _calib.badChannelsPHA(ieta,  iphi  );
		else
		  bpmt0 = _calib.badChannelsPHA(ieta, 2*iphi  );

		float energy( 0.0 );
		if ( !bpmt0 ) energy = pmt0;

		if (energy>0) {
		  time=float(ihit_tdc_pha[ieta][iphi][ih])-t0-slew/sqrt(energy)+float(HTDC_Offset);
		} else {
		  time=float(ihit_tdc_pha[ieta][iphi][ih])-t0+float(HTDC_Offset);
		}
		t->setTime(cKey,int(time*10),iht);
		iht++;
	      }	    
	    } else if (cKey.detector()==CEM){
	      int iht=0;
	      t->setTDCHits(cKey,nhit_tdc_cem[ieta][iphi]);
	      // store them in reverse order which is then the real time order
	      for (int ih=nhit_tdc_cem[ieta][iphi]-1;ih>-1; ih--) {
		float time   = 0;
// 		float slew = _calib.slewTdcCEM(ieta,iphi);
// 		float t0     = _calib.t0TdcCEM(ieta,iphi);
		float slew = 0.;
		float t0 = -10000.;
		bool bpmt0(_calib.badChannelsCEM(ieta, 2*iphi));
		bool bpmt1(_calib.badChannelsCEM(ieta, 2*iphi+1));
		float energy( 0.0 );
		if ( !bpmt0 && !bpmt1 ) {
		  energy = sqrt(pmt0*pmt1);
		  if (energy == 0) energy = pmt0+pmt1; // 14-Jan-04 If one PMT is zero, use other.
		}
		else if ( !bpmt0 && bpmt1 ) energy = pmt0;
		else if ( !bpmt1 && bpmt0 ) energy = pmt1;

		if (energy>0) {
		  time=float(ihit_tdc_cem[ieta][iphi][ih])-t0-slew/sqrt(energy)+float(EMTDC_Offset);
		} else {
		  time=float(ihit_tdc_cem[ieta][iphi][ih])-t0+float(EMTDC_Offset);
		}
		t->setTime(cKey,int(time*10),iht);
		iht++;
	      }	    
	    } else if (cKey.detector()==PEM){
	      int iht=0;
	      t->setTDCHits(cKey,nhit_tdc_pem[ieta][iphi]);
	      // store them in reverse order which is then the real time order
	      for (int ih=nhit_tdc_pem[ieta][iphi]-1;ih>-1; ih--) {
		float time   = 0;
// 		float slew = _calib.slewTdcPEM(ieta,iphi);
// 		float t0     = _calib.t0TdcPEM(ieta,iphi);
		float slew = 0.;
		float t0 = -10000.;
		int TType( TOWER_TYPE[ieta] );
		bool bpmt0( false ), bpmt1( false );
		if ( TType < 5 ) {
		  bpmt0 = _calib.badChannelsPEM(ieta, 2*iphi  );
		  bpmt1 = _calib.badChannelsPEM(ieta, 2*iphi+1);
		}
		else if ( TType == 5 )
		  bpmt0 = _calib.badChannelsPEM(ieta,  iphi  );
		else
		  bpmt0 = _calib.badChannelsPEM(ieta, 2*iphi  );

		float energy( 0.0 );
		if ( !bpmt0 && !bpmt1 )  {
		  energy = sqrt(pmt0*pmt1);
		  if (energy == 0) energy = pmt0+pmt1; // 14-Jan-04 If one PMT is zero, use other.
		}
		else if ( !bpmt0 && bpmt1 ) energy = pmt0;
		else if ( !bpmt1 && bpmt0 ) energy = pmt1;
	
		if (energy>0) {
		  time=float(ihit_tdc_pem[ieta][iphi][ih])-t0-slew/sqrt(energy)+float(EMTDC_Offset);
		} else {
		  time=float(ihit_tdc_pem[ieta][iphi][ih])-t0+float(EMTDC_Offset);
		}
		t->setTime(cKey,int(time*10),iht);
		iht++;
	      }
	    }
	  }
	} 
      }
    }
  }


  //store bad chanel info. in validcode
  for (int ieta=0; ieta<TENETA; ieta++) {
    for (int iphi=0; iphi<TENPHI; iphi++) {
      // get detector
      for(int detindex=0; detindex < 5; detindex++){
	
	Detector Det=Detector(detindex);
	cKey = CellKey(ieta,iphi,Det);
      	// check whether tower exists already. 
	CalTower* t = _calData->tower(ieta,iphi);
	if (t>0) {
	  
	  //cem
	  if (cKey.detector()==CEM){
	    bool bpmt0(_calib.badChannelsCEM(ieta, 2*iphi));
	    bool bpmt1(_calib.badChannelsCEM(ieta, 2*iphi+1));
	    if( bpmt0)  t->setBad(CEMPMT0BAD);
	    if( bpmt1)  t->setBad(CEMPMT1BAD);
	  }

	  //pem
	  if (cKey.detector()==PEM){
	    bool bpmt0( false ), bpmt1( false );
	    int TType( TOWER_TYPE[ieta] );
	    if ( TType < 5 ) {
	      bpmt0 = _calib.badChannelsPEM(ieta, 2*iphi  );
	      bpmt1 = _calib.badChannelsPEM(ieta, 2*iphi+1);
	    }
	    else if ( TType == 5 )
	      bpmt0 = _calib.badChannelsPEM(ieta,  iphi  );
	    else
	      bpmt0 = _calib.badChannelsPEM(ieta, 2*iphi  );
	    if( bpmt0)  t->setBad(PEMPMT0BAD);
	    if( bpmt1)  t->setBad(PEMPMT1BAD);
	    
	  }
	  
	  
	  // hadronic
	  for (int ihad = 0; ihad != cKey.numHadDepth(); ++ihad) {
	    if (cKey.hadDetector(ihad)==CHA && cKey.detector()==CHA){
	      bool bpmt0(_calib.badChannelsCHA(ieta, 2*iphi));
	      bool bpmt1(_calib.badChannelsCHA(ieta, 2*iphi+1));
	      if( bpmt0)  t->setBad(CHAPMT0BAD);
	      if( bpmt1)  t->setBad(CHAPMT1BAD);
	    }

	    if (cKey.hadDetector(ihad)==WHA && cKey.detector()==WHA) {
	      bool bpmt0(_calib.badChannelsWHA(ieta, 2*iphi));
	      bool bpmt1(_calib.badChannelsWHA(ieta, 2*iphi+1));
	      if( bpmt0)  t->setBad(WHAPMT0BAD);
	      if( bpmt1)  t->setBad(WHAPMT1BAD);
	    }

	    if (cKey.hadDetector(ihad)==PHA && cKey.detector()==PHA) {
    	      bool bpmt0( false ), bpmt1( false );
	      int TType( TOWER_TYPE[ieta] );
	      if ( TType < 5 ) {
		bpmt0 = _calib.badChannelsPHA(ieta, 2*iphi  );
		bpmt1 = _calib.badChannelsPHA(ieta, 2*iphi+1);
	      }
	      else if ( TType == 5 )
		bpmt0 = _calib.badChannelsPHA(ieta,  iphi  );
	      else
		bpmt0 = _calib.badChannelsPHA(ieta, 2*iphi  );
	      if( bpmt0)  t->setBad(PHAPMT0BAD);
	      if( bpmt1)  t->setBad(PHAPMT1BAD);
	    }

	  }
	}
      }
    }
  }



     //   }
  return 1;
  
}

//-------------------------------------------------------------------
//   AddSimEnergy:
//   add simulation energy in a tower 
//
//-------------------------------------------------------------------

void CalDataMaker::addSimEnergy(AbsEvent* anEvent,TowerKey& key, Detector tag, float energy0, float energy1){

  switch ( tag ) {
  case CEM:
      
      if((energy0 > CEMMIN || energy1 >CEMMIN) && !(key.iPhi()==5 && key.iEta()>33 && key.iEta()<36 )){
	
	// Put the energy into Caldata:
	CalTower* t = _calData->towerMake(key);
	uint2 data0 = RawQieStorableBank::make_packedData(_calib.pedCem()+(energy0 / _calib.detCem())); 
	uint2 data1 = RawQieStorableBank::make_packedData(_calib.pedCem()+(energy1 / _calib.detCem())); 
	t->addPmtEnergy(CalorDetInfo::CEM, 0, energy0, data0);
	t->addPmtEnergy(CalorDetInfo::CEM, 1, energy1, data1);
      }
      break;
    case PEM:
      if (energy0 > PEMMIN || energy1 > PEMMIN) {
	// Put the energy into Caldata:
	CalTower* t = _calData->towerMake(key);
	uint2 data0 = RawQieStorableBank::make_packedData(_calib.pedPem()+(energy0 / _calib.detPem()));
 	uint2 data1 = RawQieStorableBank::make_packedData(_calib.pedPem()+(energy1 / _calib.detPem()));
	int ieta=int(key.iEta());
	if (TOWER_TYPE[ieta]==3 || TOWER_TYPE[ieta]==4) {	// Towers 0,1,2,3 
	  t->addPmtEnergy(CalorDetInfo::PEM, 0, energy0, data0);
	  t->addPmtEnergy(CalorDetInfo::PEM, 1, energy1, data1);
	} else {
	  t->addPmtEnergy(CalorDetInfo::PEM, 0, energy0, data0);
	}
      }     
      break;
    case CHA:
		
      if((energy0 > CHAMIN || energy1 > CHAMIN) && !(key.iPhi()==5 && key.iEta()>33 && key.iEta()<36 )){
	
	// Put the energy into Caldata:
	CalTower* t = _calData->towerMake(key);
	uint2 data0 = RawQieStorableBank::make_packedData(_calib.pedCha()+(energy0 / _calib.detCha())); 
	uint2 data1 = RawQieStorableBank::make_packedData(_calib.pedCha()+(energy1 / _calib.detCha())); 
	t->addPmtEnergy(CalorDetInfo::CHA, 0, energy0, data0);
	t->addPmtEnergy(CalorDetInfo::CHA, 1, energy1, data1);
      }     
      break;
    case WHA:
		
      if(energy0 > WHAMIN || energy1 > WHAMIN){
	
	// Put the energy into Caldata:
	CalTower* t = _calData->towerMake(key);
	uint2 data0 = RawQieStorableBank::make_packedData(_calib.pedWha()+(energy0 / _calib.detWha())); 
	uint2 data1 = RawQieStorableBank::make_packedData(_calib.pedWha()+(energy1 / _calib.detWha())); 
	t->addPmtEnergy(CalorDetInfo::WHA, 0, energy0, data0);
	t->addPmtEnergy(CalorDetInfo::WHA, 1, energy1, data1);
      }     
      break;
    case PHA:
		
      if (energy0 > PHAMIN || energy1 > PHAMIN) {	        // WKS:21Jun03
	// Put the energy into Caldata:
	CalTower* t = _calData->towerMake(key);
	uint2 data0 = RawQieStorableBank::make_packedData(_calib.pedPha()+(energy0 / _calib.detPha()));
 	uint2 data1 = RawQieStorableBank::make_packedData(_calib.pedPha()+(energy1 / _calib.detPha()));
	int ieta=int(key.iEta());
	if (TOWER_TYPE[ieta]==4) { 			        // Towers 2,3
	  t->addPmtEnergy(CalorDetInfo::PHA, 0, energy0, data0);
	  t->addPmtEnergy(CalorDetInfo::PHA, 1, energy1, data1);
	} else {
	  t->addPmtEnergy(CalorDetInfo::PHA, 0, energy0, data0);
	}
      }     
      break;
    case PPR:

      if (energy0 > 0 || energy1 > 0) {			        // WKS:21Jun03
	// Put the energy into Caldata:
	CalTower* t = _calData->towerMake(key);
	uint2 data0 = RawQieStorableBank::make_packedData(_calib.pedPpr()+(energy0 / _calib.detPpr()));
 	uint2 data1 = RawQieStorableBank::make_packedData(_calib.pedPpr()+(energy1 / _calib.detPpr()));
	int ieta=int(key.iEta());
	if (TOWER_TYPE[ieta]==3 || TOWER_TYPE[ieta]==4) {       // Towers 0,1,2,3
	  t->addPmtEnergy(CalorDetInfo::PPR, 0, energy0, data0);
	  t->addPmtEnergy(CalorDetInfo::PPR, 1, energy1, data1);
	} else {
	  t->addPmtEnergy(CalorDetInfo::PPR, 0, energy0, data0);
	}
      }     
      break;
    case MPA:

      if(energy0 > MPAMIN){
        // Put the energy into Caldata:
        CalTower* t = _calData->towerMake(key);
        uint2 data0 = RawQieStorableBank::make_packedData(_calib.pedMpa() +(energy0 / _calib.detMpa()));
        t->addPmtEnergy(CalorDetInfo::MPA, 0, energy0, data0);
      }
      break;
      //    default:
      //          ERRLOG(ELerror,"CalDataMaker: ") << " Detector tag unknown: " << tag << " " << key <<endmsg;
    }
}


// Make D banks

int CalDataMaker::makeDBanks(AbsEvent* anEvent,bool verbose) {

  CalData_ch caldata;

  if(_calData.is_null()){
  
    //  get CalData from the event. 
    
    CalData::Error result = CalData::find(caldata);
  
    if(result != CalData::OK){
       ERRLOG(ELerror,"CalDataMaker: ") << "could not find CalData in the event" << endmsg;
      return 0 ;
    }
  }
  else{
    caldata.assign(&(*_calData));
  }

  //
  // Create a the D Bank classes:
  Handle<CEMD_StorableBank> cemd(new CEMD_StorableBank) ;
  int cemChanCount = 0;
  Handle<PEMD_StorableBank> pemd(new PEMD_StorableBank) ;
  int pemChanCount = 0;
  Handle<CHAD_StorableBank> chad(new CHAD_StorableBank) ;
  int chaChanCount = 0;
  Handle<WHAD_StorableBank> whad(new WHAD_StorableBank) ;
  int whaChanCount = 0;
  Handle<PHAD_StorableBank> phad(new PHAD_StorableBank) ;
  int phaChanCount = 0;
  Handle<PPRD_StorableBank> pprd(new PPRD_StorableBank) ;
  int pprChanCount = 0;
  Handle<MPAD_StorableBank> mpad(new MPAD_StorableBank) ;
  int mpaChanCount = 0;
  // this will be a list of calorimeters for a given tower
  DetList caloList;

  //
  // Loop over the towers in CalData
  //
 for(CalData::ConstIterator tower= caldata->begin(); tower != caldata->end(); ++tower){ 
 
    tower->detectors(caloList);    // get detector list

    // Loop over detectors (calorimeters) in the tower
    for(DetList::iterator calIter = caloList.begin(); calIter != caloList.end(); ++calIter){

        if(*calIter != MPA){    
	CellKey cellKey(tower->iEta(),tower->iPhi(),*calIter);
	// loop over photo-tubes
	int theMaxPmt=1;
	if(cellKey.towerType()<=4) theMaxPmt=2;
	for (int ipmt=0; ipmt<theMaxPmt; ++ipmt){
	  CalorKey calKey;
	  calKey = cellKey.makeCalorKey(cellKey,ipmt);
	  if ( !calKey.isValid() && cellKey.towerType()!=8 ) {
	    std::cout << "CalDataMaker::makeDBanks() CalorKey is inValid:" << std::endl;
	    calKey.print();
	    std::cout << "CellKey print:" << std::endl;
	    std::cout << cellKey ;
	    std::cout << "detector " << *calIter << " ipmt " << ipmt << std::endl;
	  }
	  else{

	    // CEM
	    
	    if (calKey.detector()==CalorDetInfo::CEM){	    
	      uint2 data = tower->pmtData(cellKey,ipmt);
	      // add the channel id and data to bank 
	      if(data > 0){
		uint2 chid = cemd->make_chid(calKey); 
		if (cemChanCount<cemd->get_max_size()/2) {
		  cemd->add_datum( chid, data ); 
		  cemChanCount++;
		}
		else {
		      ERRLOG(ELsevere,"CalDataMaker: ") << "makeDBanks more CEMD data than maximum number of channels allowed" << endmsg;
		}
	      }
	    }   // endif CEM
	 
	    // PEM
	    
	    if (calKey.detector()==CalorDetInfo::PEM) {
	      
	      uint2 data = tower->pmtData(cellKey,ipmt);
	      if(data > 0){
		// add the channel id and data to bank 
		uint2 chid = pemd->make_chid(calKey); 
		if (pemChanCount<pemd->get_max_size()/2) {
		pemd->add_datum( chid, data ); 
		pemChanCount++;
		}
		else {
		  ERRLOG(ELsevere,"CalDataMaker: ") << "makeDBanks : more PEMD data than maximum number of channels allowed" << endmsg;
		}
	      }
	    }     // endif PEM

	    // PPR
	    
	    if (calKey.detector()==5) {
	      
	      uint2 data = tower->pmtData(cellKey,ipmt);
	      if(data > 0){
		// add the channel id and data to bank 
		uint2 chid = pprd->make_chid(calKey); 
		if (pprChanCount<pprd->get_max_size()/2) {
		pprd->add_datum( chid, data ); 
		pprChanCount++;
		}
		else {
		  ERRLOG(ELsevere,"CalDataMaker: ") << "makeDBanks : more PPRD data than maximum number of channels allowed" << endmsg;

		}
	      }
	    }     // endif PPR
	    
	    // CHA
  
	    if (calKey.detector()==CalorDetInfo::CHA){
	      uint2 data = tower->pmtData(cellKey,ipmt);
	      if(data > 0){
		// add the channel id and data to bank 
		uint2 chid = chad->make_chid(calKey); 
		if (chaChanCount<chad->get_max_size()/2) {
		  chad->add_datum( chid, data ); 
		  chaChanCount++;
		}
		else {
		  ERRLOG(ELsevere,"CalDataMaker: ") << "makeDBanks : more PPRD data than maximum number of channels allowed" << endmsg;
		}
	      }	  
	    }
	    
	    // WHA

	    if (calKey.detector()==CalorDetInfo::WHA){
	    
	      uint2 data= tower->pmtData(cellKey,ipmt);
	      if(data > 0){    
		// add the channel id and data to bank 
		uint2 chid = whad->make_chid(calKey); 
		if (whaChanCount<whad->get_max_size()/2) {
		  whad->add_datum( chid, data ); 
		  whaChanCount++;
		}
		else {
		  ERRLOG(ELsevere,"CalDataMaker: ") << "makeDBanks : more WHAD data than maximum number of channels allowed" << endmsg;
		}
	      }
	    }                           // endif WHA
	  
	    // PHA

	    if (calKey.detector()==CalorDetInfo::PHA){
	      uint2 data = tower->pmtData(cellKey,ipmt);
	      if(data > 0){
		// add the channel id and data to bank 
		uint2 chid = phad->make_chid(calKey); 
		if (phaChanCount<phad->get_max_size()/2) {
		  phad->add_datum( chid, data ); 
		  phaChanCount++;
		}
		else {
		  ERRLOG(ELsevere,"CalDataMaker: ") << "makeDBanks : more PHAD data than maximum number of channels allowed" << endmsg;
		}
	      }
	    }                           // endif PHA
	  }
	}           // end loop over pmts                        
      }        // end if MPA condition
        if(*calIter == MPA){ 
	    // MPA
	CellKey cellKey(tower->iEta(),tower->iPhi(),*calIter);
	// loop over photo-tubes
	int theMaxPmt=1;
	for (int ipmt=0; ipmt<theMaxPmt; ++ipmt){
	  CalorKey calKey;
	  calKey = cellKey.makeCalorKey(cellKey,ipmt);
	    if (calKey.detector()==CalorDetInfo::MPA){
	      uint2 data = tower->pmtData(cellKey,ipmt);
	      if(data > 0){
		// add the channel id and data to bank 
		uint2 chid = mpad->make_chid(calKey); 
		if (mpaChanCount<mpad->get_max_size()/2) {
		  mpad->add_datum( chid, data ); 
		  mpaChanCount++;
		}
		else {
		  ERRLOG(ELsevere,"CalDataMaker: ") << "makeDBanks : more MPAD data than maximum number of channels allowed" << endmsg;
		}
	      }
	    }                           // endif MPA
	}
	}
    }  // end loop over detectors
  }  // end loop over towers

  //
  // transfer data buffers to the banks.

  if (cemChanCount>0) cemd->commit();
  if (pemChanCount>0) pemd->commit();
  if (chaChanCount>0) chad->commit();
  if (whaChanCount>0) whad->commit();
  if (phaChanCount>0) phad->commit();
  if (pprChanCount>0) pprd->commit();
  if (mpaChanCount>0) mpad->commit();

  if (verbose) {
    cemd->print();
    chad->print();
    whad->print();
    pprd->print();
    pemd->print();
    phad->print();
    mpad->print();
  }
  if (cemChanCount>0){ 
    Handle<CEMD_StorableBank> cemd_handle(cemd);
    anEvent->append(cemd_handle);
  } 
  if (pemChanCount>0){ 
    anEvent->append(pemd);
  }  
  if (chaChanCount>0){ 
    anEvent->append(chad);
  }
  if (whaChanCount>0){ 
    anEvent->append(whad);
  }
  if (phaChanCount>0){ 
    anEvent->append(phad);
  }
  if (pprChanCount>0){ 
    anEvent->append(pprd);
  }
  if (mpaChanCount>0){ 
    anEvent->append(mpad);
  }
  return 1;
}



// } // namespace calor