#include <SiPixelLorentzAngle.h>

Inheritance diagram for SiPixelLorentzAngle:

Public Member Functions
virtual void	analyze (const edm::Event &e, const edm::EventSetup &c)
virtual void	beginJob ()
virtual void	endJob ()
float	getLorentzAngle (const uint32_t &) const
const std::map< unsigned int, float > &	getLorentzAngles () const
void	putLorentsAngles (std::map< unsigned int, float > &LA)
bool	putLorentzAngle (const uint32_t &, float &)
	SiPixelLorentzAngle ()
	SiPixelLorentzAngle (const edm::ParameterSet &conf)
virtual	~SiPixelLorentzAngle ()
	~SiPixelLorentzAngle ()
Private Member Functions
void	fillPix (const SiPixelCluster &LocPix, const PixelTopology *topol, Pixinfo &pixinfo)
void	findMean (int i, int i_ring)
Private Attributes
std::map< int, TH2F * >	_h_drift_depth_
std::map< int, TH2F * >	_h_drift_depth_adc2_
std::map< int, TH2F * >	_h_drift_depth_adc_
std::map< int, TH2F * >	_h_drift_depth_noadc_
std::map< int, TH1F * >	_h_mean_
int	bladeF_
double	chi2_
double	chi2F_
Clust	clust_
double	clustChargeMax_
Clust	clustF_
int	clustSizeYMin_
edm::ParameterSet	conf_
int	diskF_
float	eta_
float	etaF_
int	event_
int	event_counter_
int	eventF_
std::string	filename_
std::string	filenameFit_
TH2F *	h_cluster_shape_
TH2F *	h_cluster_shape_adc_
TH2F *	h_cluster_shape_adc_rot_
TH2F *	h_cluster_shape_noadc_
TH2F *	h_cluster_shape_noadc_rot_
TH2F *	h_cluster_shape_rot_
TH1F *	h_drift_depth_adc_slice_
TH1F *	h_tracks_
TFile *	hFile_
int	hist_depth_
int	hist_drift_
int	hist_x_
int	hist_y_
int	hitCounter_
int	isflipped_
int	ladder_
int	layer_
std::map< unsigned int, float >	m_LA
const MagneticField *	magfield
double	max_depth_
double	max_drift_
double	max_x_
double	max_y_
double	min_depth_
double	min_drift_
double	min_x_
double	min_y_
int	module_
int	moduleF_
double	ndof_
double	ndofF_
double	normChi2Max_
int	panelF_
float	phi_
float	phiF_
int	pixelTracksCounter_
Pixinfo	pixinfo_
Pixinfo	pixinfoF_
float	pt_
float	ptF_
double	ptmin_
Rechit	rechit_
Rechit	rechitF_
double	residualMax_
const TransientTrackingRecHitBuilder *	RHBuilder
int	run_
int	runF_
int	sideF_
bool	simData_
Hit	simhit_
Hit	simhitF_
TTree *	SiPixelLorentzAngleTree_
TTree *	SiPixelLorentzAngleTreeForward_
Chi2MeasurementEstimator *	theEstimator
const KFTrajectoryFitter *	theFitter
PropagatorWithMaterial *	thePropagator
PropagatorWithMaterial *	thePropagatorOp
const KFTrajectorySmoother *	theSmoother
KFUpdator *	theUpdator
const TrackerGeometry *	tracker
int	trackEventsCounter_
Hit	trackhit_
Hit	trackhitF_
TrajectoryStateTransform	tsTransform
int	usedHitCounter_
double	width_

Detailed Description

Definition at line 85 of file SiPixelLorentzAngle.h.

Constructor & Destructor Documentation

SiPixelLorentzAngle::SiPixelLorentzAngle ( const edm::ParameterSet & conf ) [explicit]

Definition at line 36 of file SiPixelLorentzAngle.cc.

References hist_x_, hist_y_, max_depth_, max_drift_, max_x_, max_y_, min_depth_, min_drift_, min_x_, min_y_, and width_.

                                                                    : 
  conf_(conf), filename_(conf.getParameter<std::string>("fileName")), filenameFit_(conf.getParameter<std::string>("fileNameFit")), ptmin_(conf.getParameter<double>("ptMin")), simData_(conf.getParameter<bool>("simData")),    normChi2Max_(conf.getParameter<double>("normChi2Max")), clustSizeYMin_(conf.getParameter<int>("clustSizeYMin")), residualMax_(conf.getParameter<double>("residualMax")), clustChargeMax_(conf.getParameter<double>("clustChargeMax")),hist_depth_(conf.getParameter<int>("binsDepth")), hist_drift_(conf.getParameter<int>("binsDrift"))
{
  //    anglefinder_=new  TrackLocalAngle(conf);
  hist_x_ = 50;
  hist_y_ = 100;
  min_x_ = -500.;
  max_x_ = 500.;
  min_y_ = -1500.;
  max_y_ = 500.;
  width_ = 0.0285;
  min_depth_ = -100.;
  max_depth_ = 400.;
  min_drift_ = -1000.; //-200.;(conf.getParameter<double>("residualMax"))
  max_drift_ = 1000.; //400.;

}

SiPixelLorentzAngle::~SiPixelLorentzAngle ( ) [virtual]

Definition at line 55 of file SiPixelLorentzAngle.cc.

{  }

SiPixelLorentzAngle::SiPixelLorentzAngle ( ) [inline]

Definition at line 14 of file SiPixelLorentzAngle.h.

{};

SiPixelLorentzAngle::~SiPixelLorentzAngle ( ) [inline]

