/home/bes3soft/bes3soft/Boss/7.0.2/dist/7.0.2/Reconstruction/TrkReco/TrkReco-00-08-59-patch4-slc6tag/src/TBuilderCurl.cxx File Reference

#include "TrkReco/TBuilderCurl.h"
#include "TrkReco/TMLink.h"
#include "TrkReco/TTrack.h"
#include "TrackUtil/Lpav.h"
#include "TrkReco/TMLine.h"
#include "TrkReco/TRobustLineFitter.h"
#include "GaudiKernel/StatusCode.h"
#include "GaudiKernel/IInterface.h"
#include "GaudiKernel/Kernel.h"
#include "GaudiKernel/Service.h"
#include "GaudiKernel/ISvcLocator.h"
#include "GaudiKernel/SvcFactory.h"
#include "GaudiKernel/IDataProviderSvc.h"
#include "GaudiKernel/Bootstrap.h"
#include "GaudiKernel/MsgStream.h"
#include "GaudiKernel/SmartDataPtr.h"
#include "GaudiKernel/IMessageSvc.h"

Go to the source code of this file.

Defines

#define DEBUG_CURL_DUMP   0
#define DEBUG_CURL_GNUPLOT   0
#define DEBUG_CURL_MC   0

Functions

int doLineFit (AList< TMLink > &points, double &m_a, double &m_b, double &chi2, double &nhits, int ipC=1)
int checkBorder (AList< TMLink > &layer0, AList< TMLink > &layer1, AList< TMLink > &layer2)
int checkBorder (AList< TMLink > &layer0, AList< TMLink > &layer1, AList< TMLink > &layer2, AList< TMLink > &layer3)
int offsetBorder (TMLink *l)
void makeList (AList< TMLink > &layer, AList< TMLink > &list, double q, int border, int checkB, TMLink *layer0)
unsigned findMaxLocalId (unsigned layerId)
unsigned isIsolation (unsigned localId, unsigned maxLocalId, unsigned layerId, int lr, const AList< TMLink > &allStereoList)
void findTwoHits (AList< TMLink > &twoOnLayer, const AList< TMLink > &hitsOnLayer, const AList< TMLink > &allStereoList)
void setLR (AList< TMLink > &hitsOnLayer, unsigned LR=0)
bool moveLR (AList< TMLink > &hitsOnLayer)
void selectGoodWires (const AList< TMLink > &allWires, AList< TMLink > &goodWires)
void calVirtualCircle (const TMLink &hit, const TTrack &track, const int LR, HepPoint3D &center, double &radius)
void moveLR (AList< TMLink > &hits, const AList< TMLink > &hitsOnLayerOrg, const TTrack &track)

Define Documentation

#define DEBUG_CURL_DUMP   0

Definition at line 29 of file TBuilderCurl.cxx.

Referenced by TBuilderCurl::buildStereo(), TCurlFinder::checkExceptionalSegmentsType01(), TCurlFinder::checkExceptionalSegmentsType02(), TCurlFinder::checkExceptionalSegmentsType03(), TCurlFinder::make3DTrack(), and TCurlFinder::makeCurlTracks().

#define DEBUG_CURL_GNUPLOT   0

Definition at line 30 of file TBuilderCurl.cxx.

#define DEBUG_CURL_MC   0

Definition at line 31 of file TBuilderCurl.cxx.

Function Documentation

void calVirtualCircle ( const TMLink &  hit,
const TTrack &  track,
const int  LR,
HepPoint3D &  center,
double &  radius 
)

Definition at line 1472 of file TBuilderCurl.cxx.

References abs, Helix::center(), TTrack::charge(), TMDCWireHit::drift(), TTrack::helix(), TMLink::hit(), TTrack::radius(), unit, TMLink::wire(), and TMDCWire::xyPosition().

Referenced by TBuilderCurl::makeLine(), and moveLR().

01473                                                     {
01474   if(abs(LR) != 1)return;
01475   double Q = track.charge();
01476   int isOuter = 1;
01477   if(Q > 0. && LR == 1)isOuter = -1; // Inner
01478   else if(Q < 0. && LR == -1)isOuter = -1; // Inner
01479   radius = fabs(track.radius());
01480   center = track.helix().center();    
01481   HepPoint3D wire = hit.wire()->xyPosition();
01482   center.setZ(0.);
01483   wire.setZ(0.);
01484   // new center(virtual)
01485   center = wire + 
01486     (radius+((double)isOuter)*(hit.hit()->drift()))*(center-wire).unit();
01487 }

