HoughView Class Reference

#include <HoughView.h>

Inheritance diagram for HoughView:
Page

List of all members.

Public Member Functions

 HoughView ()
 ~HoughView ()
virtual void BuildDisplay (TCanvas *c)
virtual void DrawEvent (ParticleObjectHolder *poh=0, NtpStRecord *ntp=0, int ntpEvt=-1)
void DrawChains (int view)
void DrawTrueVertex (int view)
void DrawVertex (int view)
void DrawClusters (int view, TH2 *h)
void DrawHough (int view, TH2 *h)

Private Member Functions

void ClearHit (TH2D *his, double below_val, int curx, int cury)
void SaveHitGroup (TH2D *his, TH2D *saver, double save_val, double with_val, int curx, int cury)
void GetPeakAreaAverage (double &x, double &y, double &val, int &cnt, int curx, int cury, TH2D *hist)

Private Attributes

TPad * padU
TPad * padV
TH2D * viewU
TH2D * viewV
TH2D * houghviewU
TH2D * houghviewV
ParticleObjectHoldermypoh

Detailed Description

Definition at line 13 of file HoughView.h.


Constructor & Destructor Documentation

HoughView::HoughView (  ) 

Definition at line 24 of file HoughView.cxx.

References houghviewU, houghviewV, padU, padV, viewU, and viewV.

00024                     :Page()
00025 {
00026         viewU=0;
00027         viewV=0;
00028         padU=0;
00029         padV=0;
00030         houghviewU=0;
00031         houghviewV=0;
00032 
00033 }

HoughView::~HoughView (  ) 

Definition at line 36 of file HoughView.cxx.

00037 {}


Member Function Documentation

void HoughView::BuildDisplay ( TCanvas *  c  )  [virtual]

Reimplemented from Page.

Definition at line 40 of file HoughView.cxx.

References Page::myCanvas, padU, and padV.

Referenced by ParticleDisplay::BuildDisplay().

00041 {
00042         myCanvas=c;
00043         myCanvas->cd();
00044         padU=new TPad("houghview_padu","houghview_padu",0,0.5,0.7,1.0);
00045         padV=new TPad("houghview_padv","houghview_padv",0,0,0.7,0.5);
00046 
00047 
00048         padU->Draw();
00049         padV->Draw();
00050         
00051 
00052         
00053 }

void HoughView::ClearHit ( TH2D *  his,
double  below_val,
int  curx,
int  cury 
) [private]

Definition at line 932 of file HoughView.cxx.

00933 {
00934         if(curx<1 || cury < 1 || curx> his->GetNbinsX() || cury >his->GetNbinsY())return;
00935         
00936 //      printf("looking in %d %d\n",curx,cury);
00937         
00938         double thisval = his->GetBinContent(curx,cury);
00939         if(thisval>=below_val)return;
00940         
00941         
00942 //      printf("val %f >= %f ?\n",thisval,below_val);
00943         ClearHit(his,below_val,curx,cury+1);
00944         ClearHit(his,below_val,curx+1,cury);
00945         ClearHit(his,below_val,curx,cury-1);
00946         ClearHit(his,below_val,curx-1,cury);
00947         ClearHit(his,below_val,curx+1,cury+1);
00948         ClearHit(his,below_val,curx-1,cury-1);
00949         ClearHit(his,below_val,curx+1,cury-1);
00950         ClearHit(his,below_val,curx-1,cury+1);
00951 
00952 //      printf("recursion done\n");
00953 
00954         his->SetBinContent(curx,cury,0);        
00955         
00956 }

void HoughView::DrawChains ( int  view  ) 

Definition at line 174 of file HoughView.cxx.

References ParticleObjectHolder::chu, ParticleObjectHolder::chv, copy(), Chain::end_t, Chain::end_z, ChainHelper::FindMaxPath(), ChainHelper::finished, ChainHelper::GetAllChildren(), ChainHelper::GetChain(), Chain::myId, mypoh, Chain::parentChain, Chain::t, and Chain::z.