Definition at line 15 of file SiPixelLorentzAngle.h.

{};

Member Function Documentation

void SiPixelLorentzAngle::analyze	(	const edm::Event &	e,
		const edm::EventSetup &	c
	)		`[virtual]`

Implements edm::EDAnalyzer.

Definition at line 150 of file SiPixelLorentzAngle.cc.

References _h_drift_depth_adc2_, _h_drift_depth_adc_, _h_drift_depth_noadc_, Pixinfo::adc, Hit::alpha, TrackerHitAssociator::associateHit(), Hit::beta, PXFDetId::blade(), bladeF_, Clust::charge, chi2_, Trajectory::chiSquared(), clust_, clustChargeMax_, SiPixelRecHit::cluster(), clustF_, clustSizeYMin_, Pixinfo::col, conf_, gather_cfg::cout, PXFDetId::disk(), diskF_, shallow::drift(), reco::TrackBase::eta(), eta_, edm::EventID::event(), event_, event_counter_, HcalObjRepresent::Fill(), fillPix(), Hit::gamma, edm::EventSetup::get(), edm::Event::getByLabel(), edm::ParameterSet::getParameter(), h_cluster_shape_adc_, h_cluster_shape_adc_rot_, h_cluster_shape_noadc_, h_cluster_shape_noadc_rot_, h_tracks_, hitCounter_, edm::EventBase::id(), TrackerGeometry::idToDet(), isflipped_, TrajectoryStateOnSurface::isValid(), j, PXBDetId::ladder(), ladder_, PXBDetId::layer(), layer_, TrajectoryStateOnSurface::localDirection(), BaseSiTrackerRecHit2DLocalPos::localPosition(), TrajectoryStateOnSurface::localPosition(), Clust::maxPixelCol, Clust::maxPixelRow, Trajectory::measurements(), Clust::minPixelCol, Clust::minPixelRow, PXFDetId::module(), PXBDetId::module(), module_, moduleF_, Trajectory::ndof(), ndof_, normChi2Max_, Pixinfo::npix, PXFDetId::panel(), panelF_, PV3DBase< T, PVType, FrameType >::perp(), reco::TrackBase::phi(), phi_, PixelSubdetector::PixelBarrel, PixelSubdetector::PixelEndcap, pixelTracksCounter_, pixinfo_, pixinfoF_, reco::TrackBase::pt(), pt_, ptmin_, rechit_, rechitF_, residualMax_, Pixinfo::row, edm::EventID::run(), run_, PXFDetId::side(), sideF_, simData_, simhit_, simhitF_, SiPixelLorentzAngleTree_, SiPixelLorentzAngleTreeForward_, Clust::size_x, Clust::size_y, PixelGeomDetUnit::specificTopology(), GeomDet::surface(), funct::tan(), Surface::toGlobal(), tracker, align::Tracker, trackEventsCounter_, trackhit_, trackhitF_, usedHitCounter_, width_, Rechit::x, PV3DBase< T, PVType, FrameType >::x(), Pixinfo::x, Hit::x, x, Clust::x, Rechit::y, PV3DBase< T, PVType, FrameType >::y(), Pixinfo::y, Hit::y, Clust::y, detailsBasic3DVector::y, and PV3DBase< T, PVType, FrameType >::z().