int checkBorder ( AList< TMLink > &  layer0,
AList< TMLink > &  layer1,
AList< TMLink > &  layer2,
AList< TMLink > &  layer3 
)

Definition at line 925 of file TBuilderCurl.cxx.

References genRecEmupikp::i, and delete_small_size::size.

00928                                   {
00929   AList<TMLink> list = layer0;
00930   list.append(layer1);
00931   list.append(layer2);
00932   list.append(layer3);
00933   int size = list.length();
00934   if(size <= 1)return 0;
00935   int layerId = list[0]->hit()->wire()->layerId();
00936   int maxLocalId = 79;
00937   if(layerId >= 15)maxLocalId = 127;
00938   if(layerId >= 23)maxLocalId = 159;
00939   if(layerId >= 32)maxLocalId = 207;
00940   if(layerId >= 41)maxLocalId = 255;
00941   int HId = (int)(0.8*(maxLocalId+1));
00942   int LId = (int)(0.2*(maxLocalId+1));
00943   int low  = 0;
00944   int high = 0;
00945   for(int i=0;i>size;++i){
00946     if(list[i]->hit()->wire()->localId() < LId)low = 1;
00947     else if(list[i]->hit()->wire()->localId() > HId)high = 1;
00948     if(low == 1 && high == 1)return 1;
00949   }
00950   return 0;
00951 }

int checkBorder ( AList< TMLink > &  layer0,
AList< TMLink > &  layer1,
AList< TMLink > &  layer2 
)

Definition at line 898 of file TBuilderCurl.cxx.

References genRecEmupikp::i, and delete_small_size::size.

00900                                   {
00901   AList<TMLink> list = layer0;
00902   list.append(layer1);
00903   list.append(layer2);
00904   int size = list.length();
00905   if(size <= 1)return 0;
00906   int layerId = list[0]->hit()->wire()->layerId();
00907   int maxLocalId = 79;
00908   if(layerId >= 15)maxLocalId = 127;
00909   if(layerId >= 23)maxLocalId = 159;
00910   if(layerId >= 32)maxLocalId = 207;
00911   if(layerId >= 41)maxLocalId = 255;
00912   int HId = (int)(0.8*(maxLocalId+1));
00913   int LId = (int)(0.2*(maxLocalId+1));
00914   int low  = 0;
00915   int high = 0;
00916   for(int i=0;i>size;++i){
00917     if(list[i]->hit()->wire()->localId() < LId)low = 1;
00918     else if(list[i]->hit()->wire()->localId() > HId)high = 1;
00919     if(low == 1 && high == 1)return 1;
00920   }
00921   return 0;
00922 }

int doLineFit ( AList< TMLink > &  points,
double &  m_a,
double &  m_b,
double &  chi2,
double &  nhits,
int  ipC = 1 
)

Definition at line 783 of file TBuilderCurl.cxx.

References genRecEmupikp::i, and x.

Referenced by TBuilderCurl::fitLine().

00786 {
00787     m_a = m_b = nhits = 0.;
00788     chi2 = 1.e+10;
00789     unsigned n = points.length();
00790     double sum = double(n);
00791     double sumX = 0., sumY = 0., sumX2 = 0., sumXY = 0., sumY2 = 0.;
00792     for (unsigned i = 0; i < n; i++) {
00793       TMLink & l = * points[i];
00794       
00795       double x = l.position().x();
00796       double y = l.position().y();
00797       sumX  += x;
00798       sumY  += y;
00799       sumX2 += x * x;
00800       sumXY += x * y;
00801       sumY2 += y * y;
00802     }
00803     if(ipC != 0)sum += 1.0;// IP Constraint
00804     nhits = sum;
00805     
00806     double m_det = sum * sumX2 - sumX * sumX;
00807     if(m_det == 0. && sum != 2.){
00808       return -1;
00809     }else if(m_det == 0. && sum == 2.){
00810       double x0 = points[0]->position().x();
00811       double y0 = points[0]->position().y();
00812       double x1 = points[1]->position().x();
00813       double y1 = points[1]->position().y();
00814       if(x0 == x1)return -1;
00815       m_a = (y0-y1)/(x0-x1);
00816       m_b = -m_a*x1 + y1;
00817       chi2 = 0.;
00818       return 0;
00819     }
00820     chi2 = 0.;
00821     m_a = (sumXY * sum - sumX * sumY) / m_det;
00822     m_b = (sumX2 * sumY - sumX * sumXY) / m_det;
00823 
00824     for(unsigned i = 0; i < n; i++) {
00825       TMLink & l = * points[i];
00826       
00827       double x = l.position().x();
00828       double y = l.position().y();
00829       double d =  y-(x*m_a+m_b);
00830       chi2 += d*d;
00831     }
00832 
00833     return 0;
00834 }