00175 {       
00176         
00177         ChainHelper * ch = 0;
00178     if(view==2)
00179         ch=&mypoh->chu;
00180     else if(view==3)
00181         ch=&mypoh->chv;
00182     else return;    
00183 
00184 
00185     std::vector<int> path =ch->FindMaxPath();
00186 
00187         for(unsigned int i=0;i<ch->finished.size();i++)
00188     {
00189             if(ch->finished[i].parentChain>-1 || ch->finished[i].entries < 1)continue;
00190         int parentcount=2;
00191         Chain *todraw=ch->GetChain(ch->finished[i].myId);
00192 
00193                 //draw the parent
00194             int dsize = todraw->t.size();
00195         int needconnect=0;
00196 
00197             if (dsize==1)  //chain with only 1 cluster needs to draw line to parent explicitly... if not a head chain
00198                                                 //other wise, make the single cluster, head chain, no children noticeable
00199          {
00200          //      if(todraw->parentChain<0)continue;
00201                  if(todraw->parentChain>-1)
00202              {
00203                      dsize=2;
00204                  needconnect=1;
00205              }
00206          }
00207 
00208          double * tt = new double[todraw->t.size()];
00209          double * tz = new double[todraw->t.size()];
00210          if(needconnect)
00211          {
00212                  tt[0]=ch->GetChain(todraw->parentChain)->end_t;
00213              tz[0]=ch->GetChain(todraw->parentChain)->end_z;
00214              tz[1]=todraw->z[0];
00215              tt[1]=todraw->t[0];
00216          }
00217          else
00218          {
00219              std::copy(todraw->t.begin(),todraw->t.end(),tt);
00220              std::copy(todraw->z.begin(),todraw->z.end(),tz);
00221                 
00222              if(dsize==1)
00223              {
00224                      tt[1]=tt[0];
00225                  tz[1]=tz[0]-0.03;
00226                  dsize=2;
00227 
00228              }
00229 
00230          }
00231 
00232          TPolyLine *pl = new TPolyLine(dsize,tz,tt );
00233 
00234          //int linecol =2+todraw->level;
00235 
00236         //      pl->SetLineColor(linecol);
00237         //      pl->SetLineColor(parentcount);
00238                  pl->SetLineColor(2);
00239 
00240          int width=1;
00241          //should make this better if it is to be really used
00242 
00243          for(unsigned int j=0;j<path.size();j++)
00244          {
00245                 if(path[j]==todraw->myId)width+=1;
00246          }
00247 
00248          pl->SetLineWidth(width);
00249          //    pl->Draw();
00250          pl->Draw("same");
00251                         
00252          for (unsigned int k=0;k<ch->finished[i].children.size();k++)
00253          {
00254                 parentcount++;
00255                 
00256             std::vector<int> c = ch->GetAllChildren(ch->GetChain(ch->finished[i].children[k]));
00257                         
00258             printf("parent %d children ", ch->finished[i].myId);
00259                        
00260             for(unsigned int j=0;j<c.size();j++)
00261             {               
00262                                 
00263                 Chain *todraw=ch->GetChain(c[j]);
00264 
00265                 printf("%d ",todraw->myId);
00266 
00267                 
00268                 int dsize = todraw->t.size();
00269                 int needconnect=0;
00270 
00271                 if (dsize==1)  //chain with only 1 cluster needs to draw line to parent explicitly... if not a head chain
00272                                                 //other wise, make the single cluster, head chain, no children noticeable 
00273                 {
00274                 //      if(todraw->parentChain<0)continue;
00275                 
00276                         if(todraw->parentChain>-1) 
00277                         {
00278                                 dsize=2;
00279                                 needconnect=1;
00280                         }
00281                 }
00282 
00283                 double * tt = new double[todraw->t.size()];
00284                 double * tz = new double[todraw->t.size()];
00285                 if(needconnect)
00286                 {
00287                         tt[0]=ch->GetChain(todraw->parentChain)->end_t;
00288                         tz[0]=ch->GetChain(todraw->parentChain)->end_z;
00289                         tz[1]=todraw->z[0];
00290                         tt[1]=todraw->t[0];
00291                 }
00292                 else
00293                 {
00294                         std::copy(todraw->t.begin(),todraw->t.end(),tt);
00295                         std::copy(todraw->z.begin(),todraw->z.end(),tz);
00296 
00297                         if(dsize==1)
00298                         {
00299                                 tt[1]=tt[0];
00300                                 tz[1]=tz[0]-0.03;
00301                                 dsize=2;
00302 
00303                         }
00304 
00305                 }
00306 
00307                 TPolyLine *pl = new TPolyLine(dsize,tz,tt );
00308 
00309                 //int linecol =2+todraw->level;
00310 
00311         //      pl->SetLineColor(linecol);
00312                 pl->SetLineColor(parentcount);
00313 
00314 
00315                 int width=1;
00316                 //should make this better if it is to be really used
00317 
00318                 for(unsigned int j=0;j<path.size();j++)
00319                 {
00320                         if(path[j]==todraw->myId)width+=1;
00321                 }
00322 
00323                 pl->SetLineWidth(width);
00324         //      pl->Draw();
00325                 pl->Draw("same");
00326  
00327                 }printf("\n");
00328         }}      
00329                 
00330 }       

void HoughView::DrawClusters ( int  view,
TH2 *  h 
)

Definition at line 336 of file HoughView.cxx.

References ParticleObjectHolder::cluster_saver, ParticleObjectHolder::clusterer, Managed::ManagedCluster::e, Managed::ClusterSaver::GetCluster(), Managed::ClusterManager::GetCluster(), Managed::ClusterSaver::GetClusterMap(), Managed::ClusterManager::GetClusterMap(), Managed::ManagedCluster::GetStatus(), Munits::m, mypoh, Managed::ManagedCluster::t, and Managed::ManagedCluster::z.

Referenced by DrawEvent().