{  
  event_counter_++;
  //    if(event_counter_ % 500 == 0) cout << "event number " << event_counter_ << endl;
  cout << "event number " << event_counter_ << endl;

  edm::ESHandle<TrackerGeometry> estracker;
  es.get<TrackerDigiGeometryRecord>().get(estracker);
  tracker=&(* estracker);

  TrackerHitAssociator* associate;
  if(simData_) associate = new TrackerHitAssociator(e); else associate = 0; 
  // restet values
  module_=-1;
  layer_=-1;
  ladder_ = -1;
  isflipped_ = -1;
  pt_ = -999;
  eta_ = 999;
  phi_ = 999;
  pixinfo_.npix = 0;

  run_       = e.id().run();
  event_     = e.id().event();
        
  // get the association map between tracks and trajectories
  edm::Handle<TrajTrackAssociationCollection> trajTrackCollectionHandle;
  e.getByLabel(conf_.getParameter<std::string>("src"),trajTrackCollectionHandle);
  if(trajTrackCollectionHandle->size() >0){
    trackEventsCounter_++;
    for(TrajTrackAssociationCollection::const_iterator it = trajTrackCollectionHandle->begin(); it!=trajTrackCollectionHandle->end();++it){
      const Track&      track = *it->val;
      const Trajectory& traj  = *it->key;
                
      // get the trajectory measurements
      std::vector<TrajectoryMeasurement> tmColl = traj.measurements(); 
      //                        TrajectoryStateOnSurface tsos = tsoscomb( itTraj->forwardPredictedState(), itTraj->backwardPredictedState() );
      pt_ = track.pt();
      eta_ = track.eta();
      phi_ = track.phi();
      chi2_ = traj.chiSquared();
      ndof_ = traj.ndof();
      if(pt_ < ptmin_) continue;
      // iterate over trajectory measurements
      std::vector<PSimHit> matched;
      h_tracks_->Fill(0);
      bool pixeltrack = false;
      for(std::vector<TrajectoryMeasurement>::const_iterator itTraj = tmColl.begin(); itTraj != tmColl.end(); itTraj++) {
        if(! itTraj->updatedState().isValid()) continue;
        TransientTrackingRecHit::ConstRecHitPointer recHit = itTraj->recHit();
        if(! recHit->isValid() || recHit->geographicalId().det() != DetId::Tracker ) continue;
        unsigned int subDetID = (recHit->geographicalId().subdetId());
        if( subDetID == PixelSubdetector::PixelBarrel || subDetID == PixelSubdetector::PixelEndcap){
          if(!pixeltrack){
            h_tracks_->Fill(1);
            pixelTracksCounter_++;
          }
          pixeltrack = true;
        }
        if( subDetID == PixelSubdetector::PixelBarrel){
                                
          hitCounter_++;
                                
          DetId detIdObj = recHit->geographicalId();
          const PixelGeomDetUnit * theGeomDet = dynamic_cast<const PixelGeomDetUnit*> ( tracker->idToDet(detIdObj) );
          if(!theGeomDet) continue;

          const PixelTopology * topol = &(theGeomDet->specificTopology());

          if(!topol) continue;
          PXBDetId pxbdetIdObj(detIdObj);
          layer_ = pxbdetIdObj.layer();
          ladder_ = pxbdetIdObj.ladder();
          module_ = pxbdetIdObj.module();
          float tmp1 = theGeomDet->surface().toGlobal(Local3DPoint(0.,0.,0.)).perp();
          float tmp2 = theGeomDet->surface().toGlobal(Local3DPoint(0.,0.,1.)).perp();
          if ( tmp2<tmp1 ) isflipped_ = 1;
          else isflipped_ = 0;
          const SiPixelRecHit * recHitPix = dynamic_cast<const SiPixelRecHit *>((*recHit).hit());
          if(!recHitPix) continue;
          rechit_.x  = recHitPix->localPosition().x();
          rechit_.y  = recHitPix->localPosition().y();
          SiPixelRecHit::ClusterRef const& cluster = recHitPix->cluster();      
        
          // fill entries in clust_
          clust_.x = (cluster)->x();
          clust_.y = (cluster)->y();
          clust_.charge = (cluster->charge())/1000.;
          clust_.size_x = cluster->sizeX();
          clust_.size_y = cluster->sizeY();
          clust_.maxPixelCol = cluster->maxPixelCol();
          clust_.maxPixelRow = cluster->maxPixelRow();
          clust_.minPixelCol = cluster->minPixelCol();
          clust_.minPixelRow = cluster->minPixelRow();
          // fill entries in pixinfo_:
          fillPix(*cluster ,topol, pixinfo_);
          // fill the trackhit info
          TrajectoryStateOnSurface tsos=itTraj->updatedState();
          if(!tsos.isValid()){
            cout << "tsos not valid" << endl;
            continue;   
          }     
          LocalVector trackdirection=tsos.localDirection();
          LocalPoint trackposition=tsos.localPosition();
        
          if(trackdirection.z()==0) continue;                           
          // the local position and direction
          trackhit_.alpha = atan2(trackdirection.z(),trackdirection.x());
          trackhit_.beta = atan2(trackdirection.z(),trackdirection.y());
          trackhit_.gamma = atan2(trackdirection.x(),trackdirection.y());
          trackhit_.x = trackposition.x();
          trackhit_.y = trackposition.y();
                                        
        
          // fill entries in simhit_:   
          if(simData_){
            matched.clear();        
            matched = associate->associateHit((*recHitPix));    
            float dr_start=9999.;
            for (std::vector<PSimHit>::iterator isim = matched.begin(); isim != matched.end(); ++isim){
              DetId simdetIdObj((*isim).detUnitId());
              if (simdetIdObj == detIdObj) {
                float sim_x1 = (*isim).entryPoint().x(); // width (row index, in col direction)
                float sim_y1 = (*isim).entryPoint().y(); // length (col index, in row direction)
                float sim_x2 = (*isim).exitPoint().x();
                float sim_y2 = (*isim).exitPoint().y();
                float sim_xpos = 0.5*(sim_x1+sim_x2);
                float sim_ypos = 0.5*(sim_y1+sim_y2);
                float sim_px = (*isim).momentumAtEntry().x();
                float sim_py = (*isim).momentumAtEntry().y();
                float sim_pz = (*isim).momentumAtEntry().z();
                                                
                float dr = (sim_xpos-(recHitPix->localPosition().x()))*(sim_xpos-recHitPix->localPosition().x()) +
                  (sim_ypos-recHitPix->localPosition().y())*(sim_ypos-recHitPix->localPosition().y());
                if(dr<dr_start) {
                  simhit_.x     = sim_xpos;
                  simhit_.y     = sim_ypos;
                  simhit_.alpha = atan2(sim_pz, sim_px);
                  simhit_.beta  = atan2(sim_pz, sim_py);
                  simhit_.gamma = atan2(sim_px, sim_py);
                  dr_start = dr;
                }
              }
            } // end of filling simhit_
          }
          // is one pixel in cluster a large pixel ? (hit will be excluded)
          bool large_pix = false;
          for (int j = 0; j <  pixinfo_.npix; j++){
            int colpos = static_cast<int>(pixinfo_.col[j]);
            if (pixinfo_.row[j] == 0 || pixinfo_.row[j] == 79 || pixinfo_.row[j] == 80 || pixinfo_.row[j] == 159 || colpos % 52 == 0 || colpos % 52 == 51 ){
              large_pix = true; 
            }
          }
                                        
          double residual = TMath::Sqrt( (trackhit_.x - rechit_.x) * (trackhit_.x - rechit_.x) + (trackhit_.y - rechit_.y) * (trackhit_.y - rechit_.y) );
        
          SiPixelLorentzAngleTree_->Fill();
          if( !large_pix && (chi2_/ndof_) < normChi2Max_ && cluster->sizeY() >= clustSizeYMin_ && residual < residualMax_ && (cluster->charge() < clustChargeMax_)){
            usedHitCounter_++;
            // iterate over pixels in hit
            for (int j = 0; j <  pixinfo_.npix; j++){
              // use trackhits
              float dx = (pixinfo_.x[j]  - (trackhit_.x - width_/2. / TMath::Tan(trackhit_.alpha))) * 10000.;
              float dy = (pixinfo_.y[j]  - (trackhit_.y - width_/2. / TMath::Tan(trackhit_.beta))) * 10000.;
              float depth = dy * tan(trackhit_.beta);
              float drift = dx - dy * tan(trackhit_.gamma);
              _h_drift_depth_adc_[module_ + (layer_ -1) * 8]->Fill(drift, depth, pixinfo_.adc[j]);
              _h_drift_depth_adc2_[module_ + (layer_ -1) * 8]->Fill(drift, depth, pixinfo_.adc[j]*pixinfo_.adc[j]);
              _h_drift_depth_noadc_[module_ + (layer_ -1) * 8]->Fill(drift, depth);             
              if( layer_ == 3 && module_==1 && isflipped_){
                float dx_rot = dx * TMath::Cos(trackhit_.gamma) + dy * TMath::Sin(trackhit_.gamma);
                float dy_rot = dy * TMath::Cos(trackhit_.gamma) - dx * TMath::Sin(trackhit_.gamma) ;
                h_cluster_shape_adc_->Fill(dx, dy, pixinfo_.adc[j]);
                h_cluster_shape_noadc_->Fill(dx, dy);
                h_cluster_shape_adc_rot_->Fill(dx_rot, dy_rot, pixinfo_.adc[j]);
                h_cluster_shape_noadc_rot_->Fill(dx_rot, dy_rot);
              }
            }
          }
        } else if (subDetID == PixelSubdetector::PixelEndcap) {
          DetId detIdObj = recHit->geographicalId();
          const PixelGeomDetUnit * theGeomDet = dynamic_cast<const PixelGeomDetUnit*> ( tracker->idToDet(detIdObj) );
          if(!theGeomDet) continue;
          
          const PixelTopology * topol = &(theGeomDet->specificTopology());

          if(!topol) continue;
          PXFDetId pxfdetIdObj(detIdObj);
          sideF_ = pxfdetIdObj.side();
          diskF_ = pxfdetIdObj.disk();
          bladeF_ = pxfdetIdObj.blade();
          panelF_ = pxfdetIdObj.panel();
          moduleF_ = pxfdetIdObj.module();
          //                                    float tmp1 = theGeomDet->surface().toGlobal(Local3DPoint(0.,0.,0.)).perp();
          //                                    float tmp2 = theGeomDet->surface().toGlobal(Local3DPoint(0.,0.,1.)).perp();
          //                                    if ( tmp2<tmp1 ) isflipped_ = 1;
          //                                    else isflipped_ = 0;
          const SiPixelRecHit * recHitPix = dynamic_cast<const SiPixelRecHit *>((*recHit).hit());
          if(!recHitPix) continue;
          rechitF_.x  = recHitPix->localPosition().x();
          rechitF_.y  = recHitPix->localPosition().y();
          SiPixelRecHit::ClusterRef const& cluster = recHitPix->cluster();      
        
          // fill entries in clust_
          clustF_.x = (cluster)->x();
          clustF_.y = (cluster)->y();
          clustF_.charge = (cluster->charge())/1000.;
          clustF_.size_x = cluster->sizeX();
          clustF_.size_y = cluster->sizeY();
          clustF_.maxPixelCol = cluster->maxPixelCol();
          clustF_.maxPixelRow = cluster->maxPixelRow();
          clustF_.minPixelCol = cluster->minPixelCol();
          clustF_.minPixelRow = cluster->minPixelRow();
          // fill entries in pixinfo_:
          fillPix(*cluster ,topol, pixinfoF_);
          // fill the trackhit info
          TrajectoryStateOnSurface tsos=itTraj->updatedState();
          if(!tsos.isValid()){
            cout << "tsos not valid" << endl;
            continue;   
          }     
          LocalVector trackdirection=tsos.localDirection();
          LocalPoint trackposition=tsos.localPosition();
        
          if(trackdirection.z()==0) continue;                           
          // the local position and direction
          trackhitF_.alpha = atan2(trackdirection.z(),trackdirection.x());
          trackhitF_.beta = atan2(trackdirection.z(),trackdirection.y());
          trackhitF_.gamma = atan2(trackdirection.x(),trackdirection.y());
          trackhitF_.x = trackposition.x();
          trackhitF_.y = trackposition.y();
                                        
        
          // fill entries in simhit_:   
          if(simData_){
            matched.clear();        
            matched = associate->associateHit((*recHitPix));    
            float dr_start=9999.;
            for (std::vector<PSimHit>::iterator isim = matched.begin(); isim != matched.end(); ++isim){
              DetId simdetIdObj((*isim).detUnitId());
              if (simdetIdObj == detIdObj) {
                float sim_x1 = (*isim).entryPoint().x(); // width (row index, in col direction)
                float sim_y1 = (*isim).entryPoint().y(); // length (col index, in row direction)
                float sim_x2 = (*isim).exitPoint().x();
                float sim_y2 = (*isim).exitPoint().y();
                float sim_xpos = 0.5*(sim_x1+sim_x2);
                float sim_ypos = 0.5*(sim_y1+sim_y2);
                float sim_px = (*isim).momentumAtEntry().x();
                float sim_py = (*isim).momentumAtEntry().y();
                float sim_pz = (*isim).momentumAtEntry().z();
                                                
                float dr = (sim_xpos-(recHitPix->localPosition().x()))*(sim_xpos-recHitPix->localPosition().x()) +
                  (sim_ypos-recHitPix->localPosition().y())*(sim_ypos-recHitPix->localPosition().y());
                if(dr<dr_start) {
                  simhitF_.x     = sim_xpos;
                  simhitF_.y     = sim_ypos;
                  simhitF_.alpha = atan2(sim_pz, sim_px);
                  simhitF_.beta  = atan2(sim_pz, sim_py);
                  simhitF_.gamma = atan2(sim_px, sim_py);
                  dr_start = dr;
                }
              }
            } // end of filling simhit_
          }
          SiPixelLorentzAngleTreeForward_->Fill();
        }
      } //end iteration over trajectory measurements
    } //end iteration over trajectories
  }

        
        
        
        
}

