00001 #include "StEEmcClusterMaker.h"
00002
00003 #include <algorithm>
00004 #include <iostream>
00005
00006 #include "TMath.h"
00007
00008
00009
00010 #include "StEvent/StEvent.h"
00011 #include "StEvent/StEmcCollection.h"
00012 #include "StEvent/StEmcDetector.h"
00013 #include "StEvent/StEmcModule.h"
00014 #include "StEvent/StEmcClusterCollection.h"
00015 #include "StEvent/StEmcCluster.h"
00016
00017
00018 #define FIXED_SMD_SEED
00019
00020
00021
00022
00023 #define DEBUG 0
00024
00025 ClassImp(StEEmcClusterMaker);
00026
00027
00028 StEEmcClusterMaker::StEEmcClusterMaker(const Char_t *name):StMaker(name)
00029 {
00030
00032 mFillStEvent = 0;
00033
00035 mSuppress = false;
00036
00039 mClusterId = 0;
00040
00042 const Float_t eseeds[] = { 0.6, 1.0/1000, 1.0/1000, 1.0/1000., 0.3/1000., 0.3/1000. };
00043 for ( Int_t i = 0; i < 6; i++ ) seedEnergy( eseeds[i], i );
00044
00045 mMaxExtent = 6;
00046 mSeedFloor = 2.0;
00047
00060
00061 StEEmcClusterVec_t t;
00062 std::vector< StEEmcClusterVec_t > layers;
00063 for ( Int_t i = 0; i < 4; i++ ) layers.push_back(t);
00064 for ( Int_t i = 0; i < 12; i++ ) mTowerClusters.push_back(layers);
00065
00066 StEEmcSmdClusterVec_t s;
00067 std::vector< StEEmcSmdClusterVec_t > planes;
00068 planes.push_back(s);
00069 planes.push_back(s);
00070 for ( Int_t i = 0; i < 12; i++ ) mSmdClusters.push_back(planes);
00071
00072 StEEmcTowerVec_t tow;
00073 std::vector< StEEmcTowerVec_t > lay;
00074 for ( Int_t i = 0; i < 4; i++ ) lay.push_back(tow);
00075 for ( Int_t i = 0; i < 12; i++ ) mSeedTowers.push_back(lay);
00076
00077 mEEtow=new EEmcGeomSimple();
00078 mEEsmd=EEmcSmdGeom::instance();
00079 mEEmap=EEmcSmdMap::instance();
00080
00081 }
00082
00083
00084 Int_t StEEmcClusterMaker::Init()
00085 {
00086 mEEanalysis=(StEEmcA2EMaker*)GetMaker(mAnalysisName);
00087 assert(mEEanalysis);
00088
00089 return StMaker::Init();
00090 }
00091
00092
00093 Int_t StEEmcClusterMaker::Make()
00094 {
00095
00097 if ( !buildTowerClusters() ) return kStWarn;
00098
00100 if ( !buildSmdClusters() ) return kStWarn;
00101
00103 if ( mFillStEvent ) fillStEvent();
00104
00106 if ( mFillStEvent )
00107 if ( !verifyStEvent() ) Warning("Make","StEvent not properly copied");
00108
00109 return kStOK;
00110 }
00111
00112
00113 void StEEmcClusterMaker::Clear( Option_t *opts )
00114 {
00115
00117 for ( Int_t sector=0; sector<12; sector++ ) {
00118 for ( Int_t layer=0; layer<4; layer++ )
00119 mTowerClusters[sector][layer].clear();
00120 for ( Int_t plane=0; plane<2; plane++ )
00121 mSmdClusters[sector][plane].clear();
00122 }
00123
00124 for ( Int_t i=0;i<6;i++ ) mNumberOfClusters[i]=0;
00125
00126 mClusterId = 0;
00127
00128 return;
00129 }
00130
00131
00132
00133 Bool_t StEEmcClusterMaker::buildTowerClusters()
00134 {
00135
00138
00139
00140
00142 for ( Int_t layer=0;layer<4;layer++)
00143 {
00144
00147 Float_t weights[720]; for ( Int_t i=0;i<720;i++ ) weights[i]=0.;
00148
00150 StEEmcClusterVec_t myClusters;
00151
00152
00154 StEEmcTowerVec_t towers = mEEanalysis -> towers( layer );
00155
00156
00157
00159 std::sort(towers.begin(),towers.end());
00161 std::reverse(towers.begin(),towers.end());
00162
00163
00166 StEEmcTowerVec_t::iterator last = towers.begin();
00167 while ( last != towers.end() ) {
00168
00169 if ( (*last).energy() < mSeedEnergy[layer] ) break;
00174
00175 #if DEBUG
00176 std::cout << "-- Seed tower ----------------------" << std::endl;
00177 (*last).print();
00178 #endif
00179 last++;
00180 }
00181
00182 StEEmcTowerVec_t::iterator iter = towers.begin();
00183 while ( iter != towers.end() ) {
00184
00186 if ( iter==last ) break;
00187
00188
00195 if ( weights[ (*iter).index() ] > 0. ) {
00196 iter++;
00197 continue;
00198 }
00199
00202 for ( Int_t in=0; in < (*iter).numberOfNeighbors(); in++ ) {
00203 StEEmcTower t=(*iter).neighbor(in);
00204 weights[ t.index() ] += (*iter).energy();
00205 }
00206
00207
00208 iter++;
00209
00210 }
00211
00212
00214 iter=towers.begin();
00215 while ( iter != towers.end() ) {
00216
00218 if ( iter==last ) break;
00219
00223 if ( weights[ (*iter).index() ] > 0. ) {
00224 iter++;
00225 continue;
00226 }
00227
00228 StEEmcCluster cluster;
00229 StEEmcTower seed=(*iter);
00230 #if DEBUG
00231 std::cout << "--- Clustering ----------------" << std::endl;
00232 seed.print();
00233 #endif
00234
00235 TVector3 momentum;
00236 cluster.add(seed,1.0);
00237 UInt_t sec,sub,eta;
00238 sec=(UInt_t)seed.sector();
00239 sub=(UInt_t)seed.subsector();
00240 eta=(UInt_t)seed.etabin();
00241 TVector3 d=mEEtow->getTowerCenter(sec,sub,eta).Unit();
00242 momentum += ( seed.energy() * d );
00243
00244 for ( Int_t in=0; in<seed.numberOfNeighbors(); in++ ) {
00245 StEEmcTower t=seed.neighbor(in);
00246 sec=(UInt_t)t.sector();
00247 sub=(UInt_t)t.subsector();
00248 eta=(UInt_t)t.etabin();
00249 d=mEEtow->getTowerCenter(sec,sub,eta).Unit();
00250 Float_t weight = seed.energy() / weights[ t.index() ];
00251 momentum += ( t.energy() * d * weight );
00252 cluster.add(t, weight);
00253 #if DEBUG
00254 std::cout << "adding " << t.name() << " E=" << t.energy() << " W=" << weight << std::endl;
00255 #endif
00256 }
00257
00258
00260 cluster.momentum( momentum );
00261 cluster.key( mClusterId++ );
00262 #if DEBUG
00263 cluster.print();
00264 #endif
00265
00266 mTowerClusters[ seed.sector() ][ layer ].push_back( cluster );
00267 mNumberOfClusters[layer]++;
00268
00269 iter++;
00270
00271 }
00272
00273 }
00274
00275
00276 return true;
00277
00278 }
00279
00280
00281 #ifdef FIXED_SMD_SEED
00282 Bool_t StEEmcClusterMaker::buildSmdClusters()
00283 {
00284
00285 Int_t max_extent = mMaxExtent;
00286
00287
00289 for ( Int_t sector=0; sector<12; sector++ )
00291 for ( Int_t plane=0; plane<2; plane++ ) {
00292
00294 Float_t floor[288]; for ( Int_t i=0; i<288; i++ ) floor[i]=0.;
00295
00297 Float_t energy[288]; for ( Int_t i=0; i<288; i++ ) energy[i]=0.;
00298
00300 StEEmcStripVec_t strips=mEEanalysis->strips(sector,plane);
00301 StEEmcStripVec_t seeds;
00303 std::sort(strips.begin(),strips.end());
00305 std::reverse(strips.begin(),strips.end());
00306
00307
00309 StEEmcStripVec_t::iterator istrip=strips.begin();
00310 while ( istrip != strips.end() ) {
00311 if ( (*istrip).stat()||(*istrip).fail() ) {
00312 istrip++;
00313 continue;
00314 }
00315
00316 energy[ (*istrip).index() ] = (*istrip).energy();
00317 istrip++;
00318 }
00319
00321 std::vector<Bool_t> seed_flags( strips.size(), false );
00322
00324 Int_t nstrip=0;
00325 Int_t nseeds=0;
00326 istrip=strips.begin();
00327 while ( istrip!=strips.end() ) {
00328
00329 Int_t index=(*istrip).index();
00330 Float_t eseed=(*istrip).energy();
00331
00333 if ( index <= 3 || index >= 283 ) {
00334 istrip++;
00335 nstrip++;
00336 continue;
00337 }
00338
00340 #ifdef SKIP_STAT
00341 if ( (*istrip).fail() || (*istrip).stat() ) {
00342 istrip++;
00343 nstrip++;
00344 continue;
00345 }
00346 #endif
00347
00350 if ( eseed < mSeedFloor*floor[ index ] + mSeedEnergy[4+plane] ) {
00351 istrip++;
00352 nstrip++;
00353 continue;
00354 }
00355
00357 seeds.push_back( (*istrip) );
00358 nstrip++;
00359 nseeds++;
00360
00361 #ifndef LOOSE_CUTS
00363 for ( Int_t i=0; i < 288; i++ ) {
00364
00365 Int_t dx=TMath::Abs(index-i);
00366
00368 if ( dx<3 )
00369 if ( eseed > floor[i] ) floor[i]=eseed;
00370
00372 if ( dx<5 )
00373 if ( 0.20 * eseed > floor[i] ) floor[i] = 0.2 * eseed;
00374
00376 if ( dx<11 )
00377 if ( 0.10 * eseed > floor[i] ) floor[i] = 0.1 * eseed;
00378
00380 if ( dx<21 )
00381 if ( 0.20*eseed > floor[i] ) floor[i] = 0.05 * eseed;
00382
00383 }
00384 #else
00386 for ( Int_t i=0; i < 288; i++ ) {
00387
00388 Int_t dx=TMath::Abs(index-i);
00389
00391 if ( dx<7 )
00392 if ( 0.05 * eseed > floor[i] ) floor[i] = 0.05 * eseed;
00393
00394 }
00395
00396 #endif
00397 istrip++;
00398
00399 }
00400
00401
00407 Bool_t owned[288]; for (Int_t i=0;i<288;i++) owned[i]=false;
00408 Bool_t xseed[288]; for (Int_t i=0;i<288;i++) xseed[i]=false;
00409
00410 StEEmcStripVec_t::iterator iseed=seeds.begin();
00411 while ( iseed != seeds.end() ) {
00412
00414 Int_t index=(*iseed).index();
00415 if ( owned[index] || (mSuppress&&xseed[index]) ) {
00416 iseed++;
00417 continue;
00418 }
00419
00421 owned[index]=true;
00422 xseed[index]=true;
00424 StEEmcSmdCluster cluster;
00426 cluster.add( (*iseed) );
00427
00430 for ( Int_t i=index+1; i<=index+max_extent; i++ ) {
00432 StEEmcStrip strip=mEEanalysis->strip(sector,plane,i);
00434 #ifdef SKIP_STAT
00435 if (strip.energy()<=0.&&!strip.fail()&&!strip.stat()) break;
00436 #else
00437 if ( strip.energy()<=0. ) break;
00438 #endif
00440 owned[ strip.index() ] = true;
00441 xseed[ strip.index() ] = true;
00443 cluster.add(strip);
00444 }
00445 for ( Int_t i=index-1; i>=index-max_extent; i-- ) {
00447 StEEmcStrip strip=mEEanalysis->strip(sector,plane,i);
00449 #ifdef SKIP_STAT
00450 if (strip.energy()<=0.&&!strip.fail()&&!strip.stat()) break;
00451 #else
00452 if (strip.energy()<=0.) break;
00453 #endif
00455 owned[ strip.index() ] = true;
00456 xseed[ strip.index() ] = true;
00458 cluster.add(strip);
00459 }
00460
00462 if ( cluster.size() >= 3 ) {
00463 cluster.key( mClusterId++ );
00464 mSmdClusters[ sector ][ plane ].push_back(cluster);
00465 mNumberOfClusters[4+plane]++;
00466
00467
00468
00469 Int_t ns=cluster.numberOfStrips();
00470 Int_t left=999,right=-999;
00471 for ( Int_t is=0;is<ns;is++ )
00472 {
00473 StEEmcStrip s=cluster.strip(is);
00474 if ( s.index()<left ) left=s.index();
00475 if ( s.index()>right ) right=s.index();
00476 }
00477
00478 xseed[left-1]=true;
00479 xseed[left-2]=true;
00480 xseed[right+1]=true;
00481 xseed[right+2]=true;
00482
00483 }
00484
00485
00486 iseed++;
00487 }
00488
00489 }
00490
00491 return true;
00492
00493 }
00494 #endif
00495 #ifdef DYNAMIC_SMD_SEED
00496
00497 Bool_t StEEmcClusterMaker::buildSmdClusters()
00498 {
00499
00500 Int_t max_extent = mMaxExtent;
00501
00502
00504 for ( Int_t sector=0; sector<12; sector++ )
00505 {
00506
00508 StEEmcClusterVec_t clist = clusters(sector,0);
00509
00511 StEEmcStripVec_t seeds;
00512
00514 for ( UInt_t iuv=0; iuv<2; iuv++ )
00515 {
00516
00518 StEEmcStripVec_t seeds;
00519
00521 std::vector<Int_t> used(288,0);
00522
00524 std::vector<Float_t> floor(288,0.);
00525
00526
00528 for ( UInt_t ic=0; ic<clist.size(); ic++ )
00529 {
00530
00532 StEEmcCluster cluster=clist[ic];
00533 Int_t sec=cluster.tower(0).sector();
00534 Int_t sub=cluster.tower(0).subsector();
00535 Int_t eta=cluster.tower(0).etabin();
00536
00539 Float_t esmd = cluster.energy() * 0.0056;
00540
00543 Float_t ethresh = esmd * 0.1;
00544
00548 Int_t smin[2], smax[2];
00549 mEEmap->getRangeU(sec,sub,eta, smin[0], smax[0] );
00550 mEEmap->getRangeV(sec,sub,eta, smin[1], smax[1] );
00551 smin[0]=TMath::Max(smin[0]-20,4);
00552 smax[0]=TMath::Min(smax[0]+20,284);
00553 smin[1]=TMath::Max(smin[1]-20,4);
00554 smax[1]=TMath::Min(smax[1]+20,284);
00555
00556
00559 for ( Int_t istrip=smin[iuv]; istrip<smax[iuv]; istrip++ )
00560 {
00561
00565 StEEmcStrip strip=mEEanalysis->strip(sector,iuv,istrip);
00566
00569 Float_t mythresh=TMath::Max(ethresh, floor[istrip]);
00570
00573 if ( strip.energy() > mythresh && !used[istrip] ) {
00574
00575 seeds.push_back( strip );
00576 used[istrip]++;
00577 floor[istrip] = strip.energy();
00578 for ( Int_t i=0; i<20; i++ ) {
00579 if ( istrip+i<288 ) floor[istrip+i]+=0.10*strip.energy();
00580 if ( istrip-i>0 ) floor[istrip-i]+=0.10*strip.energy();
00581 }
00582
00584 StEEmcSmdCluster smdc;
00585 smdc.key( mClusterId++ );
00586 Bool_t goL=1, goR=1;
00587 for ( Int_t i=0; i<max_extent; i++ ) {
00588 goL&=(istrip-i>0);
00589 goR&=(istrip+i<287);
00590
00591 if ( goL ) {
00592 StEEmcStrip s=mEEanalysis->strip(sector,iuv,istrip-i);
00593 if ( !s.stat() && !s.fail() ) {
00594 goL&=s.energy()>0.;
00595 if ( goL ) {
00596 smdc.add(s);
00597 used[istrip-i]=1;
00598 }
00599 }
00600 }
00601
00602 if ( goR ) {
00603 StEEmcStrip s=mEEanalysis->strip(sector,iuv,istrip+i);
00604 if ( !s.stat() && !s.fail() ) {
00605 goR&=s.energy()>0.;
00606 if ( goR ) {
00607 smdc.add(s);
00608 used[istrip+i]=1;
00609 }
00610 }
00611 }
00612
00613 }
00614
00615 if ( smdc.numberOfStrips() > 2 )
00616 mSmdClusters[ sector ][ iuv ].push_back( smdc );
00617
00618 }
00619
00620
00621 }
00622
00623
00624 }
00625
00626 }
00627
00628 }
00629
00630 }
00631
00632 #endif
00633
00634
00635
00636 void StEEmcClusterMaker::fillStEvent()
00637 {
00638
00639 StEvent *stevent=(StEvent*)GetInputDS("StEvent");
00640 if ( !stevent ) {
00641 Warning("fillStEvent","called, but no StEvent to be found");
00642 return;
00643 }
00644
00645 std::cout << "Adding tower clusters to StEvent at " << stevent << std::endl << std::flush;
00646
00650 StEmcDetector *detector=stevent->emcCollection()->detector(kEndcapEmcTowerId);
00651 if ( !detector )
00652 {
00653 Warning("fillStEvent","detector == NULL, MAJOR StEvent problem, continuing");
00654 return;
00655 }
00656
00657
00665
00666 if ( mNumberOfClusters[0] > 0 )
00667 {
00668
00669
00674 StEmcClusterCollection *collect = detector -> cluster();
00675 if ( !collect )
00676 {
00677
00678 collect = new StEmcClusterCollection();
00679 detector->setCluster( collect );
00680 }
00681
00682 assert(collect);
00683 collect->setDetector( kEndcapEmcTowerId );
00684 collect->setClusterFinderId( 123 );
00685 collect->setClusterFinderParamVersion( 123 );
00686
00688 for ( Int_t isector=0; isector<12; isector++ )
00689 {
00690
00692 for ( UInt_t iclust=0; iclust<mTowerClusters[isector][0].size(); iclust++ )
00693 {
00694
00695 StEEmcCluster cl=(mTowerClusters[isector][0])[iclust];
00696
00700 StEmcCluster *emccluster = new StEmcCluster();
00701 emccluster->setEta( cl.momentum().Eta() );
00702 emccluster->setPhi( cl.momentum().Phi() );
00703 emccluster->setSigmaEta(-1.);
00704 emccluster->setSigmaPhi(-1.);
00705 emccluster->setEnergy( cl.energy() );
00706 emccluster->SetUniqueID( cl.key() );
00707 #if 1
00708 for ( Int_t i=0; i< cl.numberOfTowers(); i++ )
00709 {
00710 StEmcRawHit *hit=cl.tower(i).stemc();
00711 assert( hit );
00712 emccluster->addHit( hit );
00713 }
00714 #endif
00715
00716 collect->addCluster( emccluster );
00717
00718
00720 mEtoEE[ emccluster ] = cl;
00721 cl.stemc( emccluster );
00722
00723 }
00724
00725 }
00726
00727 }
00728
00729 else
00730 {
00731
00732 detector->setCluster( NULL );
00733
00734 }
00735
00736
00737
00738
00740
00741
00745 detector=stevent->emcCollection()->detector(kEndcapEmcPreShowerId);
00746 if ( !detector )
00747 {
00748 Warning("fillStEvent","detector == NULL for pre/post, no clusters for you");
00749 }
00750 else if ( mNumberOfClusters[1] > 0 ||
00751 mNumberOfClusters[2] > 0 ||
00752 mNumberOfClusters[3] > 0 )
00753 {
00754
00755 StEmcClusterCollection *pqr = detector -> cluster();
00756 if ( !pqr )
00757 {
00758
00759 pqr = new StEmcClusterCollection();
00760 detector->setCluster( pqr );
00761 }
00762 assert(pqr);
00763 pqr -> setDetector( kEndcapEmcPreShowerId );
00764 pqr -> setClusterFinderId( 123 );
00765 pqr -> setClusterFinderParamVersion( 321 );
00766
00768 for ( Int_t isector=0; isector<12; isector++ )
00769
00771 for ( Int_t ilayer=1; ilayer<4; ilayer++ )
00772
00774 for ( UInt_t iclust=0; iclust<mTowerClusters[isector][ilayer].size(); iclust++ )
00775 {
00776
00777 StEEmcCluster cl=(mTowerClusters[isector][ilayer])[iclust];
00778
00782 StEmcCluster *emccluster = new StEmcCluster();
00783 emccluster->setEta( cl.momentum().Eta() );
00784 emccluster->setPhi( cl.momentum().Phi() );
00785 emccluster->setSigmaEta(-1.);
00786 emccluster->setSigmaPhi(-1.);
00787 emccluster->setEnergy( cl.energy() );
00788 emccluster->SetUniqueID( cl.key() );
00789 #if 1
00790 for ( Int_t i=0; i< cl.numberOfTowers(); i++ )
00791 {
00792 StEmcRawHit *hit=cl.tower(i).stemc();
00793 assert( hit );
00794 emccluster->addHit( hit );
00795 }
00796 #endif
00797
00798
00799 pqr->addCluster( emccluster );
00800
00801 mEtoEE[ emccluster ] = cl;
00802 cl.stemc( emccluster );
00803
00804 }
00805
00806 }
00807 else
00808 {
00809
00810 detector->setCluster( NULL );
00811
00812 }
00813
00814
00818 StDetectorId ids[]={ kEndcapSmdUStripId, kEndcapSmdVStripId };
00819
00820
00821 for ( Int_t iplane=0; iplane<2; iplane++ )
00822 {
00823
00824 detector=stevent->emcCollection()->detector(ids[iplane]);
00825 if ( !detector )
00826 {
00827 Warning("fillStEvent","detector == NULL for smd plane, no clusters for you");
00828 }
00829 else if ( mNumberOfClusters[4+iplane] > 0 )
00830 {
00831
00832 StEmcClusterCollection *smdc = detector -> cluster();
00833 if ( !smdc )
00834 {
00835
00836 smdc = new StEmcClusterCollection();
00837 detector->setCluster( smdc );
00838 }
00839
00840 smdc->setDetector( ids[iplane] );
00841 smdc->setClusterFinderId( 123 );
00842 smdc->setClusterFinderParamVersion( 321 );
00843
00844
00845 for ( Int_t isector=0; isector<12; isector++ )
00846 {
00847
00848 for ( UInt_t iclust=0; iclust<mSmdClusters[isector][iplane].size(); iclust++ )
00849 {
00850
00851 StEEmcSmdCluster cl = (mSmdClusters[isector][iplane])[iclust];
00852
00853
00854 StEmcCluster *emccluster = new StEmcCluster();
00855 emccluster->setEta( cl.mean() );
00856 emccluster->setPhi( (Float_t)iplane );
00857 emccluster->setSigmaEta(-1.);
00858 emccluster->setSigmaPhi(-1.);
00859 emccluster->setEnergy( cl.energy() );
00860 emccluster->SetUniqueID( cl.key() );
00861 for ( Int_t i=0; i< cl.numberOfStrips(); i++ )
00862 {
00863 StEmcRawHit *hit=cl.strip(i).stemc();
00864 assert( hit );
00865 emccluster->addHit( hit );
00866 }
00867 smdc->addCluster( emccluster );
00868
00869 mEtoEEsmd[ emccluster ] = cl;
00870 cl.stemc( emccluster );
00871
00872 }
00873
00874
00875 }
00876
00877 }
00878
00879 else
00880 {
00881 detector -> setCluster( NULL );
00882 }
00883
00884
00885
00886 }
00887
00888
00889 }
00890
00891
00892
00893
00894 Bool_t StEEmcClusterMaker::verifyStEvent()
00895 {
00896
00898 StEvent *stevent=(StEvent*)GetInputDS("StEvent");
00899 if ( !stevent ) {
00900 Warning("verifyStEvent","No StEvent found.");
00901 return true;
00902 }
00903
00904
00905 Bool_t go = true;
00906
00907 StEmcDetector *detector=stevent->emcCollection()->detector(kEndcapEmcTowerId);
00908 if ( !detector )
00909 {
00910 Warning("verifyStEvent","detector == NULL for towers");
00911 return false;
00912 }
00913
00914
00915 StEmcClusterCollection *cc=detector->cluster();
00916 if ( cc )
00917 {
00918
00919
00924
00925 Float_t emc_sum_towers = 0.;
00926 Float_t eemc_sum_towers = 0.;
00927
00928 StSPtrVecEmcCluster &emcClusters=cc->clusters();
00929 for ( UInt_t i=0; i<emcClusters.size(); i++ )
00930 {
00931 StEmcCluster *cl=emcClusters[i];
00932 assert(cl);
00933 emc_sum_towers += cl->energy();
00934 }
00935
00936 for ( Int_t sec=0; sec<12; sec++ )
00937 {
00938
00939 for ( Int_t i=0; i<numberOfClusters(sec,0); i++ )
00940 eemc_sum_towers += cluster(sec,0,i).energy();
00941 }
00942
00943 std::cout << "StEEmcClusterMaker tower checksum: ";
00944 if ( emc_sum_towers == eemc_sum_towers ) {
00945 std::cout << "passed";
00946 }
00947 else {
00948 std::cout << "FAILED"; go=false;
00949 }
00950 std::cout << std::endl;
00951
00952 }
00953 else {
00954
00955 std::cout << "StEEmcClusterMaker tower checksum: NULL collection, nclust=" << mNumberOfClusters[0] << std::endl;
00956 go &= (mNumberOfClusters[0]==0);
00957
00958 }
00959
00964 Float_t emc_sum_smdu = 0.;
00965 Float_t eemc_sum_smdu = 0.;
00966
00967 detector=stevent->emcCollection()->detector(kEndcapSmdUStripId);
00968 if ( !detector )
00969 {
00970 Warning("verifyStEvent","detector == NULL for smdu");
00971 return false;
00972 }
00973
00974 cc=detector->cluster();
00975 if ( cc )
00976 {
00977
00978 StSPtrVecEmcCluster &smduClusters=cc->clusters();
00979
00980 for ( UInt_t i=0; i<smduClusters.size(); i++ )
00981 {
00982 StEmcCluster *cl=smduClusters[i];
00983 assert(cl);
00984 emc_sum_smdu += cl->energy();
00985 }
00986
00987
00988 for ( Int_t sec=0; sec<12; sec++ )
00989 {
00990
00991 for ( Int_t i=0; i<numberOfSmdClusters(sec,0); i++ )
00992
00993 {
00994 eemc_sum_smdu += smdcluster(sec,0,i).energy();
00995
00996 }
00997
00998 }
00999
01000 std::cout << "StEEmcClusterMaker smdu checksum: ";
01001 if ( emc_sum_smdu == eemc_sum_smdu ) {
01002 std::cout << "passed";
01003 }
01004 else {
01005 std::cout << "FAILED"; go=false;
01006 }
01007 std::cout << std::endl;
01008
01009 }
01010 else
01011 {
01012 std::cout << "StEEmcClusterMaker smdu checksum: NULL collection, nclust=" << mNumberOfClusters[4] << std::endl;
01013 go &= (mNumberOfClusters[4]==0);
01014 }
01015
01016
01018
01019 Float_t emc_sum_smdv = 0.;
01020 Float_t eemc_sum_smdv = 0.;
01021
01022 detector=stevent->emcCollection()->detector(kEndcapSmdVStripId);
01023 if (!detector)
01024 {
01025 Warning("verifyStEvent","detector == NULL for smdv");
01026 return false;
01027 }
01028
01029 cc=detector->cluster();
01030
01031 if ( cc )
01032 {
01033
01034 StSPtrVecEmcCluster &smdvClusters=cc->clusters();
01035
01036 for ( UInt_t i=0; i<smdvClusters.size(); i++ )
01037 {
01038 StEmcCluster *cl=smdvClusters[i];
01039 assert(cl);
01040 emc_sum_smdv += cl->energy();
01041 }
01042
01043
01044 for ( Int_t sec=0; sec<12; sec++ )
01045 {
01046
01047 for ( Int_t i=0; i<numberOfSmdClusters(sec,1); i++ )
01048
01049 {
01050 eemc_sum_smdv += smdcluster(sec,1,i).energy();
01051
01052 }
01053
01054 }
01055
01056 std::cout << "StEEmcClusterMaker smdv checksum: ";
01057 if ( emc_sum_smdv == eemc_sum_smdv ) {
01058 std::cout << "passed";
01059 }
01060 else {
01061 std::cout << "FAILED"; go=false;
01062 }
01063 std::cout << std::endl;
01064
01065 }
01066 else
01067 {
01068 std::cout << "StEEmcClusterMaker smdv checksum: NULL collection, nclust=" << mNumberOfClusters[5] << std::endl;
01069 go &= (mNumberOfClusters[5]==0);
01070 }
01071
01072
01073
01074 return go;
01075 }
01076
01077
01078
01079
01080 void StEEmcClusterMaker::print()
01081 {
01082
01083 std::cout << "StEEmcClusterMaker::print()" << std::endl;
01084 const Char_t *names[] = { "tower", "pre1", "pre2", "post", "smdu", "smdv" };
01085 for ( Int_t i=0;i<6;i++ )
01086 {
01087
01088 std::cout << "Number of " << names[i]
01089 << " clusters = " << mNumberOfClusters[i]
01090 << std::endl;
01091
01092 }
01093
01094 std::cout << "printout of tower clusters follows:" << std::endl;
01095 for ( Int_t sec=0;sec<12;sec++)
01096 for ( Int_t i=0;i<numberOfClusters(sec,0);i++ )
01097 {
01098 StEEmcCluster clust=cluster(sec,0,i);
01099 clust.print();
01100
01101
01102
01103
01104 }
01105
01106
01107
01108 }