00337 {
00338 
00339          std::map<double, std::map<double, int>  >::iterator p_iterr;
00340      std::map<double, int>::iterator s_iterr;
00341         
00342 
00343         Managed::ClusterManager *cluh = &mypoh->clusterer;
00344 /*
00345              for(p_iterr=cluh->GetClusterMap(view)->begin();p_iterr!=cluh->GetClusterMap(view)->end(); p_iterr++)
00346          {
00347          for(s_iterr=p_iterr->second.begin();s_iterr!=p_iterr->second.end(); s_iterr++)
00348          {
00349          
00350                 
00351                 h->Fill( p_iterr->first,s_iterr->first, cluh->GetCluster(s_iterr->second)->e);  
00352          }
00353          }
00354                 h->SetStats(0);
00355                 h->Draw("colz");
00356 */
00357 
00358   //Managed::ManagedCluster *largest=0;
00359   //double last_plane=0;
00360     for(p_iterr=cluh->GetClusterMap(view)->begin();p_iterr!=cluh->GetClusterMap(view)->end(); p_iterr++)
00361     {   
00362 
00363         for(s_iterr=p_iterr->second.begin();s_iterr!=p_iterr->second.end(); s_iterr++)
00364         {
00365                         //Managed::ManagedCluster * clus = cluh->GetCluster(s_iterr->second);
00366                         h->Fill( p_iterr->first,s_iterr->first, cluh->GetCluster(s_iterr->second)->e); 
00367 
00368                 }
00369         }
00370 
00371 
00372         
00373         Managed::ClusterSaver * cluhs = &mypoh->cluster_saver;
00374     for(p_iterr=cluhs->GetClusterMap(view)->begin();p_iterr!=cluhs->GetClusterMap(view)->end(); p_iterr++)
00375     {   
00376 
00377         for(s_iterr=p_iterr->second.begin();s_iterr!=p_iterr->second.end(); s_iterr++)
00378         {
00379                         Managed::ManagedCluster * clus = cluhs->GetCluster(s_iterr->second);
00380                         
00381                         //printf("%f %f %d\n",clus->z,clus->t,clus->GetStatus());
00382                         if(clus->GetStatus()==1)continue;//long muon clusters
00383                         
00384                         h->Fill( p_iterr->first,s_iterr->first, cluhs->GetCluster(s_iterr->second)->e); 
00385 
00386                 }
00387         }
00388         h->SetStats(0);
00389         h->Draw("colz");       
00390         
00391         
00392          
00393 /*
00394     for(p_iterr=cluh->GetClusterMap(view)->begin();p_iterr!=cluh->GetClusterMap(view)->end(); p_iterr++)
00395     {   
00396         //new plane and a largest cluster in the previous plane?
00397                 if(fabs(last_plane-p_iterr->first)>0.03 && largest)
00398                 {
00399                         TMarker *m=new TMarker(largest->z,largest->t, 25);
00400                         m->Draw("same");
00401                         largest=0;
00402                         last_plane=p_iterr->first;
00403                 }
00404         
00405         for(s_iterr=p_iterr->second.begin();s_iterr!=p_iterr->second.end(); s_iterr++)
00406         {
00407                         Managed::ManagedCluster * clus = cluh->GetCluster(s_iterr->second);
00408 
00409                         if(!largest)
00410                         {
00411                                 largest=clus;
00412                                 continue;
00413                         }
00414                         if(largest->e<clus->e)largest=clus;
00415                 }
00416         }
00417         if(fabs(last_plane-largest->z)>0.03 && largest)
00418         {
00419                 TMarker *m=new TMarker(largest->z,largest->t, 25);
00420                 m->Draw("same");
00421                 largest=0;
00422         }    
00423 */
00424 
00425 //to take all within xxx% of the max hit in the event!
00426 
00427 //find max value!
00428 double maxvale=0;
00429     for(p_iterr=cluh->GetClusterMap(view)->begin();p_iterr!=cluh->GetClusterMap(view)->end(); p_iterr++)
00430     {   
00431         for(s_iterr=p_iterr->second.begin();s_iterr!=p_iterr->second.end(); s_iterr++)
00432         {
00433                         Managed::ManagedCluster * clus = cluh->GetCluster(s_iterr->second);
00434                         if(clus->e>maxvale)
00435                         {
00436                                 maxvale=clus->e;
00437                         }
00438                 }
00439         }
00440 
00441 
00442     for(p_iterr=cluh->GetClusterMap(view)->begin();p_iterr!=cluh->GetClusterMap(view)->end(); p_iterr++)
00443     {   
00444         for(s_iterr=p_iterr->second.begin();s_iterr!=p_iterr->second.end(); s_iterr++)
00445         {
00446                         Managed::ManagedCluster * clus = cluh->GetCluster(s_iterr->second);
00447 
00448                 //      if(clus->e>maxvale*0.1)
00449         //              {
00450                                                 if(clus->GetStatus()==1)continue;//long muon clusters
00451                 TMarker *m=new TMarker(clus->z,clus->t, 25);
00452                 m->Draw("same");
00453                         
00454                 //      }
00455                 }
00456         }
00457         
00458         
00459 /*
00460     for(p_iterr=cluhs->GetClusterMap(view)->begin();p_iterr!=cluhs->GetClusterMap(view)->end(); p_iterr++)
00461     {   
00462         for(s_iterr=p_iterr->second.begin();s_iterr!=p_iterr->second.end(); s_iterr++)
00463         {
00464                         Managed::ManagedCluster * clus = cluhs->GetCluster(s_iterr->second);
00465 
00466                         if(clus->GetStatus()==1)continue;//long muon clusters
00467                 TMarker *m=new TMarker(clus->z,clus->t, 27);
00468                 m->Draw("same");
00469                         
00470 
00471                 }
00472         }
00473 */      
00474         
00475         
00476 
00477 }

void HoughView::DrawEvent ( ParticleObjectHolder poh = 0,
NtpStRecord ntp = 0,
int  ntpEvt = -1 
) [virtual]

Reimplemented from Page.

Definition at line 56 of file HoughView.cxx.

References DrawClusters(), DrawHough(), DrawTrueVertex(), DrawVertex(), ParticleObjectHolder::event, RecRecordImp< T >::GetHeader(), VldContext::GetSimFlag(), RecHeader::GetVldContext(), houghviewU, houghviewV, isMC, SimFlag::kMC, ParticleEvent::maxu, ParticleEvent::maxv, ParticleEvent::maxz, ParticleObjectHolder::mctrue, ParticleEvent::minu, ParticleEvent::minv, ParticleEvent::minz, Page::myCanvas, mypoh, padU, padV, viewU, viewV, MCTrue::vtx_u, MCTrue::vtx_v, and MCTrue::vtx_z.

Referenced by ParticleDisplay::UpdateDisplay().