void SiPixelLorentzAngle::beginJob ( void ) [virtual]

Reimplemented from edm::EDAnalyzer.

Definition at line 57 of file SiPixelLorentzAngle.cc.

References _h_drift_depth_, _h_drift_depth_adc2_, _h_drift_depth_adc_, _h_drift_depth_noadc_, _h_mean_, Pixinfo::adc, bladeF_, chi2_, clust_, clustF_, Pixinfo::col, diskF_, eta_, event_, event_counter_, filename_, h_cluster_shape_, h_cluster_shape_adc_, h_cluster_shape_adc_rot_, h_cluster_shape_noadc_, h_cluster_shape_noadc_rot_, h_cluster_shape_rot_, h_tracks_, hFile_, hist_depth_, hist_drift_, hist_x_, hist_y_, hitCounter_, isflipped_, ladder_, layer_, max_depth_, max_drift_, max_x_, max_y_, min_depth_, min_drift_, min_x_, min_y_, module_, moduleF_, mergeVDriftHistosByStation::name, ndof_, Pixinfo::npix, panelF_, phi_, pixelTracksCounter_, pixinfo_, pixinfoF_, pt_, rechit_, rechitF_, Pixinfo::row, run_, sideF_, simhit_, simhitF_, SiPixelLorentzAngleTree_, SiPixelLorentzAngleTreeForward_, trackEventsCounter_, trackhit_, trackhitF_, usedHitCounter_, Pixinfo::x, and Pixinfo::y.