unsigned findMaxLocalId ( unsigned  layerId  ) 

Definition at line 1261 of file TBuilderCurl.cxx.

Referenced by findTwoHits().

01261                                 {
01262 /*  unsigned maxLocalId = 79;
01263   if(superLayerId == 3)maxLocalId = 127;
01264   else if(superLayerId == 5)maxLocalId = 159;
01265   else if(superLayerId == 7)maxLocalId = 207;
01266   else if(superLayerId == 9)maxLocalId = 255;
01267   return maxLocalId;
01268 */
01269 //changed to BESiii, Liuqg
01270   unsigned maxLocalId = 39;
01271   if(layerId == 1)maxLocalId = 43;
01272   else if(layerId == 2)maxLocalId = 47;
01273   else if(layerId == 3)maxLocalId = 55;
01274   else if(layerId == 4)maxLocalId = 63;
01275   else if(layerId == 5)maxLocalId = 71;
01276   else if(layerId == 6 || layerId == 7 )maxLocalId = 79;
01277   else if(layerId == 20 || layerId == 21 || layerId == 22 || layerId == 23)maxLocalId = 159;
01278   else if(layerId == 24 || layerId == 25 || layerId == 26 || layerId == 27)maxLocalId = 175;
01279   else if(layerId == 28 || layerId == 29 || layerId == 30 || layerId == 31)maxLocalId = 207;
01280   else if(layerId == 32 || layerId == 33 || layerId == 34 || layerId == 35)maxLocalId = 239;
01281 
01282   return maxLocalId;
01283 }

void findTwoHits ( AList< TMLink > &  twoOnLayer,
const AList< TMLink > &  hitsOnLayer,
const AList< TMLink > &  allStereoList 
)

Definition at line 1313 of file TBuilderCurl.cxx.

References findMaxLocalId(), isIsolation(), and rb::R().

Referenced by TBuilderCurl::makeLine().