00057 {
00058         myCanvas->cd();
00059 
00060         padU->Clear();
00061         padV->Clear();
00062 
00063         padU->Divide(2,1);
00064         padV->Divide(2,1);
00065 
00066         mypoh=poh;
00067         
00068         int isMC=0;
00069         if(poh->GetHeader().GetVldContext().GetSimFlag() == SimFlag::kMC)isMC=1;
00070         
00071         double minz= mypoh->event.minz;
00072         double minu= mypoh->event.minu;
00073         double minv= mypoh->event.minv;
00074         double maxz= mypoh->event.maxz;
00075         double maxu= mypoh->event.maxu;
00076         double maxv= mypoh->event.maxv;
00077         
00078 
00079         
00080         if(isMC)
00081         {
00082                 if(mypoh->mctrue.vtx_z>=0)minz=minz<mypoh->mctrue.vtx_z ? minz : mypoh->mctrue.vtx_z;
00083                 maxz=maxz>mypoh->mctrue.vtx_z ? maxz : mypoh->mctrue.vtx_z;
00084                 minu=minu<mypoh->mctrue.vtx_u ? minu : mypoh->mctrue.vtx_u;
00085                 maxu=maxu>mypoh->mctrue.vtx_u ? maxu : mypoh->mctrue.vtx_u;
00086                 minv=minv<mypoh->mctrue.vtx_v ? minv : mypoh->mctrue.vtx_v;
00087                 maxv=maxv>mypoh->mctrue.vtx_v ? maxv : mypoh->mctrue.vtx_v;
00088         }
00089         
00090 
00091 
00092 
00093         minz -= 0.1;
00094         minu -= 0.1;
00095         minv -= 0.1;
00096         maxz += 0.1;
00097         maxu += 0.1;
00098         maxv += 0.1;
00099         
00100         int nu=(int)((maxu-minu)/0.0412);
00101         int nv=(int)((maxv-minv)/0.0412);
00102         int nz=(int)((maxz-minz)/0.0354);
00103 
00104 
00105         if(viewU)delete viewU;
00106         if(viewV)delete viewV;
00107 
00108         if(houghviewU)delete houghviewU;
00109         if(houghviewV)delete houghviewV;
00110 
00111         viewU = new TH2D("hu","Clusters U", nz, minz,maxz,nu,minu, maxu);
00112         viewV = new TH2D("hv","Clusters V", nz, minz,maxz,nv,minv, maxv);
00113         houghviewU = new TH2D("hhu","Hough U", 80, 0,3.141592,50,-1,1);
00114         houghviewV = new TH2D("hhv","Hough V", 80, 0,3.141592,50,-1,1);
00115         viewU->SetStats(0);
00116         viewV->SetStats(0);
00117         viewU->SetDirectory(0);
00118         viewV->SetDirectory(0);
00119 
00120         houghviewU->SetStats(0);
00121         houghviewV->SetStats(0);
00122         houghviewU->SetDirectory(0);
00123         houghviewV->SetDirectory(0);
00124 
00125         if(mypoh)
00126         {       
00127         
00128 
00129                 //psf.FindPrimaryShower(&mypoh->clusterer);
00130         
00131         
00132                 padU->cd(1);    
00133                 
00134                 viewU->Draw();
00135                 DrawClusters(2,viewU);
00136                 DrawVertex(2);
00137                 if(isMC)DrawTrueVertex(2);
00138         //      DrawChains(2);
00139 
00140                 padU->cd(2);
00141                 DrawHough(2,houghviewU);
00142 
00143                 
00144 
00145                 padV->cd(1);
00146                 viewV->Draw();
00147                 DrawClusters(3,viewV);
00148                 DrawVertex(3);
00149                 if(isMC)DrawTrueVertex(3);
00150         //      DrawChains(3);  
00151                 padV->cd(2);
00152                 DrawHough(3,houghviewV);
00153 
00154 
00155 
00156 
00157         }
00158 
00159         
00160         double maxz_u = viewU->GetMaximum();
00161         double maxz_v = viewV->GetMaximum();
00162         double maxz_both = maxz_u>maxz_v?maxz_u:maxz_v;
00163         viewU->GetZaxis()->SetRangeUser(0, maxz_both);
00164         viewV->GetZaxis()->SetRangeUser(0, maxz_both);
00165         
00166         padU->Update();
00167         padV->Update();
00168 
00169         myCanvas->Update();
00170         
00171 }

void HoughView::DrawHough ( int  view,
TH2 *  h 
)

the old code....

Definition at line 518 of file HoughView.cxx.

References PrimaryShowerFinder::FindFirstIntersection(), PrimaryShowerFinder::GetHoughLines(), PrimaryShowerFinder::intU, PrimaryShowerFinder::intV, Munits::m, max, mypoh, padU, padV, ParticleObjectHolder::psf, and PrimaryShowerFinder::ran.

Referenced by DrawEvent().