{

  //    cout << "started SiPixelLorentzAngle" << endl;
  hFile_ = new TFile (filename_.c_str(), "RECREATE" );
  int bufsize = 64000;
  // create tree structure
  SiPixelLorentzAngleTree_ = new TTree("SiPixelLorentzAngleTree_","SiPixel LorentzAngle tree", bufsize);
  SiPixelLorentzAngleTree_->Branch("run", &run_, "run/I", bufsize);
  SiPixelLorentzAngleTree_->Branch("event", &event_, "event/I", bufsize);
  SiPixelLorentzAngleTree_->Branch("module", &module_, "module/I", bufsize);
  SiPixelLorentzAngleTree_->Branch("ladder", &ladder_, "ladder/I", bufsize);
  SiPixelLorentzAngleTree_->Branch("layer", &layer_, "layer/I", bufsize);
  SiPixelLorentzAngleTree_->Branch("isflipped", &isflipped_, "isflipped/I", bufsize);
  SiPixelLorentzAngleTree_->Branch("pt", &pt_, "pt/F", bufsize);
  SiPixelLorentzAngleTree_->Branch("eta", &eta_, "eta/F", bufsize);
  SiPixelLorentzAngleTree_->Branch("phi", &phi_, "phi/F", bufsize);
  SiPixelLorentzAngleTree_->Branch("chi2", &chi2_, "chi2/D", bufsize);
  SiPixelLorentzAngleTree_->Branch("ndof", &ndof_, "ndof/D", bufsize);
  SiPixelLorentzAngleTree_->Branch("trackhit", &trackhit_, "x/F:y/F:alpha/D:beta/D:gamma_/D", bufsize);
  SiPixelLorentzAngleTree_->Branch("simhit", &simhit_, "x/F:y/F:alpha/D:beta/D:gamma_/D", bufsize);
  SiPixelLorentzAngleTree_->Branch("npix", &pixinfo_.npix, "npix/I", bufsize);
  SiPixelLorentzAngleTree_->Branch("rowpix", pixinfo_.row, "row[npix]/F", bufsize);
  SiPixelLorentzAngleTree_->Branch("colpix", pixinfo_.col, "col[npix]/F", bufsize);
  SiPixelLorentzAngleTree_->Branch("adc", pixinfo_.adc, "adc[npix]/F", bufsize);
  SiPixelLorentzAngleTree_->Branch("xpix", pixinfo_.x, "x[npix]/F", bufsize);
  SiPixelLorentzAngleTree_->Branch("ypix", pixinfo_.y, "y[npix]/F", bufsize);
  SiPixelLorentzAngleTree_->Branch("clust", &clust_, "x/F:y/F:charge/F:size_x/I:size_y/I:maxPixelCol/I:maxPixelRow:minPixelCol/I:minPixelRow/I", bufsize);
  SiPixelLorentzAngleTree_->Branch("rechit", &rechit_, "x/F:y/F", bufsize);
        
  SiPixelLorentzAngleTreeForward_ = new TTree("SiPixelLorentzAngleTreeForward_","SiPixel LorentzAngle tree forward", bufsize);
  SiPixelLorentzAngleTreeForward_->Branch("run", &run_, "run/I", bufsize);
  SiPixelLorentzAngleTreeForward_->Branch("event", &event_, "event/I", bufsize);
  SiPixelLorentzAngleTreeForward_->Branch("side", &sideF_, "side/I", bufsize);
  SiPixelLorentzAngleTreeForward_->Branch("disk", &diskF_, "disk/I", bufsize);
  SiPixelLorentzAngleTreeForward_->Branch("blade", &bladeF_, "blade/I", bufsize);
  SiPixelLorentzAngleTreeForward_->Branch("panel", &panelF_, "panel/I", bufsize);
  SiPixelLorentzAngleTreeForward_->Branch("module", &moduleF_, "module/I", bufsize);
  SiPixelLorentzAngleTreeForward_->Branch("pt", &pt_, "pt/F", bufsize);
  SiPixelLorentzAngleTreeForward_->Branch("eta", &eta_, "eta/F", bufsize);
  SiPixelLorentzAngleTreeForward_->Branch("phi", &phi_, "phi/F", bufsize);
  SiPixelLorentzAngleTreeForward_->Branch("chi2", &chi2_, "chi2/D", bufsize);
  SiPixelLorentzAngleTreeForward_->Branch("ndof", &ndof_, "ndof/D", bufsize);
  SiPixelLorentzAngleTreeForward_->Branch("trackhit", &trackhitF_, "x/F:y/F:alpha/D:beta/D:gamma_/D", bufsize);
  SiPixelLorentzAngleTreeForward_->Branch("simhit", &simhitF_, "x/F:y/F:alpha/D:beta/D:gamma_/D", bufsize);
  SiPixelLorentzAngleTreeForward_->Branch("npix", &pixinfoF_.npix, "npix/I", bufsize);
  SiPixelLorentzAngleTreeForward_->Branch("rowpix", pixinfoF_.row, "row[npix]/F", bufsize);
  SiPixelLorentzAngleTreeForward_->Branch("colpix", pixinfoF_.col, "col[npix]/F", bufsize);
  SiPixelLorentzAngleTreeForward_->Branch("adc", pixinfoF_.adc, "adc[npix]/F", bufsize);
  SiPixelLorentzAngleTreeForward_->Branch("xpix", pixinfoF_.x, "x[npix]/F", bufsize);
  SiPixelLorentzAngleTreeForward_->Branch("ypix", pixinfoF_.y, "y[npix]/F", bufsize);
  SiPixelLorentzAngleTreeForward_->Branch("clust", &clustF_, "x/F:y/F:charge/F:size_x/I:size_y/I:maxPixelCol/I:maxPixelRow:minPixelCol/I:minPixelRow/I", bufsize);
  SiPixelLorentzAngleTreeForward_->Branch("rechit", &rechitF_, "x/F:y/F", bufsize);
        
        
  //book histograms
  char name[128];
  for(int i_module = 1; i_module<=8; i_module++){
    for(int i_layer = 1; i_layer<=3; i_layer++){
      sprintf(name, "h_drift_depth_adc_layer%i_module%i", i_layer, i_module); 
      _h_drift_depth_adc_[i_module + (i_layer -1) * 8] = new TH2F(name,name,hist_drift_ , min_drift_, max_drift_, hist_depth_, min_depth_, max_depth_);
      sprintf(name, "h_drift_depth_adc2_layer%i_module%i", i_layer, i_module); 
      _h_drift_depth_adc2_[i_module + (i_layer -1) * 8] = new TH2F(name,name,hist_drift_ , min_drift_, max_drift_, hist_depth_, min_depth_, max_depth_);
      sprintf(name, "h_drift_depth_noadc_layer%i_module%i", i_layer, i_module); 
      _h_drift_depth_noadc_[i_module + (i_layer -1) * 8] = new TH2F(name,name,hist_drift_ , min_drift_, max_drift_, hist_depth_, min_depth_, max_depth_);
      sprintf(name, "h_drift_depth_layer%i_module%i", i_layer, i_module); 
      _h_drift_depth_[i_module + (i_layer -1) * 8] = new TH2F(name,name,hist_drift_ , min_drift_, max_drift_, hist_depth_, min_depth_, max_depth_);
      sprintf(name, "h_mean_layer%i_module%i", i_layer, i_module); 
      _h_mean_[i_module + (i_layer -1) * 8] = new TH1F(name,name,hist_depth_, min_depth_, max_depth_);
    }
  }
        
  // just for some expaining plots
  h_cluster_shape_adc_  = new TH2F("h_cluster_shape_adc","cluster shape with adc weight", hist_x_, min_x_, max_x_, hist_y_, min_y_, max_y_);
  h_cluster_shape_noadc_  = new TH2F("h_cluster_shape_noadc","cluster shape without adc weight", hist_x_, min_x_, max_x_, hist_y_, min_y_, max_y_);
  h_cluster_shape_  = new TH2F("h_cluster_shape","cluster shape", hist_x_, min_x_, max_x_, hist_y_, min_y_, max_y_);
  h_cluster_shape_adc_rot_  = new TH2F("h_cluster_shape_adc_rot","cluster shape with adc weight", hist_x_, min_x_, max_x_, hist_y_, -max_y_, -min_y_);
  h_cluster_shape_noadc_rot_  = new TH2F("h_cluster_shape_noadc_rot","cluster shape without adc weight", hist_x_, min_x_, max_x_, hist_y_, -max_y_, -min_y_);
  h_cluster_shape_rot_  = new TH2F("h_cluster_shape_rot","cluster shape", hist_x_, min_x_, max_x_, hist_y_, -max_y_, -min_y_);
  h_tracks_ = new TH1F("h_tracks","h_tracks",2,0.,2.);
  event_counter_ = 0;
  trackEventsCounter_ = 0;
  //    trackcounter_ = 0;
  hitCounter_ = 0;
  usedHitCounter_ = 0;
  pixelTracksCounter_ = 0;
//   edm::ESHandle<TrackerGeometry> estracker; //this block should not be in beginJob()
//   c.get<TrackerDigiGeometryRecord>().get(estracker);
//   tracker=&(* estracker);
}

