00001 #include <iostream> 00002 #include "treesearch.h" 00003 #include <math.h> 00004 using namespace std; 00005 00006 //__________________________________ 00007 static int _hashgen = 1; 00008 static int hashgen(void) { 00009 _hashgen += 2; 00010 _hashgen &= 0x7ffffff; 00011 return _hashgen; 00012 } 00013 extern int tlayers; 00014 extern int trelinanz; 00015 extern TreeLine *trelin; 00016 extern Options opt; 00017 //__________________________________ 00018 static int has_hits[TLAYERS]; 00019 //__________________________________ 00102 treesearch::treesearch(){ 00103 tlayers = TLAYERS; 00104 } 00105 treesearch::~treesearch(){ 00106 00107 } 00108 /*---------------------------------------------------------------------------*\ 00109 00110 wireselection() - this function steps through the hits from the unprimed 00111 and primed planes for a tree-plane to decide whether 00112 hits should or should not be paired together when the 00113 hit pattern for the tree-plane is constructed. 00114 00115 inputs: (1) Hit **x - pointer to the hit in the linked list of 00116 hits for the unprimed plane at which the 00117 scan is to start. 00118 (2) Hit **X - pointer to the hit in the linked list of 00119 hits for the primed plane at which the 00120 scan is to start. 00121 (3) double maxdist - maximum separation between hits on the 00122 primed and unprimed hits for these hits 00123 to paired together when making the bit 00124 pattern is formed. 00125 00126 outputs: (1) Hit **x - pointer to the hit in the linked list of 00127 for the unprimed plane which should be 00128 used for generating the bit pattern; =0 if 00129 the hit on the unprimed plane should not 00130 be used at all. 00131 (2) Hit **X - pointer to the hit in the linked list of 00132 for the primed plane which should be 00133 used for generating the bit pattern; =0 if 00134 the hit on the primed plane should not 00135 be used at all. 00136 (3) Hit **xn - pointer to the next hit to consider for 00137 the unprimed plane 00138 (4) Hit **Xn - pointer to the next hit to consider for 00139 the primed plane 00140 00141 \*---------------------------------------------------------------------------*/ 00142 00143 void treesearch::wireselection( Hit **x, Hit **X, Hit **xn, Hit**Xn, double maxdist) 00144 { 00145 cerr << "THIS FUNCTION NEEDS REVISION" << endl; 00146 double wireDistance1; 00147 double wireDistance2; 00148 00149 if( *x && *X) { /* there are hits on both primed and unprimed planes */ 00150 00151 /* ---- compute the distance between the hits on the two planes ---- */ 00152 00153 wireDistance2 = (*x)->rPos2 - (*X)->rPos1; /* unprimed left-primed right */ 00154 wireDistance1 = (*x)->rPos1 - (*X)->rPos2; /* unprimed right-primed left */ 00155 00156 if( wireDistance2 < -maxdist) { 00157 00158 /* ---- CASE 1: the unprimed hit cannot be paired with any of the 00159 primed hits. So, it needs to be considered by 00160 itself as the bit pattern is constructed. Also, 00161 the next scan should begin with the next hit on 00162 this unprimed plane and see if it can be paired 00163 with the current hit on the primed plane. ---- */ 00164 00165 *Xn = *X; /* keep current point as starting point for 00166 the primed hits in the next scan */ 00167 *xn = (*x)->nextdet; /* use next hit as the starting point for the 00168 unprimed hits in the next scan */ 00169 *X = 0; /* no primed match, so return the unprimed hit 00170 for use as the bit pattern is generated. */ 00171 } else if( wireDistance1 > maxdist) { 00172 00173 /* ---- CASE 2: the primed hit cannot be paired with any of the 00174 unprimed hits. So, it needs to be considered by 00175 itself as the bit pattern is constructed. Also, 00176 the next scan should begin with the next hit on 00177 this primed plane and see if it can be paired 00178 with the current hit on the unprimed plane. ---- */ 00179 00180 *xn = *x; /* keep current point as starting point for 00181 the unprimed hits in the next scan */ 00182 *Xn = (*X)->nextdet; /* use next hit as the starting point for the 00183 primed hits in the next scan */ 00184 *x = 0; /* no unprimed match, so return the primed hit 00185 for use as the bit pattern is generated. */ 00186 } else { 00187 00188 /* ---- CASE 3: the primed hit and the unprimed hit are paired. 00189 So, both need to be considered as the hit pattern 00190 is constructed. Also, the next scan should begin 00191 with the next hit on each of these planes. ---- */ 00192 00193 *xn = (*x)->nextdet; /* use next hit as the starting point for the 00194 unprimed hits in the next scan */ 00195 *Xn = (*X)->nextdet; /* use next hit as the starting point for the 00196 primed hits in the next scan */ 00197 } 00198 } else if( *x ) { /* only hits on primed plane */ 00199 *Xn = 0; /* no unprimed plane for the next scan */ 00200 *xn = (*x)->nextdet; /* use next hit as the starting point for the 00201 primed hits in the next scan */ 00202 } else if( *X ) { /* only hits on unprimed plane */ 00203 *xn = 0; /* no primed plane for the next scan */ 00204 *Xn = (*X)->nextdet; /* use next hit as the starting point for the 00205 primed hits in the next scan */ 00206 } else /* no hits on either plane */ 00207 *xn = *Xn = 0; /* no next scan! */ 00208 00209 return; 00210 } 00211 //________________________________________________________________________ 00212 /*---------------------------------------------------------------------------*\ 00213 00214 _setpoints() - this function sets the bins in a hit pattern for a 00215 range of positions. The range of hit patterns is specified 00216 by a start and a stop position in the detector. This 00217 function turns on the bins in the hit pattern for each 00218 level of the bin-division used in the treesearch algorithm. 00219 00220 inputs: (1) double posStart - position (in cm) of the start of the 00221 range for which bins are to be turned 00222 on 00223 (2) double posEnd - position (in cm) of the stop of the 00224 range for which bins are to be turned 00225 on 00226 (3) double detectorwidth - width of the tree-detector (in cm) 00227 (4) unsigned binwidth - width of a bin at the deepest level 00228 of bin-division in the treesearch. 00229 (5) char *pattern - pointer to the hit pattern for this 00230 tree-detector 00231 (6) int *hash - pointer to ??? 00232 00233 outputs: (1) char *pattern - pointer to the hit pattern for this 00234 tree-detector 00235 (2) int *hash - pointer to ??? 00236 00237 \*---------------------------------------------------------------------------*/ 00238 00239 void treesearch::_setpoints(double posStart,double posEnd,double detectorwidth, 00240 unsigned binwidth,char *pattern,int *hash) 00241 { 00242 int i,j; 00243 int ia, ie, hashint = hashgen(); 00244 unsigned oldwidth = binwidth; 00245 char *oldpattern = pattern; 00246 00247 /* ---- compute the first bin in the deepest tree level to turn on ---- */ 00248 00249 ia = (int)floor(posStart/detectorwidth * binwidth); 00250 00251 /* ---- compute the last bin in the deepest tree level to turn on ---- */ 00252 00253 ie = (int)floor(posEnd /detectorwidth * binwidth); 00254 00255 /* ---- step through each of the bins at the deepest bin-division 00256 level in the hit pattern and turn on the bits in the 00257 pattern at all the bin-division levels for this bin. ---- */ 00258 //cerr << "(" << ia << "," << ie << "," << posStart<< "," << posEnd << "," << detectorwidth<< "," << binwidth << ")" << endl; 00259 00260 for( j = ia; j <= ie; j ++ ) { /* loop over the bins to be set */ 00261 00262 i = j; /* remember the bin at the deepest level which is being turn on */ 00263 00264 binwidth = oldwidth; /* the size of the bit pattern for the detector 00265 at the deepest level of bin-division. */ 00266 pattern = oldpattern; /* pointer to start of the bit pattern */ 00267 00268 /* ---- check if the bin is inside the detector ---- */ 00269 //I added "(signed int)" to the 00270 //following line to prevent the warning from comparing signed and 00271 //unsigned integers 00272 if( i >= (signed int)binwidth ) 00273 return; 00274 if( i < 0 ) 00275 continue; 00276 /* ---- compute some hash value ??? ---- */ 00277 hash[i] = ((hash[i]<<1) + hashint)|1 ; 00278 /* ---- set the bits in the bit pattern at all depths of bin-division ---- */ 00279 00280 if( pattern ) { 00281 while( binwidth ) { /* starting at maximum depth, loop over 00282 each depth of the treesearch */ 00283 pattern[i] = 1; /* turn on the bit in this bin */ 00284 pattern += binwidth; /* set ahead to the part of the bit 00285 pattern in which the bits for the 00286 next higher bin-division are stored */ 00287 i >>= 1; /* go up one depth of the depth */ 00288 binwidth >>= 1; /* size of the bit pattern for the next 00289 higher bin-division is half of the 00290 size of the bit pattern for a given 00291 level of bin-division. */ 00292 00293 } /* end of loop over the depths of the treesearch */ 00294 } 00295 } /* end of loop over the bins to be set */ 00296 } 00297 //________________________________________________________________________ 00298 /*---------------------------------------------------------------------------*\ 00299 00300 _setpoint() - this function sets the bins in the hit pattern for a 00301 range of positions around a central point within a specified 00302 distance/resolution by calling the setpoint() function 00303 described above. 00304 00305 inputs: (1) double position - position (in cm) of the center point 00306 around which the bins are to be turned 00307 on 00308 (2) double resolution - distance (in cm) around the central 00309 point for which bins are to be turned 00310 on 00311 (3) double detectorwidth - width of the tree-detector (in cm) 00312 (4) unsigned binwidth - width of a bin at the deepest level 00313 of bin-division in the treesearch. 00314 (5) char *pattern - pointer to the hit pattern for this 00315 tree-detector 00316 (6) int *hash - ??? 00317 00318 outputs: (1) char *pattern - pointer to the hit pattern for this 00319 tree-detector 00320 (2) int *hash - pointer to ??? 00321 00322 \*---------------------------------------------------------------------------*/ 00323 00324 void treesearch::_setpoint(double position, double resolution,double detectorwidth, 00325 unsigned binwidth, char *pattern, int *hash) 00326 { 00327 _setpoints(position-resolution,position+resolution, 00328 detectorwidth,binwidth,pattern,hash); 00329 } 00330 //________________________________________________________________________ 00331 /*---------------------------------------------------------------------------*\ 00332 00333 setpoint() - this function sets the bins in the hit pattern for a range 00334 of positions specified by a center point and a half-distance 00335 around the center point by calling the _setpoint() function 00336 described above. 00337 00338 If hit patterns are supplied to this function for both the 00339 primed and unprimed planes, the two planes are treated as 00340 separate tree-detectors. In this case, the hit patterns 00341 are updated for each plane individually. This feature is 00342 used for the front region when the VCs are not included 00343 in the fit because, for these fits, each FC plane is treated 00344 as a single tree-detector. 00345 00346 Otherwise, the two planes are treated as a single 00347 tree-detector and a single hit pattern is formed from the 00348 hits seen in the two planes. 00349 00350 inputs: (01) double off - alignment offset of the 00351 tree-detector (in cm) 00352 (02) double h1 - center point of hit range on 00353 plane 1 (in cm) 00354 (03) double res1 - width of hit range on plane 1 (in cm) 00355 (04) double h2 - center point of hit range on 00356 plane 2 (in cm) 00357 (05) double res2 - width of hit range on plane 2 (in cm) 00358 (06) double width - width of the tree-detector (in cm) 00359 (07) unsigned binwidth - width of a bin at the deepest level 00360 of bin-division in the treesearch. 00361 (08) char *pa - pointer to the hit pattern for 00362 tree-detector 1 00363 (09) char *pb - pointer to the hit pattern for 00364 tree-detector 2, =0 if there is no 00365 tree-detector 2 because the planes are 00366 being combined into one tree-detector 00367 (10) int *hasha - pointer to ??? 00368 (11) int *hashb - pointer to ??? 00369 00370 outputs: (01) char *pa - pointer to the updated hit pattern 00371 for tree-detector 1 00372 (02) char *pb - pointer to the updated hit pattern for 00373 tree-detector 2, =0 if there is no 00374 tree-detector 2 because the planes are 00375 being combined into one tree-detector 00376 (03) int *hasha - pointer to ??? 00377 (04) int *hashb - pointer to ??? 00378 00379 \*---------------------------------------------------------------------------*/ 00380 00381 void treesearch::setpoint( double off,double h1, double res1,double h2, double res2, 00382 double width, unsigned bwidth, char *pa, char *pb, 00383 int *hasha, int *hashb) 00384 { 00385 00386 00387 if( pb ) { /* NOT combining two planes into one tree-detector, so 00388 set bins in the hit pattern for each plane separately */ 00389 _setpoint( off+h1, res1, width, bwidth, pa, hasha ); 00390 _setpoint( off+h2, res2, width, bwidth, pb, hashb ); 00391 } else /* Combining two planes into one tree-detector, so just 00392 set the bins for the overlap between the hit on the 00393 two planes. Here, the resolutions are added linearly 00394 instead of in quadrature to save cputime */ 00395 _setpoint( off+0.5*(h1+h2), 0.5*(res1+res2), width, bwidth, pa, hasha ); 00396 } 00397 //________________________________________________________________________ 00398 /*---------------------------------------------------------------------------*\ 00399 00400 TsSetPoint() - this function creates the bit patterns for two partner 00401 planes (planes with the like-pitched wires in the 00402 same chamber. 00403 00404 This function can treat the two planes in two different 00405 fashions, specifically: 00406 00407 (1) the planes can be treated as one tree-detector. In 00408 this case, a check is performed to see if a hit on 00409 one of the planes can be matched to a hit on the 00410 other plane (within the allowances of the maximum 00411 track slope specified by the HRCset parameter 00412 MaxSlope). If the hits are paired together by this 00413 search, the bit pattern is set around the allowable 00414 tracks through both hits. If a hit is not paired 00415 in this search, the bit pattern is set for the region 00416 around the hit. (This method is the standard for 00417 the back chambers since they have eight planes per 00418 treeline). 00419 00420 (2) the planes are treated as individual tree-detectors. 00421 In this case, a separate hit pattern is created for 00422 both planes and the above described searching for 00423 paired hits is not employed. (This method is the 00424 standard for the front chambers since they have only 00425 four planes per treeline). 00426 00427 inputs: (01) double detectorwidth - width of the tree-detector (in cm) 00428 (02) double zdistance - distance between the two planes 00429 forming this tree-detector (in cm) 00430 (03) Hit *Ha - linked list of hits for the unprimed 00431 plane of this tree-detector 00432 (04) Hit *Hb - linked list of hits for the primed 00433 plane of this tree-detector 00434 (05) unsigned binwidth - width of a bin (in cm) at the deepest 00435 bin-division in the treesearch. 00436 00437 outputs: (01) int TsSetPoint() - the number of single and paired hits 00438 processed by this return 00439 (02) char *patterna - pointer to the hit pattern for 00440 tree-detector 00441 (03) char *patternb - pointer to the hit pattern for 00442 the optional second tree-detector 00443 (04) int *hasha - pointer to ??? 00444 (05) int *hashb - pointer to ??? 00445 00446 \*---------------------------------------------------------------------------*/ 00447 // This TsSetPoint version is designed for setting a pattern for 1 hit at a time 00448 int treesearch::TsSetPoint(double detectorwidth, Hit *H,char *pattern, int *hash, unsigned binwidth) 00449 { 00450 double dw2 = (detectorwidth/2.0); /* half-width of the detector (in cm) */ 00451 _setpoints(dw2 - H->rPos1 - H->rTrackResolution,//set the points on the front side of the wire 00452 dw2 - H->rPos1 + H->rTrackResolution, 00453 detectorwidth,binwidth,pattern,hash); 00454 _setpoints(dw2 + H->rPos1 - H->rTrackResolution,//set the points on the back side of the wire 00455 dw2 + H->rPos1 + H->rTrackResolution, 00456 detectorwidth,binwidth,pattern,hash); 00457 00458 //for(int i=0;i<8;i++)cerr << hash[i] << " ";cerr << endl; 00459 /* 00460 for(int i=0;i<15;i++){ 00461 cerr << "<" << pattern[i] << ">" ; 00462 } 00463 cerr << endl; 00464 */ 00465 return 1; 00466 00467 00468 00469 } 00470 /*---------------------------------------------------------------------------*/ 00471 00472 int treesearch::TsSetPoint(double detectorwidth, double zdistance, Hit *Ha, Hit *Hb, 00473 char *patterna, char *patternb, int *hasha, int *hashb, 00474 unsigned binwidth) 00475 { 00476 cerr << "TsSetPoint called" << endl; 00477 int num = 0; 00478 00479 00480 double dw2 = (detectorwidth/2.0); /* half-width of the detector (in cm) */ 00481 Hit *Hanext, *Hbnext; 00482 00483 /* ---- Compute the maximum separation between hits on the two 00484 planes such that these hits are paired together when 00485 forming the hit pattern. If the hits are separated by 00486 more than this distance, the hits are treated as unmatched 00487 hits when the hit pattern is formed. ---- */ 00488 00489 double rcSETrMaxSlope = 0.60;//THIS VALUE COMES FROM hrcset.h 00490 cerr << "rcSET not defined in the code. Error. CODE NEEDS REPLACING" << endl; 00491 00492 00493 /*double maxdistance = rcSET.rMaxSlope * zdistance; */ /* maximum separation 00494 for hits (in cm) */ 00495 double maxdistance = rcSETrMaxSlope * zdistance;//I inserted this 00496 //line to fill in the gap due to the rcSET object not being 00497 //defined. 00498 00499 /* ---- Go through the hit lists for each plane and set the hit 00500 pattern. For each set of paired hits, the hit pattern 00501 is updated with the overlap between the hits. For the 00502 unpaired hits, the hit pattern is updated with just the 00503 possible tracks through the hit. ---- */ 00504 00505 while( Ha || Hb) { /* for each hit in Ha and Hb */ 00506 num++; /* Keep count of how many hits are processed */ 00507 00508 /* ---- Call wireselection() to pick out a hit pair from the 00509 linked lists of hits for the two planes (*Ha and *Hb) 00510 00511 This function earches for hits in plane A (stored in the 00512 linked list of hits) starting with the hit pointed to by 00513 *Ha and plane B starting with the hit pointed to by *Hb. 00514 Its returns: (1) pointers to the next hits to consider 00515 for pairing (*Hanext and *Hbnext) and (2) links to the 00516 hits (in *Ha and *Hb) to be inserted into the hit 00517 pattern (for unpair hits, either *Ha or *Hb will be zero ---- */ 00518 00519 wireselection( &Ha, &Hb, &Hanext, &Hbnext, maxdistance); 00520 00521 /* ---- Update the hit pattern to include these hits ---- */ 00522 00523 if( Ha && Hb ) { 00524 00525 /* ---- CASE 1: If the planes have hits within the HRCset Maxslope, 00526 then resolve the left-right ambiguity by checking 00527 the four combinations for matching. When a 00528 combination matches, the hit pattern is updated by 00529 treating the hits together. If no matching 00530 combinations are found, the hit pattern is updated 00531 by treated the hits separately. ---- */ 00532 00533 int found = 0; /* clear "I found a pairing" flag */ 00534 if( fabs(Ha->rPos1 - Hb->rPos1) < maxdistance ) { /* left A w/left B */ 00535 found = 1; /* set "I found a pairing" flag */ 00536 setpoint(dw2, 00537 Ha->rPos1 , Ha->rTrackResolution, 00538 Hb->rPos1 , Hb->rTrackResolution, 00539 detectorwidth, binwidth , patterna, patternb, 00540 hasha, hashb); 00541 } 00542 if( fabs(Ha->rPos2 - Hb->rPos1) < maxdistance ) { /* rght A w/left B */ 00543 found = 1; /* set "I found a pairing" flag */ 00544 setpoint(dw2, 00545 Ha->rPos2, Ha->rTrackResolution, 00546 Hb->rPos1, Hb->rTrackResolution, 00547 detectorwidth, binwidth , patterna, patternb, 00548 hasha, hashb); 00549 } 00550 if( fabs(Ha->rPos1 - Hb->rPos2) < maxdistance ) { /* left B w/rght B */ 00551 found = 1; /* set "I found a pairing flag */ 00552 setpoint(dw2, 00553 Ha->rPos1, Ha->rTrackResolution, 00554 Hb->rPos2, Hb->rTrackResolution, 00555 detectorwidth, binwidth , patterna, patternb, 00556 hasha, hashb); 00557 } 00558 /* 270996 - correct Hb to Ha below - finally ,-) */ 00559 if( fabs(Ha->rPos2 - Hb->rPos2) < maxdistance ) { /* rght A w/rght B */ 00560 found = 1; /* set "I found a pairing" flag */ 00561 setpoint(dw2, 00562 Ha->rPos2, Ha->rTrackResolution, 00563 Hb->rPos2, Hb->rTrackResolution, 00564 detectorwidth, binwidth , patterna, patternb, 00565 hasha, hashb); 00566 } 00567 00568 /* ---- CASE 2: The search by wireselection() paired these two 00569 hits, but the left-right ambigiuity check failed 00570 find a combination trough which an allowable track 00571 could traverse. So, these hits are added into 00572 bit pattern as unpaired hits. ---- */ 00573 00574 if( !found ) { 00575 if( patternb ) { /* treating the two planes as individual 00576 tree-planes, so put hit into both patterns */ 00577 setpoint(dw2, 00578 Ha->rPos1, Ha->rTrackResolution, 00579 Hb->rPos1, Hb->rTrackResolution, 00580 detectorwidth, binwidth, patterna, patternb, 00581 hasha, hashb); 00582 setpoint(dw2, 00583 Ha->rPos2, Ha->rTrackResolution, 00584 Hb->rPos2, Hb->rTrackResolution, 00585 detectorwidth, binwidth, patterna, patternb, 00586 hasha, hashb); 00587 } else { /* treating the two planes as one tree-plane, 00588 so just put it into the one pattern */ 00589 _setpoints(Ha->rPos1+dw2 - 0.5 * maxdistance - Ha->rTrackResolution, 00590 Ha->rPos1+dw2 + 0.5 * maxdistance + Ha->rTrackResolution, 00591 detectorwidth, binwidth, patterna, hasha); 00592 _setpoints(Ha->rPos2+dw2 - 0.5 * maxdistance - Ha->rTrackResolution, 00593 Ha->rPos2+dw2 + 0.5 * maxdistance + Ha->rTrackResolution, 00594 detectorwidth, binwidth, patterna, hasha); 00595 _setpoints(Hb->rPos1+dw2 - 0.5 * maxdistance - Hb->rTrackResolution, 00596 Hb->rPos1+dw2 + 0.5 * maxdistance + Hb->rTrackResolution, 00597 detectorwidth, binwidth, patterna, hasha); 00598 _setpoints(Hb->rPos2+dw2 - 0.5 * maxdistance - Hb->rTrackResolution, 00599 Hb->rPos2+dw2 + 0.5 * maxdistance + Hb->rTrackResolution, 00600 detectorwidth, binwidth, patterna, hasha); 00601 } 00602 } 00603 00604 /* ---- CASE 3: If, after the search, a hit on a plane was not 00605 paired with a hit on the other plane (because the 00606 other plane didn't have a hit within the HRCset 00607 MaxSlope), this hit is treated individually when 00608 being included in the bit pattern. ---- */ 00609 00610 00611 } else if ( Ha ) { /* this hit's on plane A */ 00612 if( patternb ) { /* treating the two planes as two tree-planes, 00613 so insert the hit into the bit pattern for 00614 this plane */ 00615 _setpoint( dw2 + Ha->rPos1, Ha->rTrackResolution, 00616 detectorwidth, binwidth, patterna, hasha); 00617 _setpoint( dw2 + Ha->rPos2, Ha->rTrackResolution, 00618 detectorwidth, binwidth, patterna, hasha); 00619 } else { /* treating the two planes as one tree-plane, 00620 so insert this hit into the bit pattern for 00621 this single tree-plane */ 00622 _setpoints(Ha->rPos1+dw2 - 0.5 * maxdistance - Ha->rTrackResolution, 00623 Ha->rPos1+dw2 + 0.5 * maxdistance + Ha->rTrackResolution, 00624 detectorwidth, binwidth, patterna, hasha); 00625 _setpoints(Ha->rPos2+dw2 - 0.5 * maxdistance - Ha->rTrackResolution, 00626 Ha->rPos2+dw2 + 0.5 * maxdistance + Ha->rTrackResolution, 00627 detectorwidth, binwidth, patterna, hasha); 00628 } 00629 } else { /* this hit's on plane B */ 00630 if( patternb ) { /* treating the two planes as two tree-planes, 00631 so insert the hit into the bit pattern for 00632 this plane */ 00633 _setpoint( dw2 + Hb->rPos1, Hb->rTrackResolution, 00634 detectorwidth, binwidth, patternb, hashb); 00635 _setpoint( dw2 + Hb->rPos2, Hb->rTrackResolution, 00636 detectorwidth, binwidth, patternb, hashb); 00637 } else { /* treating the two planes as one tree-planes, 00638 so insert the hit into the hit pattern for 00639 this single tree-plane */ 00640 _setpoints(Hb->rPos1+dw2 - 0.5 * maxdistance - Hb->rTrackResolution, 00641 Hb->rPos1+dw2 + 0.5 * maxdistance + Hb->rTrackResolution, 00642 detectorwidth, binwidth, patterna, hasha); 00643 _setpoints(Hb->rPos2+dw2 - 0.5 * maxdistance - Hb->rTrackResolution, 00644 Hb->rPos2+dw2 + 0.5 * maxdistance + Hb->rTrackResolution, 00645 detectorwidth, binwidth, patterna, hasha); 00646 } 00647 } 00648 Ha = Hanext; /* Continue with the next set of hits */ 00649 Hb = Hbnext; 00650 } 00651 00652 return num; 00653 } 00654 //________________________________________________________________________ 00655 00656 /*---------------------------------------------------------------------------*\ 00657 00658 exists() - this function searches through the link-list of valid treelines 00659 to see if the bit pattern for the specified treenode has already 00660 been accepted as a valid treeline. 00661 00662 inputs: (1) int *newa - 00663 (2) int front - 00664 (3) int back - 00665 (4) TreeLine *treeline - pointer to the linked-list of treelines 00666 00667 outputs: (1) int exists() - =0 if a treeline is not located 00668 =1 otherwise 00669 00670 \*---------------------------------------------------------------------------*/ 00671 00672 int treesearch::exists( int *newa, int front, int back, TreeLine *treeline) 00673 00674 { 00675 int *olda; 00676 int i, newmiss, oldmiss, diff; 00677 TreeLine *tl; 00678 int over; 00679 00680 newmiss = 0; 00681 for( i = 0; i < tlayers; i++ ) 00682 if( !newa[i] ) 00683 newmiss ++; 00684 00685 for( tl = treeline ; tl; tl = tl->next ) { /* loop over the treelines */ 00686 00687 if( tl->isvoid ) /* if the treeline has been voided, go onto next one */ 00688 continue; 00689 over = 0; 00690 if( tl->a_beg <= front && front <= tl->a_end ) 00691 over ++; 00692 if( tl->b_beg <= back && back <= tl->b_end ) 00693 over ++; 00694 00695 if( over == 2 ){ 00696 //cerr << "over = 2" << endl; 00697 return 1; 00698 } 00699 if( over == 0 ) 00700 continue; 00701 00702 olda = tl->hasharray; 00703 oldmiss = 0; 00704 diff = 0; 00705 00706 for( i = 0; i < tlayers; i++ ) { 00707 if( !olda[i] ) { 00708 oldmiss ++; 00709 } else { 00710 if( newa[i] && olda[i] != newa[i] ) { 00711 diff = 1; 00712 break; 00713 } 00714 } 00715 } 00716 if( !diff ) { 00717 if( (newmiss == 0 && oldmiss == 0) || 00718 (!newa[0] && !olda[0]) || 00719 (!newa[tlayers-1] && !olda[tlayers-1]) || 00720 (newmiss && !oldmiss && (!newa[0] || !newa[tlayers-1]))|| 00721 (oldmiss && !newmiss && (!olda[0] || !olda[tlayers-1])) 00722 ) { 00723 if( tl->a_beg > front ) 00724 tl->a_beg = front; 00725 if( tl->a_end < front ) 00726 tl->a_end = front; 00727 if( tl->b_beg > back ) 00728 tl->b_beg = back; 00729 if( tl->b_end < back ) 00730 tl->b_end = back; 00731 //cerr << "!diff" << endl; 00732 return 1; 00733 } 00734 } 00735 } 00736 return 0; 00737 } 00738 00739 //________________________________________________________________________ 00740 /*---------------------------------------------------------------------------*\ 00741 00742 _TsSearch() - this highly recursive function implements the treesearch 00743 algorithm. For a specified list of nodenodes, this function 00744 examines the attached treenode. If the bit pattern in the 00745 treenode does not match the bit pattern from the event, the 00746 function looks at the next nodenode. Otherwise, the function 00747 will call itself to see if any of the sons of this treenode 00748 at the next level of bin-division match the bit pattern from 00749 the event. This recursive calling will continue until a 00750 treenode at the deepest level of bin-division is located 00751 inside the bit pattern from the event. Since the pattern in 00752 this treenode is represents a valid treeline for the event, a 00753 treeline is constructed from the treenode and then appended 00754 to the linked list of treelines being accumulated by the 00755 treesearch. 00756 00757 inputs: (1) shortnode *node - pointer the first nodenode in the linked 00758 list of nodenodes to be checked. 00759 (2) int level - depth of the treesearch, i.e. the 00760 level of the bin-division 00761 (3) int offset - the offset of the treenode within the 00762 pattern 00763 (4) int reverse - the "pattern is flipped" flag 00764 00765 global 00766 inputs: (1) char **static_pattern - the bit pattern from the event 00767 (2) char **static_hash - 00768 (3) int static_maxlevel - 00769 (4) int static_front - 00770 00771 outputs: there are no explicit output from this function. 00772 00773 global 00774 outputs: Treelin *trelin - pointer to the link-list of valid treelines. 00775 00776 \*---------------------------------------------------------------------------*/ 00777 void treesearch::_TsSearch(shortnode *node,int level,int offset,int row_offset,int reverse,int numWires) 00778 { 00779 unsigned long pattern_start = (unsigned long)0xffffffffL; 00780 /* to add to the offset */ 00781 unsigned long pattern_offset; 00782 shorttree *tree; /* for searching in children of node*/ 00783 shortnode **cnode; 00784 shortnode **matchsons,**matchsonsnext; 00785 TreeLine *lineptr; /* evt. append to the treeline */ 00786 int *tree_pattern; 00787 int hashpat[TLAYERS]; 00788 int off, rev, off2, nlevel = level+1, i, bin,x; 00789 int miss, matched,nullhits; 00790 int frontbin, backbin; 00791 int firstwire=100,lastwire=-1; 00792 //int front_miss = 0; 00793 00794 if( trelinanz > TREELINABORT ) 00795 return; 00796 00797 00798 /* ---- Compute the offset in the bit pattern for the start 00799 position of this level of bin-division ---- */ 00800 pattern_start <<= level+1; 00801 pattern_start &= (unsigned long)0xffffffffL>>(32-static_maxlevel); 00802 //cout << "pattern start = " << pattern_start << endl; 00803 if(level == 0){ 00804 for(int u=nullhits = 0;u<tlayers;u++){ 00805 if(static_pattern[u+row_offset][pattern_start]) 00806 has_hits[u]=1; 00807 else{ 00808 has_hits[u]=0; 00809 nullhits++; 00810 } 00811 } 00812 } 00813 /* 00814 for(int s =0;s<tlayers;s++){ 00815 for(int t=0;t<15;t++){ 00816 if(static_pattern[s+row_offset][t]) 00817 cout << "1"; 00818 else 00819 cout << "0"; 00820 cout << ","; 00821 } 00822 cout << endl; 00823 } 00824 */ 00825 /* ---- Look at the trees attached to each nodenode on the 00826 specified nodenode linked list ---- */ 00827 if(numWires){ 00828 while( node) { /* search in all nodes */ 00829 00830 tree = node->tree; 00831 //tree->print(); 00832 00833 /* ---- Is the hit pattern in this treenode valid for this level 00834 of the treesearch? ---- */ 00835 //cout << "minlevel = " << tree->minlevel << endl; 00836 //cout << "row_offset = " << row_offset << endl; 00837 if( tree->minlevel > level+1) { /* check for level boundaries */ 00838 cerr << "hrm..." << endl; 00839 node = node->next; /* no, so look at the next nodenode */ 00840 continue; 00841 } 00842 00843 /* ---- Match the hit pattern to this treenode by checking the 00844 hit pattern for each tree-plane to see if the bits 00845 specified in the treenode are on ---- */ 00846 pattern_offset = pattern_start + offset; 00847 tree_pattern = tree->bit; 00848 if( reverse ) { 00849 //cout << "reversed..." << endl; 00850 for(i = matched = 0; i < tlayers; i++) { /* loop over tree-planes */ 00851 x = (*tree_pattern++); 00852 cerr << "ERROR : reversed patterns need checking/debugging" << endl; 00853 if(static_pattern[i+row_offset][pattern_offset - x]) { 00854 matched++; /* number of matched tree-planes */ 00855 if(i<firstwire)firstwire=i; 00856 if(i>lastwire)lastwire=i; 00857 } 00858 } 00859 } else { 00860 //cout << numWires << "," << tlayers << endl; 00861 for(i = matched = 0; i < tlayers; i++) { /* loop over tree-planes */ 00862 x = (*tree_pattern++); 00863 if(static_pattern[i+row_offset][pattern_offset + x]) { 00864 matched++; /* number of matched tree-planes */ 00865 if(i<firstwire)firstwire=i; 00866 if(i>lastwire && x!=0)lastwire=i; 00867 00868 } 00869 else if(x == 0 && has_hits[i]==0){ 00870 matched++; //matching null hits which are allowed in these patterns 00871 } 00872 } 00873 } 00874 //cout << "matched = " << matched << endl; 00875 /* ---- Check if there was the treenode is match now that the 00876 matching has been completely tested. ---- */ 00877 00878 if( matched == tlayers && nullhits <5) { 00879 //cout << "match found " << endl; 00880 //cout << "sons are " << endl; 00881 //cout << "-----------" << endl; 00882 matchsons = tree->son; 00883 /* 00884 for(int a=0;a<4;a++){ 00885 if(matchsons){ 00886 (*matchsons)->tree->print(); 00887 *matchsonsnext = (*matchsons)->next; 00888 matchsons = matchsonsnext; 00889 } 00890 else 00891 break; 00892 } 00893 cout << "-----------"<< endl; 00894 */ 00895 /* ---- Yes, there is a match, so now check if all the levels 00896 of the treesearch have been done. If so, then we have 00897 found a valid treeline. ---- */ 00898 00899 if( level == static_maxlevel-1 ) { 00900 //cerr << "inserting treeline..." << endl; 00901 /* all level done -> insert treeline */ 00902 /* ---- ---- */ 00903 backbin = reverse ? 00904 offset - tree->bit[tlayers-1] : offset + tree->bit[tlayers-1]; 00905 //cerr << "back =" << backbin << endl; 00906 if(reverse){ 00907 backbin = 0; 00908 for(i=0;i<tlayers;i++){ 00909 if(offset-tree->bit[i] > backbin) 00910 backbin = offset-tree->bit[i]; 00911 } 00912 00913 } 00914 else{ 00915 backbin = 0; 00916 for(i=0;i<tlayers;i++){ 00917 if(offset+tree->bit[i] > backbin) 00918 backbin = offset+tree->bit[i]; 00919 } 00920 } 00921 //cerr << "back =" << backbin << endl; 00922 frontbin = reverse ? 00923 offset - tree->bit[0] : offset + tree->bit[0]; 00924 //backbin = reverse ? 00925 // offset - tree->bit[tlayers-1] : offset + tree->bit[tlayers-1]; 00926 for( i = 0; i < tlayers; i++ ) { 00927 bin = reverse ? 00928 offset - tree->bit[i] : offset + tree->bit[i]; 00929 hashpat[i] = static_hash[i][bin]; 00930 } 00931 00932 miss = 0; 00933 if( static_pattern[0+row_offset][frontbin] == 0) 00934 miss = 1; 00935 else if( static_pattern[tlayers-1+row_offset][backbin] == 0) 00936 miss = 1; 00937 if( !exists( hashpat, frontbin, backbin,trelin) ) { 00938 if(opt.showMatchingPatterns)tree->print(); 00939 /* if track is unknown */ 00940 lineptr = (TreeLine*)malloc( sizeof(TreeLine)); /* create new one */ 00941 assert(lineptr); 00942 trelinanz ++; /* number of treelines found */ 00943 memcpy( lineptr->hasharray, hashpat, /* make space for this */ 00944 sizeof(int)*tlayers); /* treeline */ 00945 lineptr->Used = false; /* this treeline isn't used yet */ 00946 //lineptr->method = treelinemethod; /* remember the method */ 00947 lineptr->a_beg = frontbin; /* */ 00948 lineptr->a_end = frontbin; 00949 lineptr->b_beg = backbin; /* */ 00950 lineptr->b_end = backbin; 00951 lineptr->chi = 0.0; /* no chi2 value yet */ 00952 lineptr->isvoid = false; /* the treeline is voided yet */ 00953 lineptr->nummiss = miss; /* only used until TreeLineSort */ 00954 lineptr->next = trelin; /* put this treeline into the */ 00955 trelin = lineptr; /* linked-list of treelines */ 00956 lineptr->r3offset = row_offset; 00957 lineptr->firstwire = firstwire; 00958 lineptr->lastwire = lastwire; 00959 //Statist[treelinemethod]. /* update HRC statistics */ 00960 //TreeLinesGenerated[viswer][vispar][visdir-u_dir] ++; 00961 //cerr << lineptr->a_beg << "," << lineptr->a_end << "," <<endl; 00962 } 00963 } else { /* check son patterns */ 00964 for( rev = 0; rev < 4; rev += 2) { 00965 cnode = tree->son + rev; 00966 if( rev ^ reverse) { 00967 off2 = (offset << 1) + 1; 00968 for( off = 0; off < 2; off++) 00969 _TsSearch( *cnode++, nlevel, off2 - off, row_offset,2,numWires); 00970 } else { 00971 off2 = offset << 1; 00972 for( off = 0; off < 2; off++) 00973 _TsSearch( *cnode++, nlevel, off2 + off,row_offset, 0,numWires); 00974 } 00975 } /* highly optimized - time critical */ 00976 //cout << "sons done" << endl; 00977 } 00978 } 00979 node = node->next; /* ok, there wasn't a match, so go onto the 00980 next nodenode */ 00981 } /* end of loop over nodenodes */ 00982 } 00983 else{ 00984 while( node) { /* search in all nodes */ 00985 tree = node->tree; 00986 00987 /* ---- Is the hit pattern in this treenode valid for this level 00988 of the treesearch? ---- */ 00989 00990 if( tree->minlevel >= level) { /* check for level boundaries */ 00991 node = node->next; /* no, so look at the next nodenode */ 00992 continue; 00993 } 00994 00995 /* ---- Match the hit pattern to this treenode by checking the 00996 hit pattern for each tree-plane to see if the bits 00997 specified in the treenode are on ---- */ 00998 00999 pattern_offset = pattern_start + offset; 01000 tree_pattern = tree->bit; 01001 01002 if( reverse ) { 01003 for(i = matched = 0; i < tlayers; i++) { /* loop over tree-planes */ 01004 if(static_pattern[i][pattern_offset - *tree_pattern++]) { 01005 matched++; /* number of matched tree-planes */ 01006 } 01007 } 01008 } else { 01009 01010 for(i = matched = 0; i < tlayers; i++) { /* loop over tree-planes */ 01011 if(static_pattern[i][pattern_offset + *tree_pattern++]) { 01012 matched++; /* number of matched tree-planes */ 01013 } 01014 } 01015 } 01016 01017 /* ---- Check if there was the treenode is match now that the 01018 matching has been completely tested. ---- */ 01019 01020 if( matched == tlayers ) { 01021 01022 /* ---- Yes, there is a match, so now check if all the levels 01023 of the treesearch have been done. If so, then we have 01024 found a valid treeline. ---- */ 01025 01026 if( level == static_maxlevel-1 ) { 01027 /* all level done -> insert treeline */ 01028 01029 01030 /* ---- ---- */ 01031 01032 frontbin = reverse ? 01033 offset - tree->bit[0] : offset + tree->bit[0]; 01034 backbin = reverse ? 01035 offset - tree->bit[tlayers-1] : offset + tree->bit[tlayers-1]; 01036 for( i = 0; i < tlayers; i++ ) { 01037 bin = reverse ? 01038 offset - tree->bit[i] : offset + tree->bit[i]; 01039 hashpat[i] = static_hash[i][bin]; 01040 } 01041 01042 miss = 0; 01043 if( static_pattern[0][frontbin] == 0) 01044 miss = 1; 01045 else if( static_pattern[tlayers-1][backbin] == 0) 01046 miss = 1; 01047 01048 if( !exists( hashpat, frontbin, backbin,trelin) ) { 01049 /* if track is unknown */ 01050 //THE FOLLOWING LINE USED QMALLOC, BUT IT WAS 01051 //REPLACED DUE TO IT NOT EXISTING AT THE TIME 01052 //REVISED THIS PROGRAM. 01053 lineptr = (TreeLine*)malloc( sizeof(TreeLine)); /* create new one */ 01054 assert(lineptr); 01055 trelinanz ++; /* number of treelines found */ 01056 memcpy( lineptr->hasharray, hashpat, /* make space for this */ 01057 sizeof(int)*tlayers); /* treeline */ 01058 lineptr->Used = false; /* this treeline isn't used yet */ 01059 //lineptr->method = treelinemethod; /* remember the method */ 01060 lineptr->a_beg = frontbin; /* */ 01061 lineptr->a_end = frontbin; 01062 lineptr->b_beg = backbin; /* */ 01063 lineptr->b_end = backbin; 01064 lineptr->chi = 0.0; /* no chi2 value yet */ 01065 lineptr->isvoid = false; /* the treeline is voided yet */ 01066 lineptr->nummiss = miss; /* only used until TreeLineSort */ 01067 lineptr->next = trelin; /* put this treeline into the */ 01068 trelin = lineptr; /* linked-list of treelines */ 01069 //Statist[treelinemethod]. /* update HRC statistics */ 01070 //TreeLinesGenerated[viswer][vispar][visdir-u_dir] ++; 01071 } 01072 } else { /* check son patterns */ 01073 for( rev = 0; rev < 4; rev += 2) { 01074 cnode = tree->son + rev; 01075 if( rev ^ reverse) { 01076 off2 = (offset << 1) + 1; 01077 for( off = 0; off < 2; off++) 01078 _TsSearch( *cnode++, nlevel, off2 - off,0, 2,0); 01079 } else { 01080 off2 = offset << 1; 01081 for( off = 0; off < 2; off++) 01082 _TsSearch( *cnode++, nlevel, off2 + off,0, 0,0); 01083 } 01084 } /* highly optimized - time critical */ 01085 } 01086 } 01087 node = node->next; /* ok, there wasn't a match, so go onto the 01088 next nodenode */ 01089 } /* end of loop over nodenodes */ 01090 } 01091 return; 01092 } 01093 //________________________________________________________________________ 01094 /*---------------------------------------------------------------------------*\ 01095 01096 TsSearch() - this function initiates the treesearch for a set of tree- 01097 detectors by calling the _TsSearch() function described 01098 above. 01099 01100 inputs: (1) shortnode *node - 01101 (2) char *pattern[4] - 01102 (3) int *hashpat[4] - 01103 (4) int maxlevel - 01104 (5) int numWires - 01105 01106 outputs: there are no explicit output from this function. However, 01107 implicitly, it creates the linked list of valid treelines. 01108 01109 \*---------------------------------------------------------------------------*/ 01110 01111 void treesearch::TsSearch(shortnode *node, char *pattern[16], int *hashpat[16], 01112 int maxlevel, int numWires){ 01113 static_maxlevel = maxlevel; 01114 static_pattern = pattern; 01115 static_hash = hashpat; 01116 /*static_front = front;*/ 01117 if(numWires){//The region 3 version of tssearch 01118 for(int i=0;i<=numWires-tlayers;i++){ 01119 _TsSearch( node, 0, 0,i, 0,numWires); 01120 } 01121 } 01122 else _TsSearch(node,1,0,0,0,0); 01123 } 01124 //________________________________________________________________________ 01125 01126 01127 01128 01129 01130
1.4.6