00519 {
00520 
00521         if(view==2)padU->cd(2);
00522         if(view==3)padV->cd(2);
00523         
00524         double t=0;
00525         double z=0;
00526         
00527         if(!mypoh->psf.ran)return;
00528         int cp = mypoh->psf.FindFirstIntersection(view,t,z);
00529         if(cp)
00530         {
00531                 TMarker *m=new TMarker(z,t, 29);
00532                 m->Draw("same");
00533                 printf("closest intersection z %f t %f\n",z,t);
00534         }
00535 
00536 
00537 
00538         TH2D * inth = view==2? mypoh->psf.intU:mypoh->psf.intV;
00539         if(inth)
00540         {
00541                 inth->SetStats(0);
00542                 inth->Draw("colz");
00543         }
00544         
00545 /*      
00546 
00547     TH2D * hh = mypoh->psf.GetHoughMap(view);           
00548         hh->SetStats(0);
00549         hh->Draw("colz");
00550 */
00551 
00552         if(view==2)padU->cd(1);
00553         if(view==3)padV->cd(1);
00554 
00555 
00556 
00557 
00558         //double view_xmin = view==2?viewU->GetXaxis()->GetXmin():viewV->GetXaxis()->GetXmin();
00559         //double view_xmax = view==2?viewU->GetXaxis()->GetXmax():viewV->GetXaxis()->GetXmax();
00560         
00561 
00562         
00563         
00564 
00565 
00566         std::vector<HoughLine>*houghlines = mypoh->psf.GetHoughLines(view);
00567         
00568         printf("# hough lines %d\n",(int)houghlines->size());
00569         
00570         int max =-1;
00571         if((*houghlines).size()>50)max = (*houghlines).size()-50;
00572         for(int i=(*houghlines).size()-1;i>max;i--)
00573         {
00574                 if((*houghlines)[i].ncluster<1)continue;
00575 //              printf("houghline with e %f clusters %d (t,z) (%f,%f) (%f,%f) chi2/ndf %f  offz %f expt %f phi %f\n", (*houghlines)[i].sum_e, (*houghlines)[i].ncluster, (*houghlines)[i].start_t, (*houghlines)[i].start_z, (*houghlines)[i].end_t, (*houghlines)[i].end_z, (*houghlines)[i].chi2/(*houghlines)[i].ncluster,(*houghlines)[i].offset_z,(*houghlines)[i].GetExpectedT((*houghlines)[i].offset_z),(*houghlines)[i].phi);
00576                 
00577                 /*
00578                         double ym = (*houghlines)[i].GetExpectedT(view_xmin);
00579                         double yp = (*houghlines)[i].GetExpectedT(view_xmax);
00580                         
00581                                                 
00582                         TLine *l=new TLine(view_xmin,ym,view_xmax,yp);
00583                 */      
00584                 
00585                         TLine *l=new TLine((*houghlines)[i].start_z,(*houghlines)[i].start_t,(*houghlines)[i].end_z,(*houghlines)[i].end_t);
00586                 
00587                         if((*houghlines)[i].primary)l->SetLineWidth(2);
00588                         l->Draw();
00589                                                 
00590                 //      printf("fit line  --%f %f - %f %f -- theta r %f %f\n",view_xmin,ym,view_xmax,yp,houghlines[i].theta,houghlines[i].r);
00591                                                 
00592         }
00593         
00594 
00595         
00596 
00597 
00598 
00599 
00600 
00601 
00602 
00603 
00604 return;
00605 
00607 /*
00608 
00609          std::map<double, std::map<double, int>  >::iterator p_iterr;
00610      std::map<double, int>::iterator s_iterr;
00611         
00612 printf("----view %d----\n\n",view);
00613          Managed::ClusterManager *cluh = view==2?&mypoh->clusterer:&mypoh->clusterer;
00614 
00615 
00616                 double off_z = cluh->GetClusterMap(view)->begin()->first;
00617                 double off_t = cluh->GetClusterMap(view)->begin()->second.begin()->first;
00618 
00619 
00620 / *
00621              for(p_iterr=cluh->cluster_map.begin();p_iterr!=cluh->cluster_map.end(); p_iterr++)
00622          {
00623          for(s_iterr=p_iterr->second.begin();s_iterr!=p_iterr->second.end(); s_iterr++)
00624          {
00625          
00626                 for(double i=0;i<2*3.141592;i+=3.141592*2/99)
00627                 {
00628                         h->Fill(i, (p_iterr->first-off_z) * cos(i) + (s_iterr->first-off_t) * sin(i), 1cluh->GetCluster(s_iterr->second)->e);  
00629                 }
00630          }
00631          }
00632          
00633   *//*
00634   
00635   Managed::ManagedCluster *largest=0;
00636   double last_plane=0;
00637   
00638 //to take the largest in each plane  
00639 / *  
00640     for(p_iterr=cluh->GetClusterMap(view)->begin();p_iterr!=cluh->GetClusterMap(view)->end(); p_iterr++)
00641     {   
00642         //new plane and a largest cluster in the previous plane?
00643                 if(fabs(last_plane-p_iterr->first)>0.03 && largest)
00644                 {
00645                         for(double i=0;i<3.141592;i+=3.141592/(h->GetNbinsX()-1))
00646                 {
00647                         //double val = largest->e;
00648                                 double val=1;
00649                                 h->Fill(i, (largest->z-off_z) * cos(i) + (largest->t-off_t) * sin(i), val);  
00650 
00651 //                              h->Fill(i, (largest->z-off_z) * cos(i) + (largest->t-off_t) * sin(i)+2./(h->GetNbinsY())., val*.5); 
00652 //                              h->Fill(i, (largest->z-off_z) * cos(i) + (largest->t-off_t) * sin(i)+4./(h->GetNbinsY())., val*.5*.5);
00653 //                              h->Fill(i, (largest->z-off_z) * cos(i) + (largest->t-off_t) * sin(i)-2./(h->GetNbinsY())., val*.5); 
00654 //                              h->Fill(i, (largest->z-off_z) * cos(i) + (largest->t-off_t) * sin(i)-4./(h->GetNbinsY())., val*.5*.5);      
00655                         }
00656                         largest=0;
00657                         last_plane=p_iterr->first;
00658                 }
00659         
00660         for(s_iterr=p_iterr->second.begin();s_iterr!=p_iterr->second.end(); s_iterr++)
00661         {
00662                         Managed::ManagedCluster * clus = cluh->GetCluster(s_iterr->second);
00663                         if(!largest)
00664                         {
00665                                 largest=clus;
00666                                 continue;
00667                         }
00668                         if(largest->e<clus->e)largest=clus;
00669                 }
00670         }
00671         if(fabs(last_plane-largest->z)>0.03 && largest)
00672         {
00673                 for(double i=0;i<3.141592;i+=3.141592/(h->GetNbinsX()-1))
00674                 {
00675                         //double val = largest->e;
00676                         double val=1;
00677                         h->Fill(i, (largest->z-off_z) * cos(i) + (largest->t-off_t) * sin(i),  val);  
00678 //                              h->Fill(i, (largest->z-off_z) * cos(i) + (largest->t-off_t) * sin(i)+2./(h->GetNbinsY())., val*.5); 
00679 //                              h->Fill(i, (largest->z-off_z) * cos(i) + (largest->t-off_t) * sin(i)+4./(h->GetNbinsY())., val*.5*.5);
00680 //                              h->Fill(i, (largest->z-off_z) * cos(i) + (largest->t-off_t) * sin(i)-2./(h->GetNbinsY())., val*.5); 
00681 //                              h->Fill(i, (largest->z-off_z) * cos(i) + (largest->t-off_t) * sin(i)-4./(h->GetNbinsY())., val*.5*.5);      
00682                 
00683                 }
00684 
00685                 largest=0;
00686         }
00687 */       /*
00688          
00689          
00690         TH2D copy;
00691         h->Copy(copy);    
00692 
00693 //to take all within xxx% of the max hit in the event!
00694 
00695 //find max value!
00696 double maxvale=0;
00697     for(p_iterr=cluh->GetClusterMap(view)->begin();p_iterr!=cluh->GetClusterMap(view)->end(); p_iterr++)
00698     {   
00699         for(s_iterr=p_iterr->second.begin();s_iterr!=p_iterr->second.end(); s_iterr++)
00700         {
00701                         Managed::ManagedCluster * clus = cluh->GetCluster(s_iterr->second);
00702                         if(clus->e>maxvale)
00703                         {
00704                                 maxvale=clus->e;
00705                         }
00706                 }
00707         }
00708 
00709 
00710     for(p_iterr=cluh->GetClusterMap(view)->begin();p_iterr!=cluh->GetClusterMap(view)->end(); p_iterr++)
00711     {   
00712         for(s_iterr=p_iterr->second.begin();s_iterr!=p_iterr->second.end(); s_iterr++)
00713         {
00714                         Managed::ManagedCluster * clus = cluh->GetCluster(s_iterr->second);
00715 
00716 //                      if(clus->e>maxvale*0.1)
00717 //                      {
00718                 for(double i=0;i<3.141592;i+=3.141592/(h->GetNbinsX()-1))
00719                 {
00720                         //double val = clus->e;
00721                         double val=1;
00722                         h->Fill(i, (clus->z-off_z) * cos(i) + (clus->t-off_t) * sin(i),  val);  
00723                         copy.Fill(i, (clus->z-off_z) * cos(i) + (clus->t-off_t) * sin(i), clus->e);
00724 //                              h->Fill(i, (largest->z-off_z) * cos(i) + (largest->t-off_t) * sin(i)+2./(h->GetNbinsY())., val*.5); 
00725 //                              h->Fill(i, (largest->z-off_z) * cos(i) + (largest->t-off_t) * sin(i)+4./(h->GetNbinsY())., val*.5*.5);
00726 //                              h->Fill(i, (largest->z-off_z) * cos(i) + (largest->t-off_t) * sin(i)-2./(h->GetNbinsY())., val*.5); 
00727 //                              h->Fill(i, (largest->z-off_z) * cos(i) + (largest->t-off_t) * sin(i)-4./(h->GetNbinsY())., val*.5*.5);      
00728                         
00729 }                       
00730 //                      }
00731                 }
00732         }
00733         
00734         
00735         / *
00736         h->SetStats(0);
00737         h->Draw("colz");
00738          
00739          return;
00740   */
00741    //peak filter the hough map!
00742    /*
00743    TH2D ch;
00744    h->Copy(ch);
00745         
00746 
00747 h->Reset();
00748         
00749         int nx=ch.GetNbinsX();
00750         int ny=ch.GetNbinsY();   
00751 
00752         while(ch.GetMaximum()>1)
00753         {
00754                 int x;
00755                 int y;
00756                 int z;
00757                 int m = ch.GetMaximumBin();
00758                 ch.GetBinXYZ(m,x,y,z);
00759                 
00760 / *             printf("clearing hits around %d %d\n",x,y);
00761                 ClearHit(&ch, ch.GetMaximum(), x+1,y);
00762                 ClearHit(&ch, ch.GetMaximum(), x,y+1);
00763                 ClearHit(&ch, ch.GetMaximum(), x,y-1);
00764                 ClearHit(&ch, ch.GetMaximum(), x-1,y);
00765                 ClearHit(&ch, ch.GetMaximum(), x+1,y+1);
00766                 ClearHit(&ch, ch.GetMaximum(), x-1,y-1);
00767                 ClearHit(&ch, ch.GetMaximum(), x+1,y-1);
00768                 ClearHit(&ch, ch.GetMaximum(), x-1,y+1);
00769 */                              /*              
00770                 SaveHitGroup(&ch,(TH2D*)h, (double)ch.GetBinContent(x,y),(double)ch.GetBinContent(x,y), x,y);
00771                 
00772                 
00773     }    
00774                 
00775         h->SetStats(0);
00776         h->Draw("colz");
00777 
00778 
00779         if(view==2)padU->cd(1);
00780         if(view==3)padV->cd(1);
00781         
00782         
00783         / *
00784         for(int i=1;i<=h->GetNbinsX();i++)
00785         for(int j=1;j<=h->GetNbinsY();j++)
00786         {
00787                 if(h->GetBinContent(i,j)>=4)
00788                 {
00789                         
00790                         double theta = h->GetXaxis()->GetBinLowEdge(i);
00791                         double r = h->GetYaxis()->GetBinLowEdge(j);
00792                         
00793                         double ym = (- cos(theta)/sin(theta))*(viewU->GetXaxis()->GetXmin()-off_z)+r/sin(theta) + off_t;
00794                         double yp = (- cos(theta)/sin(theta))*(viewU->GetXaxis()->GetXmax()-off_z)+r/sin(theta) + off_t;
00795                                                 
00796                         TLine *l=new TLine(viewU->GetXaxis()->GetXmin(),ym,viewU->GetXaxis()->GetXmax(),yp);
00797                         l->Draw();
00798                 
00799                 }
00800         }
00801         */
00802         /*
00803         multimap<double, std::pair<double,double> > top_map;
00804         
00805         int tops=0;     
00806 / *     for(int i=1;i<=h->GetNbinsX();i++)
00807         for(int j=1;j<=h->GetNbinsY();j++)
00808         {
00809                 if(h->GetBinContent(i,j)>1)
00810                 {
00811                         top_map.insert(make_pair((double)h->GetBinContent(i,j), std::pair<int,int>(i,j) ));
00812                         tops++;
00813                         printf("topbin found x %d y %d\n",i,j);
00814                 }
00815         }
00816 */
00817 
00818 //h->Rebin2D(4,4);
00819 /*      while(h->GetMaximum()>0)
00820         {
00821                 double xv=0;
00822                 double yv=0;
00823                 double val=0;
00824                 int cnt =0;
00825                 
00826                 int x;
00827                 int y;
00828                 int z;
00829                 int m = h->GetMaximumBin();
00830                 h->GetBinXYZ(m,x,y,z);
00831                 
00832                 GetPeakAreaAverage(xv, yv,val, cnt, x, y, (TH2D*)h);
00833                 xv/=(double)cnt;
00834                 yv/=(double)cnt;
00835                 
00836                 val=copy.GetBinContent(copy.FindBin(xv,yv));
00837                 
00838                 top_map.insert(make_pair((double)val, std::pair<double,double>(xv,yv)) );
00839                 tops++;
00840         //      printf("topbin found x %f y %f\n",xv,yv);
00841 
00842         
00843         }
00844         
00845         
00846         
00847         printf("printing top %d\n",tops);
00848         
00849         multimap<double, std::pair<double,double> >::reverse_iterator it = top_map.rbegin();
00850         
00851         double view_xmin = view==2?viewU->GetXaxis()->GetXmin():viewV->GetXaxis()->GetXmin();
00852         double view_xmax = view==2?viewU->GetXaxis()->GetXmax():viewV->GetXaxis()->GetXmax();
00853         
00854         std::vector<HoughLine> houghlines;
00855         
00856         for(unsigned int k=0;k<tops;k++)
00857         {
00858                 //int i=it->second.first;
00859                 //int j=it->second.second;
00860                 
00861                         //double theta = h->GetXaxis()->GetBinCenter(i);
00862                         //double r = h->GetYaxis()->GetBinCenter(j);
00863                 double theta=it->second.first;
00864                 double r = it->second.second;   
00865                         
00866 //                      double ym = (- cos(theta)/sin(theta))*(view_xmin-off_z)+r/sin(theta) + off_t;
00867 //                      double yp = (- cos(theta)/sin(theta))*(view_xmax-off_z)+r/sin(theta) + off_t;
00868 
00869 
00870                         HoughLine hl(theta, r, off_t, off_z);
00871         
00872                         PrimaryShowerFinder::LoadCloseHits(&hl,&mypoh->clusterer,view,0.03);
00873                         
00874                         houghlines.push_back(hl);
00875 
00876                 it++;
00877         }
00878         
00879         sort(&houghlines[0],&houghlines[houghlines.size()],CompareLength);
00880         
00881         //should remove duplicates --- (caused by the same clusters being matched to different hough map peaks)
00882         std::vector<HoughLine> houghlinestemp;
00883         for(int i=0;i<houghlines.size();i++)
00884         {
00885                 HoughLine * l = &houghlines[i];
00886                 HoughLine * r=0;
00887                 if(i<houghlines.size()-1)r = &houghlines[i+1];
00888                 if(houghlines[i].ncluster>1)houghlinestemp.push_back(houghlines[i]);            
00889                 if(l && r && l->start_z == r->start_z && l->start_t == r->start_t && l ->end_z == r->end_z && l->end_t == r->end_t)
00890                 {
00891 
00892                         while (i<houghlines.size())
00893                         {
00894                                 i++;
00895                                 HoughLine * l = &houghlines[i];
00896                                 HoughLine * r = &houghlines[i+1];
00897                                 if(! (l->start_z == r->start_z && l->start_t == r->start_t && l ->end_z == r->end_z && l->end_t == r->end_t)) break;
00898                         }
00899                 }
00900         
00901         }
00902         
00903         houghlines = houghlinestemp;
00904         
00905         sort(&houghlines[0],&houghlines[houghlines.size()],CompareForwardAndClusters);
00906         
00907         int max =-1;
00908         if(houghlines.size()>10)max = houghlines.size()-10;
00909         for(int i=houghlines.size()-1;i>max;i--)
00910         {
00911                 if(houghlines[i].ncluster<1)break;
00912                 printf("houghline with e %f clusters %d (t,z) (%f,%f) (%f,%f) chi2/ndf %f\n", houghlines[i].sum_e, houghlines[i].ncluster, houghlines[i].start_t, houghlines[i].start_z, houghlines[i].end_t, houghlines[i].end_z, houghlines[i].chi2/houghlines[i].ncluster);
00913                 
00914                         double ym = houghlines[i].GetExpectedT(view_xmin);
00915                         double yp = houghlines[i].GetExpectedT(view_xmax);
00916                         
00917                                                 
00918                         TLine *l=new TLine(view_xmin,ym,view_xmax,yp);
00919                         
00920                         if(i==houghlines.size()-1)l->SetLineWidth(2);
00921                         l->Draw();
00922                                                 
00923                 //      printf("fit line  --%f %f - %f %f -- theta r %f %f\n",view_xmin,ym,view_xmax,yp,houghlines[i].theta,houghlines[i].r);
00924                                                 
00925         }
00926         
00927         
00928 */
00929 }