void SiPixelLorentzAngle::endJob ( void ) [virtual]

Reimplemented from edm::EDAnalyzer.

Definition at line 426 of file SiPixelLorentzAngle.cc.

References _h_drift_depth_, _h_drift_depth_adc2_, _h_drift_depth_adc_, _h_drift_depth_noadc_, _h_mean_, gather_cfg::cout, event_counter_, validate-o2o-wbm::f1, filenameFit_, findMean(), FitTarget::Fit, h_cluster_shape_adc_, h_cluster_shape_adc_rot_, h_cluster_shape_noadc_, h_cluster_shape_noadc_rot_, h_drift_depth_adc_slice_, h_tracks_, hFile_, hist_depth_, hist_drift_, hitCounter_, i, max_drift_, min_drift_, p1, pixelTracksCounter_, trackEventsCounter_, and usedHitCounter_.

{
  // produce histograms with the average adc counts
  for(int i_ring = 1; i_ring<=24; i_ring++){
    _h_drift_depth_[i_ring]->Divide(_h_drift_depth_adc_[i_ring], _h_drift_depth_noadc_[i_ring]);
  }
        
  h_drift_depth_adc_slice_ = new TH1F("h_drift_depth_adc_slice", "slice of adc histogram", hist_drift_ , min_drift_, max_drift_);

  TF1 *f1 = new TF1("f1","[0] + [1]*x",50., 235.); 
  f1->SetParName(0,"p0");
  f1->SetParName(1,"p1");
  f1->SetParameter(0,0);
  f1->SetParameter(1,0.4);
  ofstream fLorentzFit( filenameFit_.c_str(), ios::trunc );
  cout.precision( 4 );
  fLorentzFit << "module" << "\t" << "layer" << "\t" << "offset" << "\t" << "error" << "\t" << "slope" << "\t" << "error" << "\t" "rel.err" << "\t" "pull" << "\t" << "chi2" << "\t" << "prob" << endl;
  //loop over modlues and layers to fit the lorentz angle
  for( int i_layer = 1; i_layer<=3; i_layer++){
    for(int i_module = 1; i_module<=8; i_module++){
      //loop over bins in depth (z-local-coordinate) (in order to fit slices)
      for( int i = 1; i <= hist_depth_; i++){
        findMean(i, (i_module + (i_layer - 1) * 8));    
      }// end loop over bins in depth 
      _h_mean_[i_module + (i_layer - 1) * 8]->Fit(f1,"ERQ");
      double p0 = f1->GetParameter(0);
      double e0 = f1->GetParError(0);
      double p1 = f1->GetParameter(1);
      double e1 = f1->GetParError(1);
      double chi2 = f1->GetChisquare();
      double prob = f1->GetProb();
      fLorentzFit << setprecision( 4 ) << i_module << "\t" << i_layer << "\t" << p0 << "\t" << e0 << "\t" << p1 << setprecision( 3 ) << "\t" << e1 << "\t" << e1 / p1 *100. << "\t"<< (p1 - 0.424) / e1 << "\t"<< chi2 << "\t" << prob << endl;
    }
  } // end loop over modules and layers
  fLorentzFit.close(); 
  hFile_->cd();
  for(int i_module = 1; i_module<=8; i_module++){
    for(int i_layer = 1; i_layer<=3; i_layer++){
      _h_drift_depth_adc_[i_module + (i_layer -1) * 8]->Write();
      _h_drift_depth_adc2_[i_module + (i_layer -1) * 8]->Write();
      _h_drift_depth_noadc_[i_module + (i_layer -1) * 8]->Write();
      _h_drift_depth_[i_module + (i_layer -1) * 8]->Write();
      _h_mean_[i_module + (i_layer -1) * 8]->Write();
    }
  }
  h_cluster_shape_adc_->Write();
  h_cluster_shape_noadc_->Write();
  h_cluster_shape_adc_rot_->Write();
  h_cluster_shape_noadc_rot_->Write();
  h_tracks_->Write();
        
  hFile_->Write();
  hFile_->Close();
  cout << "events: " << event_counter_ << endl;
  cout << "events with tracks: " << trackEventsCounter_ << endl;        
  cout << "events with pixeltracks: " << pixelTracksCounter_ << endl;   
  cout << "hits in the pixel: " << hitCounter_ << endl;
  cout << "number of used Hits: " << usedHitCounter_ << endl;
}