01315                                                {
01316   //...finds "two" seq. and isolated hits.
01317   //...and then sets LR for selected hits.
01318   if(hitsOnLayer.length() == 0 ||
01319      hitsOnLayer.length() > 3)return;
01320   twoOnLayer.removeAll();
01321   if(hitsOnLayer.length() == 1){
01322     if(hitsOnLayer[0]->wire()->superLayerId() == 0)return; //origin is == 1, Liuqg 060921
01323     unsigned maxLocalId = findMaxLocalId(hitsOnLayer[0]->wire()->layerId());
01324     unsigned R = isIsolation(hitsOnLayer[0]->wire()->localId(),
01325                              maxLocalId,
01326                              hitsOnLayer[0]->wire()->layerId(),1,allStereoList);
01327     unsigned L = isIsolation(hitsOnLayer[0]->wire()->localId(),
01328                              maxLocalId,
01329                              hitsOnLayer[0]->wire()->layerId(),-1,allStereoList);
01330     if(R == 1 && L == 0){
01331       unsigned nextLocalId = hitsOnLayer[0]->wire()->localIdForPlus()+1;
01332       L = isIsolation(nextLocalId,
01333                       maxLocalId,
01334                       hitsOnLayer[0]->wire()->layerId(),-1,allStereoList);
01335       if(L == 1){ // xuox
01336         hitsOnLayer[0]->leftRight(1); // R
01337         hitsOnLayer[0]->position(hitsOnLayer[0]->arcZ(1));
01338         twoOnLayer.append(hitsOnLayer[0]);
01339       }
01340     }else if(R == 0 && L == 1){
01341       unsigned nextLocalId = hitsOnLayer[0]->wire()->localIdForMinus()+1;
01342       R = isIsolation(nextLocalId,
01343                       maxLocalId,
01344                       hitsOnLayer[0]->wire()->layerId(),1,allStereoList);
01345       if(R == 1){ // xoux
01346         hitsOnLayer[0]->leftRight(0); // L
01347         hitsOnLayer[0]->position(hitsOnLayer[0]->arcZ(0));
01348         twoOnLayer.append(hitsOnLayer[0]);
01349       }
01350     }   
01351   }
01352   if(hitsOnLayer.length() == 2){
01353     if(hitsOnLayer[0]->wire()->localIdForPlus()+1 ==
01354        hitsOnLayer[1]->wire()->localId()){
01355       unsigned maxLocalId = findMaxLocalId(hitsOnLayer[0]->wire()->layerId());
01356       unsigned R = isIsolation(hitsOnLayer[0]->wire()->localId(),
01357                                maxLocalId,
01358                                hitsOnLayer[0]->wire()->layerId(),1,allStereoList);
01359       unsigned L = isIsolation(hitsOnLayer[1]->wire()->localId(),
01360                                maxLocalId,
01361                                hitsOnLayer[1]->wire()->layerId(),-1,allStereoList);
01362       if(R == 1 && L == 1){ // xoox
01363         hitsOnLayer[0]->leftRight(1); // R
01364         hitsOnLayer[0]->position(hitsOnLayer[0]->arcZ(1));
01365         hitsOnLayer[1]->leftRight(0); // L
01366         hitsOnLayer[1]->position(hitsOnLayer[1]->arcZ(0));
01367         twoOnLayer.append(hitsOnLayer[0]);
01368         twoOnLayer.append(hitsOnLayer[1]);
01369       }
01370     }   
01371   }
01372   if(hitsOnLayer.length() == 3){
01373     if(hitsOnLayer[0]->wire()->localIdForPlus()+1 ==
01374        hitsOnLayer[1]->wire()->localId() &&
01375        hitsOnLayer[1]->wire()->localIdForPlus()+1 !=
01376        hitsOnLayer[2]->wire()->localId()){
01377       unsigned maxLocalId = findMaxLocalId(hitsOnLayer[0]->wire()->layerId());
01378       unsigned R = isIsolation(hitsOnLayer[0]->wire()->localId(),
01379                                maxLocalId,
01380                                hitsOnLayer[0]->wire()->layerId(),1,allStereoList);
01381       unsigned L = isIsolation(hitsOnLayer[1]->wire()->localId(),
01382                                maxLocalId,
01383                                hitsOnLayer[1]->wire()->layerId(),-1,allStereoList);
01384       if(R == 1 && L == 1){ // oxoox
01385         hitsOnLayer[0]->leftRight(1); // R
01386         hitsOnLayer[0]->position(hitsOnLayer[0]->arcZ(1));
01387         hitsOnLayer[1]->leftRight(0); // L
01388         hitsOnLayer[1]->position(hitsOnLayer[1]->arcZ(0));
01389         twoOnLayer.append(hitsOnLayer[0]);
01390         twoOnLayer.append(hitsOnLayer[1]);
01391       }
01392     }else if(hitsOnLayer[0]->wire()->localIdForPlus()+1 !=
01393              hitsOnLayer[1]->wire()->localId() &&
01394              hitsOnLayer[1]->wire()->localIdForPlus()+1 ==
01395              hitsOnLayer[2]->wire()->localId()){
01396       unsigned maxLocalId = findMaxLocalId(hitsOnLayer[1]->wire()->layerId());
01397       unsigned R = isIsolation(hitsOnLayer[1]->wire()->localId(),
01398                                maxLocalId,
01399                                hitsOnLayer[1]->wire()->layerId(),1,allStereoList);
01400       unsigned L = isIsolation(hitsOnLayer[2]->wire()->localId(),
01401                                maxLocalId,
01402                                hitsOnLayer[2]->wire()->layerId(),-1,allStereoList);
01403       if(R == 1 && L == 1){ // xooxo
01404         hitsOnLayer[1]->leftRight(1); // R
01405         hitsOnLayer[1]->position(hitsOnLayer[1]->arcZ(1));
01406         hitsOnLayer[2]->leftRight(0); // L
01407         hitsOnLayer[2]->position(hitsOnLayer[2]->arcZ(0));
01408         twoOnLayer.append(hitsOnLayer[1]);
01409         twoOnLayer.append(hitsOnLayer[2]);
01410       }
01411     }
01412   }
01413   /* if(twoOnLayer.length() != 0){
01414     std::cout << "TWO " << twoOnLayer.length() << std::endl;
01415     } */
01416 }