void HoughView::DrawTrueVertex ( int  view  ) 

Definition at line 502 of file HoughView.cxx.

References Munits::m, ParticleObjectHolder::mctrue, mypoh, MCTrue::vtx_u, MCTrue::vtx_v, and MCTrue::vtx_z.

Referenced by DrawEvent().

00503 {
00504 
00505         double t=0;
00506         if(view==2)
00507                 t=mypoh->mctrue.vtx_u;
00508         else if(view==3)
00509                 t=mypoh->mctrue.vtx_v;
00510         else return;
00511                                 
00512         TMarker *m = new TMarker(mypoh->mctrue.vtx_z, t, 2);
00513         m->Draw("same");
00514 
00515 }

void HoughView::DrawVertex ( int  view  ) 

Definition at line 487 of file HoughView.cxx.

References ParticleObjectHolder::event, Munits::m, mypoh, ParticleEvent::vtx_u, ParticleEvent::vtx_v, and ParticleEvent::vtx_z.

Referenced by DrawEvent().

00488 {
00489 
00490         double t=0;
00491         if(view==2)
00492                 t=mypoh->event.vtx_u;
00493         else if(view==3)
00494                 t=mypoh->event.vtx_v;
00495         else return;
00496                 
00497         TMarker *m = new TMarker(mypoh->event.vtx_z, t, 3);
00498         m->Draw("same");
00499 
00500 }