void SiPixelLorentzAngle::fillPix	(	const SiPixelCluster &	LocPix,
		const PixelTopology *	topol,
		Pixinfo &	pixinfo
	)		`[inline, private]`

Definition at line 486 of file SiPixelLorentzAngle.cc.

References Pixinfo::adc, Pixinfo::col, Topology::localPosition(), Pixinfo::npix, SiPixelCluster::pixels(), Pixinfo::row, PV3DBase< T, PVType, FrameType >::x(), Pixinfo::x, PV3DBase< T, PVType, FrameType >::y(), and Pixinfo::y.

Referenced by analyze().

{
  const std::vector<SiPixelCluster::Pixel>& pixvector = LocPix.pixels();
  pixinfo.npix = 0;
  for(std::vector<SiPixelCluster::Pixel>::const_iterator itPix = pixvector.begin(); itPix != pixvector.end(); itPix++){
    //  for(pixinfo.npix = 0; pixinfo.npix < static_cast<int>(pixvector.size()); ++pixinfo.npix) {
    pixinfo.row[pixinfo.npix] = itPix->x;
    pixinfo.col[pixinfo.npix] = itPix->y;
    pixinfo.adc[pixinfo.npix] = itPix->adc;
    LocalPoint lp = topol->localPosition(MeasurementPoint(itPix->x + 0.5, itPix->y+0.5));
    pixinfo.x[pixinfo.npix] = lp.x();
    pixinfo.y[pixinfo.npix]= lp.y();
    pixinfo.npix++;
  }
}