unsigned isIsolation ( unsigned  localId,
unsigned  maxLocalId,
unsigned  layerId,
int  lr,
const AList< TMLink > &  allStereoList 
)

Definition at line 1286 of file TBuilderCurl.cxx.

References genRecEmupikp::i.

Referenced by findTwoHits().

01290                                                  {
01291   unsigned findId;
01292   if(lr == 1){ // R : ox, Dose a wire exist at "x"?
01293     findId = maxLocalId;
01294     if(localId != 0)findId = localId-1;
01295   }else if(lr == -1){ // L : xo, Dose a wire exist at "x"?
01296     findId = 0;
01297     if(localId != maxLocalId)findId = localId+1;
01298   }else{
01299     return 1;
01300   }
01301   unsigned isolate = 1;
01302   for(int i=0;i<allStereoList.length();++i){
01303     if(allStereoList[i]->wire()->layerId() == layerId &&
01304        allStereoList[i]->wire()->localId() == findId){
01305       isolate = 0;
01306       break;     
01307     }
01308   }
01309   return isolate;
01310 }

void makeList ( AList< TMLink > &  layer,
AList< TMLink > &  list,
double  q,
int  border,
int  checkB,
TMLink *  layer0 
)

Definition at line 965 of file TBuilderCurl.cxx.

References genRecEmupikp::i, and offsetBorder().

00965                                                                                                      {
00966   int n = layer.length();
00967   if(checkB == 0){
00968     for(int i=0;i<n;++i){
00969       if(q < 0.){
00970         if(layer[i]->hit()->wire()->localId() >= border)list.append(layer[i]);
00971       }else{
00972         if(layer[i]->hit()->wire()->localId() <= border)list.append(layer[i]);
00973       }
00974     }
00975   }else{
00976     //difficult!! --> puls offset
00977     int offset = offsetBorder(layer0);
00978     if(border*2 < offset)border += offset;
00979     for(int i=0;i<n;++i){
00980       int lId = layer[i]->hit()->wire()->localId();
00981       if(lId*2 < offset)lId += offset;
00982       if(q < 0.){
00983         if(lId >= border)list.append(layer[i]);
00984       }else{
00985         if(lId <= border)list.append(layer[i]);
00986       }
00987     }
00988   }
00989 }

void moveLR ( AList< TMLink > &  hits,
const AList< TMLink > &  hitsOnLayerOrg,
const TTrack &  track 
)

Definition at line 1490 of file TBuilderCurl.cxx.

References calVirtualCircle(), TTrack::charge(), and genRecEmupikp::i.