void HoughView::GetPeakAreaAverage ( double &  x,
double &  y,
double &  val,
int &  cnt,
int  curx,
int  cury,
TH2D *  hist 
) [private]

Definition at line 994 of file HoughView.cxx.

00995 {
00996 
00997         double thisval = hist->GetBinContent(curx,cury);
00998 
00999         if(thisval==0)return;
01000         
01001         val+=thisval;
01002         x+=hist->GetXaxis()->GetBinCenter(curx);
01003         y+=hist->GetYaxis()->GetBinCenter(cury);
01004         cnt++;
01005         hist->SetBinContent(curx,cury,0);
01006 
01007         
01008         GetPeakAreaAverage(x, y, val,cnt, curx+1, cury, hist);
01009         GetPeakAreaAverage(x, y, val,cnt, curx, cury+1, hist);
01010         GetPeakAreaAverage(x, y, val,cnt, curx, cury-1, hist);
01011         GetPeakAreaAverage(x, y, val,cnt, curx-1, cury, hist);
01012         GetPeakAreaAverage(x, y, val,cnt, curx+1, cury+1, hist);
01013         GetPeakAreaAverage(x, y, val,cnt, curx-1, cury-1, hist);
01014         GetPeakAreaAverage(x, y, val,cnt, curx+1, cury-1, hist);
01015         GetPeakAreaAverage(x, y, val,cnt, curx-1, cury+1, hist);
01016 
01017                 
01018 }