void SiPixelLorentzAngle::findMean	(	int	i,
		int	i_ring
	)		`[private]`

Definition at line 503 of file SiPixelLorentzAngle.cc.

References _h_drift_depth_, _h_drift_depth_adc2_, _h_drift_depth_adc_, _h_drift_depth_noadc_, _h_mean_, error, h_drift_depth_adc_slice_, hist_drift_, j, timingPdfMaker::mean, and mathSSE::sqrt().

Referenced by endJob().

{
  double nentries = 0;
        
  h_drift_depth_adc_slice_->Reset("ICE");
        
  // determine sigma and sigma^2 of the adc counts and average adc counts
  //loop over bins in drift width
  for( int j = 1; j<= hist_drift_; j++){
    if(_h_drift_depth_noadc_[i_ring]->GetBinContent(j, i) >= 1){
      double adc_error2 = (_h_drift_depth_adc2_[i_ring]->GetBinContent(j,i) - _h_drift_depth_adc_[i_ring]->GetBinContent(j,i)*_h_drift_depth_adc_[i_ring]->GetBinContent(j, i) / _h_drift_depth_noadc_[i_ring]->GetBinContent(j, i)) /  _h_drift_depth_noadc_[i_ring]->GetBinContent(j, i);
      _h_drift_depth_adc_[i_ring]->SetBinError(j, i, sqrt(adc_error2));
      double error2 = adc_error2 / (_h_drift_depth_noadc_[i_ring]->GetBinContent(j,i) - 1.);
      _h_drift_depth_[i_ring]->SetBinError(j,i,sqrt(error2));
    } 
    else{
      _h_drift_depth_[i_ring]->SetBinError(j,i,0);
      _h_drift_depth_adc_[i_ring]->SetBinError(j, i, 0);
    }
    h_drift_depth_adc_slice_->SetBinContent(j, _h_drift_depth_adc_[i_ring]->GetBinContent(j,i));
    h_drift_depth_adc_slice_->SetBinError(j, _h_drift_depth_adc_[i_ring]->GetBinError(j,i));
    nentries += _h_drift_depth_noadc_[i_ring]->GetBinContent(j,i);      
  } // end loop over bins in drift width
                
  double mean = h_drift_depth_adc_slice_->GetMean(1); 
  double error = 0;
  if(nentries != 0){
    error = h_drift_depth_adc_slice_->GetRMS(1) / sqrt(nentries);
  }
                
  _h_mean_[i_ring]->SetBinContent(i, mean);
  _h_mean_[i_ring]->SetBinError(i, error);

}

float SiPixelLorentzAngle::getLorentzAngle ( const uint32_t & detid ) const

Definition at line 13 of file SiPixelLorentzAngle.cc.

References m_LA.

Referenced by PixelCPEBase::computeLorentzShifts(), and PixelCPEBase::driftDirection().

                                                                       {
  std::map<unsigned int,float>::const_iterator id=m_LA.find(detid);
  if(id!=m_LA.end()) return id->second;
  else {
    edm::LogError("SiPixelLorentzAngle") << "SiPixelLorentzAngle for DetID " << detid << " is not stored" << std::endl; 
  }
  return 0;
}

const std::map<unsigned int,float>& SiPixelLorentzAngle::getLorentzAngles ( ) const [inline]

Definition at line 18 of file SiPixelLorentzAngle.h.

References m_LA.

{return m_LA;}

void SiPixelLorentzAngle::putLorentsAngles ( std::map< unsigned int, float > & LA ) [inline]

Definition at line 17 of file SiPixelLorentzAngle.h.

References m_LA.

{m_LA=LA;}

bool SiPixelLorentzAngle::putLorentzAngle	(	const uint32_t &	detid,
		float &	value
	)

Definition at line 4 of file SiPixelLorentzAngle.cc.

References cond::rpcobgas::detid, m_LA, and relativeConstraints::value.

Referenced by SiPixelLorentzAngleDB::analyze().

                                                                            {
  std::map<unsigned int,float>::const_iterator id=m_LA.find(detid);
  if(id!=m_LA.end()){
    edm::LogError("SiPixelLorentzAngle") << "SiPixelLorentzAngle for DetID " << detid << " is already stored. Skippig this put" << std::endl;
    return false;
  }
  else m_LA[detid]=value;
  return true;
}