01492                            {
01493   AList<TMLink> hitsOnLayer = hitsOnLayerOrg;
01494   hitsOnLayer.remove(hits);
01495   if(hitsOnLayer.length() == 0)return;
01496 
01497   unsigned nHits = hits.length();
01498   if(nHits == 0)return;
01499 
01500   //...finds "ref" from hits.
01501   // ex) LLLL|, LLL|R, LL|RR, L|RRR, |RRRR
01502   int LR = -1; // L
01503   TMLink ref;
01504   if(hits[nHits-1]->leftRight() == 1){ // All R
01505     LR = 1; // R
01506     ref = *hits[nHits-1];
01507   }
01508   for(unsigned i=0;i<nHits;++i){
01509     if(hits[i]->leftRight() == 0){ // L
01510       ref = *hits[i];
01511       break;
01512     }
01513   }
01514 
01515   HepPoint3D center;
01516   double radius;
01517   calVirtualCircle(ref,track,LR,center,radius);
01518 
01519   double Q = track.charge();
01520   for(int i=0;i<hitsOnLayer.length();++i){
01521     int isOuter = 1;
01522     if((hitsOnLayer[i]->wire()->xyPosition()-center).mag()-radius < 0.)isOuter = -1;
01523     if(Q > 0. && isOuter == 1){
01524       hitsOnLayer[i]->position(hitsOnLayer[i]->arcZ(0)); // L
01525       hitsOnLayer[i]->leftRight(0); // L
01526     }else if(Q > 0. && isOuter == -1){
01527       hitsOnLayer[i]->position(hitsOnLayer[i]->arcZ(1)); // R
01528       hitsOnLayer[i]->leftRight(1); // R
01529     }else if(Q < 0. && isOuter == 1){
01530       hitsOnLayer[i]->position(hitsOnLayer[i]->arcZ(1)); // R
01531       hitsOnLayer[i]->leftRight(1); // R
01532     }else if(Q < 0. && isOuter == -1){
01533       hitsOnLayer[i]->position(hitsOnLayer[i]->arcZ(0)); // L
01534       hitsOnLayer[i]->leftRight(0); // L
01535     }
01536   }
01537 }

bool moveLR ( AList< TMLink > &  hitsOnLayer  ) 

Definition at line 1434 of file TBuilderCurl.cxx.

References genRecEmupikp::i.

Referenced by TBuilderCurl::makeLine().

01434                                   {
01435   unsigned nHits = hitsOnLayer.length();
01436   if(nHits == 0)return false;
01437   // ex) LLLL --> LLLR --> LLRR --> LRRR --> RRRR
01438   if(hitsOnLayer[nHits-1]->leftRight() == 1)return false; // All R
01439   for(unsigned i=0;i<nHits;++i){
01440     if(hitsOnLayer[i]->leftRight() == 0){ // L
01441       hitsOnLayer[i]->leftRight(1); // R
01442       hitsOnLayer[i]->position(hitsOnLayer[i]->arcZ(1));
01443       return true;
01444     }
01445   }
01446   return false;
01447 }

int offsetBorder ( TMLink *  l  ) 

Definition at line 954 of file TBuilderCurl.cxx.

References TMLink::hit(), TMDCWire::layerId(), and TMDCWireHit::wire().

Referenced by makeList().

00954                        {
00955   int layerId = l->hit()->wire()->layerId();
00956   int maxLocalId = 79;
00957   if(layerId >= 15)maxLocalId = 127;
00958   if(layerId >= 23)maxLocalId = 159;
00959   if(layerId >= 32)maxLocalId = 207;
00960   if(layerId >= 41)maxLocalId = 255;
00961   return maxLocalId+1;
00962 }

void selectGoodWires ( const AList< TMLink > &  allWires,
AList< TMLink > &  goodWires 
)

Definition at line 1450 of file TBuilderCurl.cxx.

References genRecEmupikp::i, and x.

Referenced by TBuilderCurl::fitLine2(), and TBuilderCurl::makeLine().

01451                                          {
01452   goodWires.removeAll();
01453   for(int i=0;i<allWires.length();++i){
01454     if(allWires[i]->position().x() != -999.){
01455       goodWires.append(allWires[i]);
01456     }
01457   }
01458 }

void setLR ( AList< TMLink > &  hitsOnLayer,
unsigned  LR = 0 
)

Definition at line 1419 of file TBuilderCurl.cxx.

References genRecEmupikp::i.

Referenced by TBuilderCurl::makeLine().

01419                                                   {
01420   // LR = 0 : L
01421   //    = 1 : R
01422   for(unsigned i=0;i<hitsOnLayer.length();++i){
01423     if(LR == 0){
01424       hitsOnLayer[i]->leftRight(0); // L
01425       hitsOnLayer[i]->position(hitsOnLayer[i]->arcZ(0));
01426     }else{
01427       hitsOnLayer[i]->leftRight(1); // R
01428       hitsOnLayer[i]->position(hitsOnLayer[i]->arcZ(1));
01429     }
01430   }
01431 }

Generated on Tue Nov 29 23:17:06 2016 for BOSS_7.0.2 by  doxygen 1.4.7