Core clustering routine that runs on a single RPC module (two gas-gaps) Strategy:
51 {
52
68
69
70 if( col.empty() ) return false;
71 std::vector<const RpcPrepData*>::const_iterator cit_begin = col.begin();
72 std::vector<const RpcPrepData*>::const_iterator cit_end = col.end();
73 if( cit_begin == cit_end ) return false;
74 if(
debug ) std::cout <<
" RPC performing clustering: " << col.size() <<
" " <<
m_rpcIdHelper->print_to_string(subid) << std::endl;
75
76
79 std::vector<const RpcPrepData*>::const_iterator cit = cit_begin;
80 const Muon::RpcPrepData* prd_first = nullptr;
81 for( ; cit!=cit_end;++cit ) {
82 const Muon::RpcPrepData* prd = *cit;
83 const Identifier&
id = prd->
identify();
86 prd_first = prd;
87 break;
88 }
89
90
91 if( !prd_first ) return false;
92
93
100
101
102 std::vector<Doublet>* channelsPtr = nullptr;
103 cit = cit_begin;
104 for( ; cit!=cit_end;++cit ) {
105
106
107 const Muon::RpcPrepData* prd = *cit;
108 const Identifier&
id = prd->
identify();
111 if(
debug ) std::cout <<
" adding prd " <<
m_rpcIdHelper->print_to_string(
id) << std::endl;
112
113
119 if( channel >= (int)channelsPtr->size() ){
120 std::cout <<
"index channels out of range: " <<
channel <<
" max " << channelsPtr->size() << std::endl;
121 continue;
122 }
123
124
126
128
129
130 int channelClusterNumber = gasgap==1 ? doublet.first : doublet.second;
131 if( channelClusterNumber != -1 ){
133 std::cout <<
" secondary hit " <<
channel <<
" " << gasgap;
134 if( measuresPhi ) std::cout << " phi " << channelClusterNumber << std::endl;
135 else std::cout << " eta " << channelClusterNumber << std::endl;
136 }
138 if( !
cluster.addSecond(prd,gasgap) ){
139 }
140 }else{
141
142 std::vector<Id> neighbours;
144 if( gasgap == 1 ) {
145 neighbours.emplace_back(2,channel );
146 if( channel != 0 ) neighbours.emplace_back(2,channel-1 );
147 if( channel < (int)channelsPtr->size()-1 ) neighbours.emplace_back(2,channel+1 );
148 }else if( gasgap == 2 ){
149 neighbours.emplace_back(1,channel );
150 if( channel != 0 ) neighbours.emplace_back(1,channel-1 );
151 if( channel < (int)channelsPtr->size()-1 ) neighbours.emplace_back(1,channel+1 );
152 }
153 }
154 if( channel != 0 ) neighbours.emplace_back(gasgap,channel-1 );
155 if( channel < (int)channelsPtr->size()-1 ) neighbours.emplace_back(gasgap,channel+1 );
156
158 std::cout <<
" adding new channel " <<
channel <<
" " << gasgap;
159 if( measuresPhi ) std::cout << " phi " << " neighbours " << neighbours.size() << std::endl;
160 else std::cout << " eta " << " neighbours " << neighbours.size() << std::endl;
161 }
162 std::vector<Id>::iterator nit = neighbours.begin();
163 std::vector<Id>::iterator nit_end = neighbours.end();
164 RpcClusterObj* currentCluster = nullptr;
165 int currentClusterId = -1;
166 for( ;nit!=nit_end;++nit ){
167
168 Doublet& doub = (*channelsPtr)[nit->ch];
169
170 int clusterNumber = nit->gp==1 ? doub.first : doub.second;
171 if( clusterNumber == -1 ) continue;
172 if(
debug ) std::cout <<
" new neighbour " << nit->gp <<
" " << nit->ch <<
" clusterid " << clusterNumber;
173
175
176
177 if( currentCluster == nullptr ){
180 if( gasgap==1 ) doublet.first = clusterNumber;
181 else doublet.second = clusterNumber;
182 currentClusterId = clusterNumber;
183 if(
debug ) std::cout <<
" adding hit to neighbour cluster with ID " << clusterNumber << std::endl;
184 }else if( clusterNumber != currentClusterId ){
185
186
189 for( ;
h!=h_end;++
h ) {
190 const Identifier& cid = (*h)->identify();
194 if( gp==1 ) doub.
first = currentClusterId;
195 else doub.second = currentClusterId;
196 }
197 if(
debug ) std::cout <<
" found cluster overlap, merging clusters " << std::endl;
198 currentCluster->merge(
cluster);
199 }else{
200 if(
debug ) std::cout <<
" cluster overlap, same cluster " << std::endl;
201 }
202 }
203
204 if( currentCluster == nullptr ){
205 if(
debug ) std::cout <<
" no neighbouring hits, creating new cluster " <<
clusters.size() << std::endl;
212 }
213 }
214 }
215
219 std::cout <<
" cluster " <<
cl.ngasgap1 <<
" " <<
cl.ngasgap2 <<
" hits " <<
cl.hitList.size() << std::endl;
220 for(
const auto *hit :
cl.hitList ){
221 std::cout <<
" " <<
m_rpcIdHelper->print_to_string(hit->identify()) << std::endl;
222 }
223 }
224 }
225
228 }
231 }
232
233 return true;
234 }
int NphiStrips() const
returns the number of phi strips
int NetaStrips() const
returns the number of eta strips
virtual const MuonGM::RpcReadoutElement * detectorElement() const override final
Returns the detector element corresponding to this PRD.
cl
print [x.__class__ for x in toList(dqregion.getSubRegions()) ]
std::vector< RpcClusterObj > clustersPhiTmp
std::vector< Doublet > channelsEta
std::vector< Doublet > channelsPhi
std::vector< RpcClusterObj > clustersEtaTmp