void HoughView::SaveHitGroup ( TH2D *  his,
TH2D *  saver,
double  save_val,
double  with_val,
int  curx,
int  cury 
) [private]

Definition at line 958 of file HoughView.cxx.

00959 {
00960         if(curx<1 || cury < 1 || curx> his->GetNbinsX() || cury >his->GetNbinsY())return;
00961         
00962 
00963         
00964         double thisval = his->GetBinContent(curx,cury);
00965 
00966 //      printf("saving %d %d thisval %f saveval %f withval %f\n",curx,cury,thisval,save_val, with_val);
00967         if(thisval==0)return;
00968         if(fabs(thisval-save_val)<0.001)
00969         {
00970                 saver->SetBinContent(curx, cury,thisval);
00971                 //thisval=0;
00972         }
00973         else if(thisval>=with_val)return;
00974 
00975 
00976         his->SetBinContent(curx,cury,0);
00977                         
00978         SaveHitGroup(his,saver,save_val,thisval,curx+1,cury);
00979         SaveHitGroup(his,saver,save_val,thisval,curx,cury+1);
00980         SaveHitGroup(his,saver,save_val,thisval,curx-1,cury);
00981         SaveHitGroup(his,saver,save_val,thisval,curx,cury-1);
00982         SaveHitGroup(his,saver,save_val,thisval,curx+1,cury+1);
00983         SaveHitGroup(his,saver,save_val,thisval,curx-1,cury-1);
00984         SaveHitGroup(his,saver,save_val,thisval,curx+1,cury-1);
00985         SaveHitGroup(his,saver,save_val,thisval,curx-1,cury+1);
00986 //      printf("save recursion done\n");
00987         
00988 
00989         //printf("saving points %d %d %f\n",curx,cury,thisval);
00990 
00991 }


Member Data Documentation

TH2D* HoughView::houghviewU [private]

Definition at line 39 of file HoughView.h.

Referenced by DrawEvent(), and HoughView().

TH2D* HoughView::houghviewV [private]

Definition at line 40 of file HoughView.h.

Referenced by DrawEvent(), and HoughView().

Definition at line 42 of file HoughView.h.

Referenced by DrawChains(), DrawClusters(), DrawEvent(), DrawHough(), DrawTrueVertex(), and DrawVertex().

TPad* HoughView::padU [private]

Definition at line 32 of file HoughView.h.

Referenced by BuildDisplay(), DrawEvent(), DrawHough(), and HoughView().

TPad* HoughView::padV [private]

Definition at line 33 of file HoughView.h.

Referenced by BuildDisplay(), DrawEvent(), DrawHough(), and HoughView().

TH2D* HoughView::viewU [private]

Definition at line 36 of file HoughView.h.

Referenced by DrawEvent(), and HoughView().

TH2D* HoughView::viewV [private]

Definition at line 37 of file HoughView.h.

Referenced by DrawEvent(), and HoughView().


The documentation for this class was generated from the following files:

Generated on 22 Nov 2017 for loon by  doxygen 1.6.1