001 #include <stdio.h>
002 #include <stdlib.h> // for exit()
003
004 #include <math.h>
005 #include <iomanip>
006
007 #include "RPCgeometry/RPCGeometry.h"
008 #include "MuonCablingTools/RPCdecoder.h"
009
010
011 static MultiChamber *T1muonGeo[MaxJtype][MaxJff*(2*MaxJzz+1)];
012
013
014 RPCGeometry* RPCGeometry::s_instance = 0;
015 bool RPCGeometry::DumpMap = false;
016 bool RPCGeometry::DumpLog = false;
017 bool RPCGeometry::DumpCon = false;
018 bool RPCGeometry::DebuGeo = false;
019
020 double RPCGeometry::s_ScaleFactor = 1.;
021 float RPCGeometry::s_Zsign = +1.;
022 double RPCGeometry::s_pi = 3.14159265358979323846;
023
024 bool RPCGeometry::s_testbeam = false;
025 int RPCGeometry::s_JobMinForEncodeProj = 9999;
026
027 std::string RPCGeometry::FileName = "LVL1conf.dump";
028
029 #ifndef LVL1_STANDALONE
030 Amdcsimrec* RPCGeometry::s_amdc = 0;
031 #endif
032
033 bool RPCGeometry::s_status = true;
034
035
036 RPCGeometry::RPCGeometry(void)
037 {
038 for (int sect=0;sect<64;++sect)
039 {
040 for(int stat=First;stat<=Extended;++stat)
041 {
042 for(int i=0;i<2;++i)
043 {
044 LogicSector[sect][stat][i] = 0;
045 }
046 }
047 }
048
049 if(s_testbeam) s_Zsign=-1.;
050
051
052
053 projMap = 0;
054
055 buildgeo();
056
057
058
059
060
061
062
063
064 setupmap();
065 setuplogic();
066
067 if(DumpCon)
068 {
069 std::ofstream dump;
070 dump.open(RPCGeometry::FileName.c_str());
071 write_geodata(dump);
072 }
073 }
074
075 RPCGeometry::~RPCGeometry(void)
076 {
077
078 for (int Jtype=1;Jtype<=MaxJtype;++Jtype)
079 {
080 for (int Jff=1;Jff<=Mff;++Jff)
081 {
082 for (int Jzz=-Mzz;Jzz<=Mzz;++Jzz)
083 {
084 MultiChamber *cham =
085 T1muonGeo[Jtype-1][(Jff-1)*(2*Mzz+1)+Jzz+Mzz];
086
087 if (cham)
088 {
089 delete cham;
090 }
091 }
092 }
093 }
094
095 }
096
097
098 const RPCGeometry*
099 RPCGeometry::instance(void)
100 {
101 if (! s_instance) {
102 s_instance = new RPCGeometry();
103 }
104
105 return s_instance;
106 }
107
108
109
110
111
112
113 void
114 RPCGeometry::myprin(char c,int i) const
115 {
116 int j;
117
118 for(j=0;j<i;j++) printf("%c",c);
119 }
120
121 void
122 RPCGeometry::cham_info(MultiChamber *cham,int stat)
123 {
124 printf(" Chamber %s, Phyisic sector = %d, Z index = %d",
125 cham->name.c_str(),cham->physics_sector+1,
126 cham->amdb_z_index);
127 if (stat > 3) printf("\n\n");
128 else printf(", Station = %d\n\n",stat + 1);
129 }
130
131 void
132 RPCGeometry::pack_info(RPC_pack *pack)
133 {
134 printf("\n Jobj = %d, Station = %d, Z index = %d\n",
135 pack->Jobj, pack->station, pack->lvl1_z_index);
136 printf(" Pack center: X=%8.2f, Y=%8.2f, Z=%8.2f\n",
137 pack->pack_coo[0],pack->pack_coo[1],pack->pack_coo[2]);
138 printf(" #z strips = %d (pitch=%5.3f cm), #phi strips = %d",
139 pack->nstrip_z, pack->stripsize_z, pack->nstrip_f);
140 printf(" (pitch=%5.2f cm)\n",pack->stripsize_f);
141 printf(" Z coordinate of first strip: Z = %8.2f\n",pack->z_proj[2]);
142 }
143
144 void
145 RPCGeometry::list_of_packs(MultiChamber *cham)
146 {
147 cham_info(cham,3);
148 printf(" has %d RPC packs:\n\n",cham->nRPC_packs);
149 RPC_pack *pack = cham->RPC_packs;
150 while(pack)
151 {
152 pack_info(pack);
153 pack = pack->next;
154 }
155 }
156
157 void
158 RPCGeometry::first_packs_info(RPC_pack *pack1,RPC_pack *pack2,int stat)
159 {
160 if (pack1)
161 {
162 printf("\n The chamber\n");
163 pack_info(pack1);
164 printf("\n has been chosen as first RPC(proj 0) of station %d\n",stat);
165 }
166 if (pack2)
167 {
168 printf("\n The chamber\n");
169 pack_info(pack2);
170 printf("\n has been chosen as first RPC(proj 1) of station %d\n",stat);
171 }
172 }
173
174 void
175 RPCGeometry::check_error(MultiChamber *cham,std::string message)
176 {
177 printf("RPCGeometry, check_links: ");
178 printf("%s\n",message.c_str());
179 cham_info(cham,10);
180 exit(1);
181 }
182
183 void
184 RPCGeometry::axis(__osstream **disp,std::string X, std::string Y)
185 {
186 *disp[6] << Y;
187 *disp[7] << "^";
188 for(int i=8;i<22;++i) *disp[i] << "|";
189 *disp[22] << "+-------------> " << X;
190 }
191
192 int
193 RPCGeometry::find_front(__osstream **disp)
194 {
195 unsigned int head = 0;
196 for (int i=0;i<22;++i)
197 {
198 if((disp[i])->COUNT > head) head = (disp[i])->COUNT;
199 }
200 return head;
201 }
202
203 void
204 RPCGeometry::align_disp(__osstream **disp,int start,int stop,int front)
205 {
206 for(int i=start;i<=stop;++i)
207 {
208 if(i<23)
209 {
210 int head_line = (disp[i])->COUNT;
211 for(int j=0;j<front - head_line;++j)
212 {
213 if(i==10) *disp[i] << "-";
214 else *disp[i] << " ";
215 }
216 }
217 }
218 }
219
220
221
222
223
224
225 void
226 RPCGeometry::dump_chamber(MultiChamber *cham)
227 {
228 __osstream *disp[23];
229
230 for(int stat=First;stat<=Extended;++stat)
231 {
232 RPC_pack* pack0 = cham->firstRPC_z[0][stat];
233 RPC_pack* pack1 = cham->firstRPC_z[1][stat];
234
235 if (pack0 || pack1)
236 {
237 #if (__GNUC__) && (__GNUC__ > 2)
238 for(int i=0;i<23;++i) disp[i] = new __osstream();
239 #else
240 for(int i=0;i<23;++i) disp[i] = new __osstream(display[i],80);
241 #endif
242
243
244 axis(disp,"Z","Phi");
245
246 float Phi = (pack0)? pack0->rotation_angle : pack1->rotation_angle;
247 *disp[10] << " Phi=";
248
249 align_disp(disp,0,22,(disp[10])->COUNT);
250
251 *disp[10] << std::setw(5) << std::setprecision(3)
252 << std::setiosflags(std::ios::fixed)
253 << Phi << "-";
254
255
256 while (pack0||pack1)
257 {
258 dump_pack_check(disp,pack1,1);
259 dump_pack_check(disp,pack0,0);
260
261 if(pack0) pack0 = pack0->next_in_z[0];
262 if(pack1) pack1 = pack1->next_in_z[1];
263
264 if(pack0 || pack1)
265 {
266 align_disp(disp,0,21,find_front(disp));
267 for(int i=0;i<10;++i) *disp[i] << " ";
268 *disp[10] << "--------";
269 for(int i=11;i<22;++i) *disp[i] << " ";
270 }
271 }
272
273 align_disp(disp,0,22,find_front(disp));
274
275
276 for(int i=0;i<23;++i) *disp[i] << " ";
277 axis(disp,"Z","R");
278 align_disp(disp,0,21,(disp[10])->COUNT);
279 for(int i=0;i<21;++i) *disp[i] << " ";
280
281 pack0 = cham->firstRPC_z[0][stat];
282 pack1 = cham->firstRPC_z[1][stat];
283
284 int head = (disp[10])->COUNT;
285 int offs = 0;
286
287
288 while (pack0||pack1)
289 {
290 int proj = (pack0)? 0 : 1;
291 RPC_pack* pack = (pack0)? pack0 : pack1;
292
293 RPC_pack * prev = pack->prev_in_z[proj];
294
295 float radius =sqrt(powf(pack->pack_coo[0],2) +
296 powf(pack->pack_coo[1],2));
297
298 float prev_r =(!prev)? 0 : sqrt(powf(prev->pack_coo[0],2) +
299 powf(prev->pack_coo[1],2));
300
301 if (prev)
302 {
303 if (pack->pack_coo[2] - pack->half_dim[1] <
304 prev->pack_coo[2] + prev->half_dim[1])
305 {
306 align_disp(disp,offs-2,offs+2,head-1);
307 } else
308 {
309 align_disp(disp,offs-2,offs+2,head+1);
310 }
311 }
312
313 if(!offs) offs = (radius > cham->Rmid)? 4 : 16;
314
315 if (radius > prev_r) --offs;
316 else if (radius < prev_r) ++offs;
317
318 if (radius == cham->Rmid) offs = 10;
319
320 *disp[offs] << "=========";
321 align_disp(disp,offs-1,offs+1,head);
322 *disp[offs-1] << " " << std::setw(6) << std::setprecision(1)
323 << std::setiosflags(std::ios::fixed)
324 << pack->gas_layer_radius[1];
325 *disp[offs+1] << " " << std::setw(6) << std::setprecision(1)
326 << std::setiosflags(std::ios::fixed)
327 << pack->gas_layer_radius[0];
328 head = (disp[offs])->COUNT;
329
330 if(pack0) pack0 = pack0->next_in_z[0];
331 if(pack1) pack1 = pack1->next_in_z[1];
332 }
333
334 *disp[10] << " ---- Rmid=" << std::setw(6) << std::setprecision(1)
335 << std::setiosflags(std::ios::fixed) << cham->Rmid << " ";
336 align_disp(disp,10,10,head);
337
338
339
340 myprin('=',79);
341 printf("\n");
342 cham_info(cham,stat);
343
344 for(int i=0;i<23;++i)
345 {
346 #if (__GNUC__) && (__GNUC__ > 2)
347 std::cout << (disp[i])->str() << std::endl;
348 #else
349 *disp[i] << std::ends;
350 printf("%s\n",display[i]);
351 #endif
352 delete disp[i];
353 }
354
355 }
356 }
357 }
358
359 void
360 RPCGeometry::dump_pack_check(__osstream **disp, RPC_pack *pack,int proj)
361 {
362 if (!pack) return;
363
364 int width = (pack->half_dim[1] > 50.)? 10 : 4;
365 float y_center = pack->pack_coo[0]*sin(-pack->rotation_angle) +
366 pack->pack_coo[1]*cos(-pack->rotation_angle);
367 if (proj == 1)
368 {
369 if (y_center - pack->half_dim[0] + 5. < 0. && !(fabsf(y_center) < 24.)
370 && !s_testbeam)
371 {
372 dump_pack_roof (disp[11],6,width);
373 dump_pack_bound (disp[10],0,width,pack->nstrip_z,0);
374 dump_pack_bound (disp[9],6,width,0,0);
375 dump_pack_bound (disp[8],6,width,0,0);
376 dump_pack_bound (disp[7],6,width,0,0);
377 }
378
379 if (fabsf(y_center) < 24. || s_testbeam)
380 {
381 dump_pack_roof (disp[10],0,width);
382 dump_pack_bound (disp[9],6,width,pack->nstrip_z,0);
383 dump_pack_bound (disp[8],6,width,0,0);
384 dump_pack_bound (disp[7],6,width,0,0);
385 }
386
387 if(y_center - pack->half_dim[0] - 5. > 0.)
388 {
389 dump_pack_roof (disp[8],6,width);
390 dump_pack_bound (disp[7],6,width,pack->nstrip_z,0);
391 dump_pack_bound (disp[6],6,width,0,0);
392
393 *disp[5] << std::setw(5) << std::setprecision(0);
394 *disp[5] << std::setiosflags(std::ios::fixed) << pack->half_dim[0];
395 dump_pack_bound (disp[5],1,width,pack->nstrip_f,1);
396 } else
397 {
398 *disp[6] << std::setw(5) << std::setprecision(0);
399 *disp[6] << std::setiosflags(std::ios::fixed) << pack->half_dim[0];
400 dump_pack_bound (disp[6],1,width,pack->nstrip_f,1);
401 dump_pack_bound (disp[5],6,width,0,0);
402 }
403
404 dump_pack_bound (disp[4],6,width,0,0);
405 dump_pack_roof (disp[3],6,width);
406
407 *disp[2] << " |";
408 for(int i=0;i<width/4;++i) *disp[2]<<" ";
409 *disp[2] << std::setw(4) << std::setprecision(0);
410 *disp[2] << std::setiosflags(std::ios::fixed) << pack->half_dim[1]*2;
411
412 *disp[1] << " |";
413
414 *disp[0] << " " << std::setw(7) << std::setprecision(0)
415 << std::setiosflags(std::ios::fixed)
416 << pack->pack_coo[2]-pack->half_dim[1];
417 }
418
419
420 if (proj == 0)
421 {
422 if (y_center + pack->half_dim[0] - 5. > 0. && !(fabsf(y_center) < 24.)
423 && !s_testbeam )
424 {
425 dump_pack_roof (disp[9],6,width);
426 dump_pack_bound (disp[10],0,width,pack->nstrip_z,0);
427 dump_pack_bound (disp[11],6,width,0,0);
428 dump_pack_bound (disp[12],6,width,0,0);
429 dump_pack_bound (disp[13],6,width,0,0);
430 }
431
432 if (fabsf(y_center) < 24. || s_testbeam)
433 {
434 dump_pack_bound (disp[11],6,width,pack->nstrip_z,0);
435 dump_pack_bound (disp[12],6,width,0,0);
436 dump_pack_bound (disp[13],6,width,0,0);
437 }
438
439 if(y_center + pack->half_dim[0] + 5. < 0.)
440 {
441 dump_pack_roof (disp[12],6,width);
442 dump_pack_bound (disp[13],6,width,pack->nstrip_z,0);
443 dump_pack_bound (disp[14],6,width,0,0);
444
445 *disp[15] << std::setw(5) << std::setprecision(0);
446 *disp[15] << std::setiosflags(std::ios::fixed) << pack->half_dim[0];
447 dump_pack_bound (disp[15],1,width,pack->nstrip_f,1);
448 } else
449 {
450 *disp[14] << std::setw(5) << std::setprecision(0);
451 *disp[14] << std::setiosflags(std::ios::fixed) << pack->half_dim[0];
452 dump_pack_bound (disp[14],1,width,pack->nstrip_f,1);
453 dump_pack_bound (disp[15],6,width,0,0);
454 }
455
456 dump_pack_bound (disp[16],6,width,0,0);
457 dump_pack_roof (disp[17],6,width);
458
459 *disp[18] << " |";
460 for(int i=0;i<width/4;++i) *disp[18]<<" ";
461 *disp[18] << std::setw(4) << std::setprecision(0);
462 *disp[18] << std::setiosflags(std::ios::fixed) << pack->half_dim[1]*2;
463
464 *disp[19] << " |";
465
466 *disp[20] << " " << std::setw(7) << std::setprecision(0)
467 << std::setiosflags(std::ios::fixed)
468 << pack->pack_coo[2]-pack->half_dim[1];
469 }
470 }
471
472 void
473 RPCGeometry::dump_pack_roof (__osstream *disp,int offset,int width)
474 {
475 for(int i=0;i<offset;++i) *disp << " ";
476 *disp << "+";
477 for(int i=0;i<width;++i) *disp << "-";
478 *disp << "+";
479 }
480
481 void
482 RPCGeometry::dump_pack_bound (__osstream *disp,
483 int offset,
484 int width,
485 int strips,
486 int xy)
487 {
488 for(int i=0;i<offset;++i) *disp << " ";
489 *disp << "|";
490 if (strips)
491 {
492 if(xy)
493 {
494 *disp << std::setw(2) << strips;
495 for(int i=0;i<width-2;++i) *disp << " ";
496 *disp << "|";
497 } else
498 {
499 for(int i=0;i<width/2-1;++i) *disp << " ";
500 *disp << std::setw(2) << strips;
501 for(int i=0;i<width/2-1;++i) *disp << " ";
502 *disp << "|";
503 }
504 } else
505 {
506 for(int i=0;i<width;++i) *disp << " ";
507 *disp << "|";
508 }
509 }
510
511
512 void
513 RPCGeometry::dump_logic_sectors(void)
514 {
515 __osstream *disp[23];
516 for(int stat=First;stat<=Extended;++stat)
517 {
518 for(int logic_sector = 62;logic_sector >= 0; logic_sector -= 2)
519 {
520 int upper_part = logic_sector;
521 int lower_part = (logic_sector%32)? logic_sector-1:logic_sector+31;
522
523 RPC_pack *next1 = LogicSector[upper_part][stat][0];
524 RPC_pack *next0 = LogicSector[lower_part][stat][0];
525 RPC_pack *pack1 = 0;
526 RPC_pack *pack0 = 0;
527 #if (__GNUC__) && (__GNUC__ > 2)
528 for(int i=0;i<23;++i) disp[i] = new __osstream();
529 #else
530 for(int i=0;i<23;++i) disp[i] = new __osstream(display[i],80);
531 #endif
532
533 int offs = (lower_part < 32)? 21 : 20;
534 int head = -1;
535
536 do
537 {
538
539 pack1 = next1;
540 pack0 = next0;
541
542 int width = 4;
543
544
545
546 if (pack1)
547 {
548 float y_center = pack1->pack_coo[0]*
549 sin(-pack1->rotation_angle) +
550 pack1->pack_coo[1]*cos(-pack1->rotation_angle);
551
552 if (y_center - pack1->half_dim[0] + 5. < 0. &&
553 !(fabsf(y_center) < 24.) && !s_testbeam)
554 {
555 dump_pack_roof (disp[11],0,width);
556 dump_pack_bound (disp[10],0,width,pack1->nstrip_z,0);
557 dump_pack_bound (disp[9],0,width,0,0);
558 dump_pack_bound (disp[8],0,width,0,0);
559 dump_pack_bound (disp[7],0,width,0,0);
560 }
561
562 if (fabsf(y_center) < 24. || s_testbeam)
563 {
564 dump_pack_roof (disp[10],0,width);
565 dump_pack_bound (disp[9],0,width,pack1->nstrip_z,0);
566 dump_pack_bound (disp[8],0,width,0,0);
567 dump_pack_bound (disp[7],0,width,0,0);
568 }
569
570 if(y_center - pack1->half_dim[0] - 5. > 0.)
571 {
572 dump_pack_roof (disp[8],0,width);
573 dump_pack_bound (disp[7],0,width,pack1->nstrip_z,0);
574 dump_pack_bound (disp[6],0,width,0,0);
575 dump_pack_bound (disp[5],0,width,pack1->nstrip_f,1);
576 } else
577 {
578 dump_pack_bound (disp[6],0,width,pack1->nstrip_f,1);
579 dump_pack_bound (disp[5],0,width,0,0);
580 }
581
582
583 dump_pack_id(disp[4],0,width,pack1->station,
584 pack1->lvl1_z_index);
585 dump_pack_roof (disp[3],0,width);
586
587 *disp[2] << " " << pack1->cham->name << pack1->cham->type;
588 }
589
590 if (pack0)
591 {
592 float y_center = pack0->pack_coo[0]*
593 sin(-pack0->rotation_angle) +
594 pack0->pack_coo[1]*cos(-pack0->rotation_angle);
595
596 if (y_center + pack0->half_dim[0] - 5. > 0. &&
597 !(fabsf(y_center) < 24.) && !s_testbeam)
598 {
599 dump_pack_roof (disp[9],0,width);
600 dump_pack_bound (disp[10],0,width,pack0->nstrip_z,0);
601 dump_pack_bound (disp[11],0,width,0,0);
602 dump_pack_bound (disp[12],0,width,0,0);
603 dump_pack_bound (disp[13],0,width,0,0);
604 }
605
606 if (fabsf(y_center) < 24. || s_testbeam)
607 {
608 dump_pack_bound (disp[11],0,width,pack0->nstrip_z,0);
609 dump_pack_bound (disp[12],0,width,0,0);
610 dump_pack_bound (disp[13],0,width,0,0);
611 }
612
613 if(y_center + pack0->half_dim[0] + 5. < 0.)
614 {
615 dump_pack_roof (disp[12],0,width);
616 dump_pack_bound (disp[13],0,width,pack0->nstrip_z,0);
617 dump_pack_bound (disp[14],0,width,0,0);
618 dump_pack_bound (disp[15],0,width,pack0->nstrip_f,1);
619 } else
620 {
621 dump_pack_bound (disp[14],0,width,pack0->nstrip_f,1);
622 dump_pack_bound (disp[15],0,width,0,0);
623 }
624
625
626 dump_pack_id(disp[16],0,width,pack0->station,
627 pack0->lvl1_z_index);
628 dump_pack_roof (disp[17],0,width);
629
630 *disp[18] << " " << pack0->cham->name << pack0->cham->type;
631
632 align_disp(disp,0,18,find_front(disp));
633
634
635
636 if(!s_testbeam)
637 {
638 RPC_pack * prev = pack0->prev_in_z[0];
639
640 float radius = give_radius(pack0);
641 float prev_r = (!prev)? 0 : give_radius(prev);
642
643 if (prev && LogicSector[lower_part][stat][0] != pack0)
644 {
645 int next_offs = give_radial_offset(radius,prev_r,offs);
646
647 if(pack0->cham != prev->cham)
648 {
649 *disp[offs] << "]";
650 align_disp(disp,19,22,(disp[17])->COUNT - 6);
651 if(same_cham(pack0->next_in_z[0],pack0))
652 *disp[next_offs] << " [====";
653 else
654 *disp[next_offs] << "[====";
655
656 } else
657 {
658 align_disp(disp,19,22,head);
659 if (pack0->pack_coo[2] - pack0->half_dim[1] <
660 prev->pack_coo[2] + prev->half_dim[1])
661 {
662 *disp[next_offs] << "====";
663
664 } else
665 {
666 *disp[next_offs] << " ====";
667 }
668 }
669
670 offs = next_offs;
671
672 } else
673 {
674 if(same_cham(pack0->next_in_z[0],pack0))
675 *disp[offs] << " [====";
676 else
677 *disp[offs] << "[====";
678 }
679
680 head = (disp[offs])->COUNT;
681
682 }
683
684 }
685
686 if (LogicSector[lower_part][stat][1] == pack0 && pack0)
687 {
688 *disp[offs] << "]";
689 }
690
691
692 if (pack1 && pack1->next_in_z[1]) next1 = pack1->next_in_z[1];
693 if (pack0 && pack0->next_in_z[0]) next0 = pack0->next_in_z[0];
694
695 } while (pack1 && pack1 != LogicSector[upper_part][stat][1] &&
696 pack0 && pack0 != LogicSector[lower_part][stat][1]);
697
698
699
700 myprin('=',21);
701 std::string sign = (lower_part < 32)? "negative" : "positive";
702 printf(" %s half barrel, Station = %1d ",sign.c_str(),stat+1);
703 myprin('=',22);
704 printf("\n");
705
706
707 myprin(' ',30);
708 printf("-------------------> Z Logic Sector = %2d\n",
709 upper_part%32);
710
711
712 for(int i=0;i<23;++i)
713 {
714 #if (__GNUC__) && (__GNUC__ > 2)
715 std::cout << (disp[i])->str() << std::endl;
716 #else
717 *disp[i] << std::ends;
718 printf("%s\n",display[i]);
719 #endif
720 delete disp[i];
721 }
722
723
724 myprin(' ',61);
725 printf("Logic Sector = %2d\n",lower_part%32);
726
727 }
728 }
729 }
730
731
732 void
733 RPCGeometry::dump_pack_id (__osstream *disp,
734 int offset,
735 int width,
736 int station,
737 int z_index)
738 {
739 for(int i=0;i<offset;++i) *disp << " ";
740 *disp << "|";
741 for(int i=0;i<(width-4)/2;++i) *disp << " ";
742 *disp << std::setw(1) << station + 1 << "*" << std::setw(2)
743 << std::setfill('0')<<z_index;
744 for(int i=0;i<(width-4)/2;++i) *disp << " ";
745 *disp << "|";
746 }
747
748
749 int
750 RPCGeometry::give_radial_offset(float current_rad, float previous_rad, int offs)
751 {
752 float delta = fabs(current_rad - previous_rad);
753 if (delta <= 1.) return offs;
754 if (current_rad > previous_rad) return(--offs);
755 else return(++offs);
756 }
757
758
759
760
761 void
762 RPCGeometry::rreset()
763 {
764 Removed = 0;
765 }
766
767
768 template<class entry> entry*
769 RPCGeometry::remove(entry* first, entry* item)
770 {
771 if(!item) return first;
772
773 entry* tmp = first;
774
775 if(tmp == item)
776 {
777 entry* newfirst = first->next;
778 item->next = 0;
779 ++Removed;
780 return newfirst;
781 }
782
783 entry* prv = tmp;
784 tmp = prv->next;
785
786 while(tmp)
787 {
788 if(tmp == item)
789 {
790 prv->next = tmp->next;
791 item->next = 0;
792 ++Removed;
793 return first;
794 }
795 prv = tmp;
796 tmp = prv->next;
797 }
798
799 printf("RPCGeometry, remove: error in removing pack!\n");
800 exit(1);
801 }
802
803
804 template<class entry> entry*
805 RPCGeometry::give_last(entry* first)
806 {
807 if(!first) return 0;
808
809 entry* last = first;
810
811 while(last->next)
812 {
813 last = last->next;
814 }
815
816 return last;
817 }
818
819
820 void
821 RPCGeometry::set_first_packs(MultiChamber* cham)
822 {
823 for (int stat=First;stat<=Extended;++stat)
824 {
825 float Zmin = 99999.;
826 RPC_pack* pack = cham->RPC_packs;
827 while(pack)
828 {
829 if(s_Zsign*pack->pack_coo[2]<=Zmin && pack->station==stat)
830 {
831 Zmin = s_Zsign*pack->pack_coo[2];
832 if(pack->f_proj[0]) cham->firstRPC_z[0][stat] = pack;
833 if(pack->f_proj[1]) cham->firstRPC_z[1][stat] = pack;
834 }
835 pack = pack->next;
836 }
837 if(cham->firstRPC_z[0][stat] == cham->firstRPC_z[1][stat])
838 {
839 cham->RPC_packs =remove(cham->RPC_packs,cham->firstRPC_z[0][stat]);
840 } else
841 {
842 cham->RPC_packs =remove(cham->RPC_packs,cham->firstRPC_z[0][stat]);
843 cham->RPC_packs =remove(cham->RPC_packs,cham->firstRPC_z[1][stat]);
844 }
845
846 if (DebuGeo) first_packs_info(cham->firstRPC_z[0][stat],
847 cham->firstRPC_z[1][stat],stat);
848 }
849 }
850
851
852 RPC_pack*
853 RPCGeometry::find_next_pack(MultiChamber* cham, RPC_pack* pack,int sect)
854 {
855 RPC_pack* it = cham->RPC_packs;
856 RPC_pack* found = 0;
857
858
859
860 float dist = 99999;
861 while (it)
862 {
863 float tmp_dist=(s_Zsign*(it->pack_coo[2]) - s_Zsign*(pack->pack_coo[2]));
864 if(tmp_dist < dist &&
865 s_Zsign*(it->pack_coo[2]) > s_Zsign*(pack->pack_coo[2]) &&
866 it->station == pack->station &&
867 it->f_proj[sect])
868 {
869 dist = tmp_dist;
870 found = it;
871 }
872 it = it->next;
873 }
874 return found;
875 }
876
877
878 void
879 RPCGeometry::connect(RPC_pack* prev, RPC_pack* next)
880 {
881 if(!prev||!next) return;
882
883 if(prev->f_proj[0] && next->f_proj[0])
884 {
885 prev->next_in_z[0] = next;
886 next->prev_in_z[0] = prev;
887 }
888
889 if(prev->f_proj[1] && next->f_proj[1])
890 {
891 prev->next_in_z[1] = next;
892 next->prev_in_z[1] = prev;
893 }
894
895
896
897
898
899
900 }
901
902
903 MultiChamber*
904 RPCGeometry::find_first_chamber(MultiChamber *cham, int stat, int proj)
905 {
906 float Zmin = 99999.;
907 MultiChamber *first = 0;
908
909 while(cham)
910 {
911
912
913
914 float z = s_Zsign*((cham->firstRPC_z[proj][stat])->pack_coo[2]);
915 if(z <= Zmin)
916 {
917 Zmin = z;
918 first = cham;
919 }
920 cham = cham->next;
921 }
922 return first;
923
924 }
925
926
927 RPC_pack*
928 RPCGeometry::give_last_pack(MultiChamber *cham, int stat, int proj)
929 {
930 RPC_pack *last = cham->firstRPC_z[proj][stat];
931 while(last->next_in_z[proj]) last = last->next_in_z[proj];
932 return last;
933 }
934
935
936 float
937 RPCGeometry::give_radius(RPC_pack *pack)
938 {
939 float radius = pack->pack_coo[0]*cos(-pack->rotation_angle) -
940 pack->pack_coo[1]*sin(-pack->rotation_angle);
941 return radius;
942 }
943
944
945 int
946 RPCGeometry::same_cham(RPC_pack *one, RPC_pack *two)
947 {
948 if(!one || !two) return 0;
949 if(one->cham == two->cham) return 1;
950 else return 0;
951 }
952
953
954
955
956 void
957 RPCGeometry::setupmap()
958 {
959
960 for (int Jtype=1;Jtype<=MaxJtype;++Jtype)
961 {
962 for (int Jff=1;Jff<=Mff;++Jff)
963 {
964 for (int Jzz=-Mzz;Jzz<=Mzz;++Jzz)
965 {
966 MultiChamber *cham =
967 T1muonGeo[Jtype-1][(Jff-1)*(2*Mzz+1)+Jzz+Mzz];
968
969 if (cham)
970 {
971 if (DebuGeo) list_of_packs(cham);
972
973 rreset();
974 set_first_packs(cham);
975 for (int stat=First;stat<=Extended;++stat)
976 {
977 if(cham->firstRPC_z[0][stat]==cham->firstRPC_z[1][stat])
978 {
979 build_InterChamS_links(cham,cham->firstRPC_z[0][stat]);
980 }else
981 {
982 build_InterChamS_links(cham,cham->firstRPC_z[0][stat]);
983 build_InterChamS_links(cham,cham->firstRPC_z[1][stat]);
984 }
985 }
986
987 if (Removed != cham->nRPC_packs)
988 {
989 std::cout << "Jtype=" << Jtype << ", Jff=" << Jff
990 << ", Jzz=" << Jzz << std::endl;
991 printf("RPCGeometry error: number of RPC packs = %d, ",
992 cham->nRPC_packs);
993 printf("removed packs = %d\n",Removed);
994 exit(1);
995 }
996
997 check_map(cham);
998
999 if(DumpMap) dump_chamber(cham);
1000 }
1001 }
1002 }
1003 }
1004 }
1005
1006
1007 void
1008 RPCGeometry::build_InterChamS_links(MultiChamber* cham, RPC_pack* pack)
1009 {
1010 if(!pack) return;
1011
1012 RPC_pack* pack1 = find_next_pack(cham,pack,0);
1013 RPC_pack* pack2 = find_next_pack(cham,pack,1);
1014
1015 if(!pack1&&!pack2) return;
1016
1017 if(pack1 == pack2)
1018 {
1019 connect(pack,pack1);
1020 if(pack->f_proj[1]) cham->RPC_packs = remove(cham->RPC_packs,pack1);
1021
1022 build_InterChamS_links(cham,pack1);
1023 } else
1024 {
1025 connect(pack,pack1);
1026 connect(pack,pack2);
1027 cham->RPC_packs = remove(cham->RPC_packs,pack1);
1028 cham->RPC_packs = remove(cham->RPC_packs,pack2);
1029
1030 build_InterChamS_links(cham,pack1);
1031 build_InterChamS_links(cham,pack2);
1032 }
1033 }
1034
1035
1036 void
1037 RPCGeometry::check_map(MultiChamber *cham)
1038 {
1039 int count = 0;
1040
1041 for (int stat=First;stat<=Extended;++stat)
1042 {
1043 int step[2] = {0,0};
1044
1045 for(int i=0;i<2;++i)
1046 {
1047 RPC_pack* pack = cham->firstRPC_z[i][stat];
1048
1049 if(pack)
1050 {
1051 if(pack->prev_in_z[0] || pack->prev_in_z[1])
1052 check_error(cham,"link for previous z pack are not empty!");
1053
1054 while (RPC_pack* next = pack->next_in_z[i])
1055 {
1056 if(pack->station!=next->station)
1057 check_error(cham,"station inconsistence!");
1058
1059 if(!pack->f_proj[i])
1060 check_error(cham,"previous(->next) logic inconsistence!");
1061 if(!next->f_proj[i])
1062 check_error(cham,"next(->previous) logic inconsistence!");
1063
1064 if(next->prev_in_z[i] != pack)
1065 check_error(cham,"next->previous link inconsistence!");
1066
1067 pack = next;
1068 ++step[i];
1069 }
1070 }
1071 }
1072
1073 if(step[0]!=step[1]) check_error(cham,"z index mismatch!");
1074
1075 RPC_pack* pack0 = cham->firstRPC_z[0][stat];
1076 RPC_pack* pack1 = cham->firstRPC_z[1][stat];
1077
1078 if (pack0&&pack1)
1079 {
1080 if(pack0==pack1) ++count;
1081 else count += 2;
1082
1083 while(step[0]--)
1084 {
1085 pack0 = pack0->next_in_z[0];
1086 pack1 = pack1->next_in_z[1];
1087 if(pack0==pack1) ++count;
1088 else count += 2;
1089 }
1090 } else
1091 {
1092 if(pack0)
1093 {
1094 count += step[0] +1;
1095 } else if (pack1)
1096 {
1097 count += step[1] +1;
1098 }
1099 }
1100
1101 }
1102
1103 if(cham->nRPC_packs != count) check_error(cham,"RPC packs lost!");
1104 }
1105
1106
1107 void
1108 RPCGeometry::setuplogic()
1109 {
1110 for(int sect=0;sect<16;++sect)
1111 {
1112 int Jff = sect/2 + 1;
1113
1114 for(int proj=0;proj<2;++proj)
1115 {
1116 for(int stat=First;stat<=Extended;++stat)
1117 {
1118 MultiChamber* cham_sector[3];
1119
1120 cham_sector[0] = 0;
1121 cham_sector[1] = 0;
1122 cham_sector[2] = 0;
1123
1124 MultiChamber* last = 0;
1125
1126 for (int Jtype=1;Jtype<=MaxJtype;++Jtype)
1127 {
1128 for (int Jzz=-Mzz;Jzz<=Mzz;++Jzz)
1129 {
1130 int z_index;
1131 if(Jzz<0)
1132 {
1133 z_index = 0;
1134 last = give_last(cham_sector[z_index]);
1135 } else if (Jzz==0)
1136 {
1137 z_index = 1;
1138 last = give_last(cham_sector[z_index]);
1139 } else
1140 {
1141 z_index = 2;
1142 last = give_last(cham_sector[z_index]);
1143 }
1144
1145 MultiChamber *cham =
1146 T1muonGeo[Jtype-1][(Jff-1)*(2*Mzz+1)+Jzz+Mzz];
1147
1148 if(cham &&
1149 cham->firstRPC_z[proj][stat] &&
1150 sect == cham->physics_sector)
1151 {
1152 if(cham_sector[z_index])
1153 {
1154 last->next = cham;
1155 } else
1156 {
1157 cham_sector[z_index] = cham;
1158 }
1159 last = cham;
1160 }
1161 }
1162 }
1163
1164 build_TransChamS_links(cham_sector,sect,stat,proj);
1165 build_RPC_coding(sect,stat,proj);
1166 }
1167 }
1168 }
1169
1170 if(!s_testbeam) check_logic();
1171
1172 if(DumpLog) dump_logic_sectors();
1173 }
1174
1175
1176 void
1177 RPCGeometry::build_TransChamS_links(MultiChamber **chambers,
1178 int sect,int stat,int proj)
1179 {
1180 if(!chambers[0] && !s_testbeam )
1181 {
1182 if (stat == Extended) return;
1183 printf("RPCGeometry, build_TransChamS_links: no chambers for Jzz < 0\n");
1184 exit(1);
1185 }
1186
1187 if(chambers[1] && (chambers[1])->next)
1188 {
1189 printf("sect = %d\n",sect);
1190 printf("RPCGeometry, build_TransChamS_links: more than 1 chamber");
1191 printf(" whit Jzz = 0\n");
1192 MultiChamber *cham = chambers[1];
1193 do {cham_info(cham,10);} while((cham = cham->next));
1194 exit(1);
1195 }
1196
1197 if(!chambers[2] && !s_testbeam )
1198 {
1199 if (stat == Extended) return;
1200 printf("RPCGeometry, build_TransChamS_links: no chambers for Jzz > 0\n");
1201 exit(1);
1202 }
1203
1204
1205 int neg_sector = sect*2 + proj - 1;
1206 if (neg_sector < 0) neg_sector =31;
1207 int z_index = 0;
1208
1209 MultiChamber *found = find_first_chamber(chambers[z_index],stat,proj);
1210 MultiChamber *last = found;
1211 RPC_pack *prev = 0;
1212 RPC_pack *next = 0;
1213
1214 LogicSector[neg_sector][stat][0] = (found)? found->firstRPC_z[proj][stat]:0;
1215 chambers[z_index] = remove(chambers[z_index],found);
1216
1217 while((found = find_first_chamber(chambers[z_index],stat,proj)))
1218 {
1219 prev = give_last_pack(last,stat,proj);
1220 next = found->firstRPC_z[proj][stat];
1221 prev->next_in_z[proj] = next;
1222 next->prev_in_z[proj] = prev;
1223 chambers[z_index] = remove(chambers[z_index],found);
1224 last = found;
1225 }
1226
1227 LogicSector[neg_sector][stat][1] = (LogicSector[neg_sector][stat][0])?
1228 give_last_pack(last,stat,proj) : 0;
1229
1230 if(chambers[1])
1231 {
1232 RPC_pack *first = (chambers[1])->firstRPC_z[proj][stat];
1233 RPC_pack *last = give_last_pack(chambers[1],stat,proj);
1234 if(first != last)
1235 {
1236 printf("RPCGeometry, build_TransChamS_links: Jzz = 0 packs ");
1237 printf("inconsistence\n");
1238 exit(1);
1239 }
1240 if(LogicSector[neg_sector][stat][1]) (LogicSector[neg_sector][stat][1])->next_in_z[proj] = first;
1241 first->prev_in_z[proj] = LogicSector[neg_sector][stat][1];
1242 LogicSector[neg_sector][stat][1] = last;
1243 (chambers[1])->next = chambers[2];
1244 }
1245
1246 z_index = (chambers[1])? 1 : 2;
1247 int pos_sector = neg_sector + 32;
1248
1249 found = find_first_chamber(chambers[z_index],stat,proj);
1250 last = found;
1251
1252 LogicSector[pos_sector][stat][0] = (found)? found->firstRPC_z[proj][stat]:0;
1253 chambers[z_index] = remove(chambers[z_index],found);
1254
1255 while((found = find_first_chamber(chambers[z_index],stat,proj)))
1256 {
1257 prev = give_last_pack(last,stat,proj);
1258 next = found->firstRPC_z[proj][stat];
1259 prev->next_in_z[proj] = next;
1260 next->prev_in_z[proj] = prev;
1261 chambers[z_index] = remove(chambers[z_index],found);
1262 last = found;
1263 }
1264
1265 LogicSector[pos_sector][stat][1] = (LogicSector[pos_sector][stat][0])?
1266 give_last_pack(last,stat,proj) : 0;
1267
1268 if((next = LogicSector[pos_sector][stat][0]) !=
1269 (prev = LogicSector[neg_sector][stat][1]))
1270 {
1271 if(prev) prev->next_in_z[proj] = next;
1272 if(next) next->prev_in_z[proj] = prev;
1273 }
1274 }
1275
1276
1277 void
1278 RPCGeometry::build_RPC_coding(int sect,int stat,int proj)
1279 {
1280 int sect_neg = sect*2 + proj - 1;
1281 if (sect_neg < 0) sect_neg = 31;
1282 int sect_pos = sect_neg + 32;
1283
1284 RPC_pack *negative_half = LogicSector[sect_neg][stat][1];
1285 RPC_pack *positive_half = LogicSector[sect_pos][stat][0];
1286
1287 bool bothExist = (negative_half && positive_half)? true : false;
1288
1289 int start = (negative_half == positive_half)? 0 : 1;
1290 int neg_index = start;
1291
1292 while(negative_half)
1293 {
1294 int code = neg_index++;
1295
1296 if(negative_half->lvl1_z_index && negative_half->lvl1_z_index != code)
1297 {
1298 printf("RPCGeometry, build_RPC_coding: wrong RPC code for ");
1299 printf("sector logic = %2d, stat = %1d\n",sect_neg,stat);
1300 printf(" old code = %d, new code = %d\n",
1301 negative_half->lvl1_z_index,code);
1302 exit(1);
1303 }
1304
1305 negative_half->lvl1_z_index = code;
1306 negative_half = negative_half->prev_in_z[proj];
1307 }
1308
1309 int pos_index = start;
1310
1311 while(positive_half)
1312 {
1313 int code = pos_index++;
1314
1315 if(positive_half->lvl1_z_index && positive_half->lvl1_z_index != code)
1316 {
1317 printf("RPCGeometry, build_RPC_coding: wrong RPC code for ");
1318 printf("sector logic = %2d, stat = %1d\n",sect_pos,stat);
1319 printf(" old code = %d, new code = %d\n",
1320 positive_half->lvl1_z_index,code);
1321 exit(1);
1322 }
1323
1324 positive_half->lvl1_z_index = code;
1325 positive_half = positive_half->next_in_z[proj];
1326 }
1327
1328 if((pos_index != neg_index) && bothExist)
1329 {
1330 printf("RPCGeometry, build_RPC_coding: different index for the two ");
1331 printf("half barrel\n");
1332 printf(" logic sector = %2d, station = %1d\n",sect_neg,stat);
1333 printf(" negative index = %2d, positive index = %2d\n",neg_index,
1334 pos_index);
1335
1336 }
1337 }
1338
1339
1340 void
1341 RPCGeometry::check_logic(void)
1342 {
1343 for(int logic_sector=0;logic_sector<64;++logic_sector)
1344 {
1345 for(int RPC_station=First;RPC_station<=Third;++RPC_station)
1346 {
1347 for(int i=0;i<2;++i)
1348 {
1349 if(!LogicSector[logic_sector][RPC_station][i])
1350 {
1351 std::string pos = (i)? "end position" : "start position";
1352 printf("RPCGeometry, check_logic: no link for %s of ",
1353 pos.c_str());
1354 printf("logic sector %2d, station = %1d\n",logic_sector,
1355 RPC_station+1);
1356 exit(1);
1357 }
1358 }
1359 }
1360 }
1361 }
1362
1363 void
1364 RPCGeometry::shift_logic(void)
1365 {
1366 int start_logic = 0;
1367
1368 for(int logic_sector=0;logic_sector<64;++logic_sector)
1369 {
1370 if(LogicSector[logic_sector][0][0])
1371 {
1372 start_logic = logic_sector;
1373 break;
1374 }
1375
1376 }
1377
1378 int sector = 0;
1379
1380 for(int logic_sector=start_logic;logic_sector<64;++logic_sector)
1381 {
1382 for(int RPC_station=First;RPC_station<=Extended;++RPC_station)
1383 {
1384 for(int i=0;i<2;++i)
1385 {
1386 LogicSector[sector][RPC_station][i] =
1387 LogicSector[logic_sector][RPC_station][i];
1388 LogicSector[logic_sector][RPC_station][i] = 0;
1389 }
1390 }
1391 ++sector;
1392 }
1393 }
1394
1395
1396
1397
1398
1399 bool
1400 RPCGeometry::RPC_packs_not_equal(RPC_pack* pack1, RPC_pack* pack2) const
1401 {
1402 if(!pack1 || !pack2) return true;
1403 if((pack1->lvl1_z_index == pack2->lvl1_z_index) &&
1404 (pack1->station == pack2->station ) &&
1405 (pack1->nstrip_z == pack2->nstrip_z ) &&
1406 (pack1->stripsize_z == pack2->stripsize_z ) &&
1407 (pack1->nstrip_f == pack2->nstrip_f ) &&
1408 (pack1->stripsize_f == pack2->stripsize_f ) ) return false;
1409 return true;
1410 }
1411
1412
1413 bool
1414 RPCGeometry::Logic_sectors_are_equal(int sector1, int sector2) const
1415 {
1416 int proj1 = (sector1+1)%2;
1417 int proj2 = (sector2+1)%2;
1418
1419 for (int RPC_station=First;RPC_station<=Extended;++RPC_station)
1420 {
1421 RPC_pack* pack1 = (sector1>=32)? LogicSector[sector1][RPC_station][0] :
1422 LogicSector[sector1][RPC_station][1];
1423 RPC_pack* pack2 = (sector2>=32)? LogicSector[sector2][RPC_station][0] :
1424 LogicSector[sector2][RPC_station][1];
1425
1426 while(pack1 && pack2)
1427 {
1428 if(RPC_packs_not_equal(pack1,pack2)) return false;
1429 pack1=(sector1>=32)? pack1->next_in_z[proj1]:pack1->prev_in_z[proj1];
1430 pack2=(sector2>=32)? pack2->next_in_z[proj2]:pack2->prev_in_z[proj2];
1431 }
1432
1433 }
1434
1435 return true;
1436 }
1437
1438
1439 void
1440 RPCGeometry::write_geodata(std::ofstream& file) const
1441 {
1442 int type = 0;
1443 int sector_map[64];
1444 int sector_typ[64];
1445 for (int i=0;i<64;++i) {sector_map[i] = 0; sector_typ[i] = 0;}
1446
1447 for (int sec=0;sec<64;++sec)
1448 {
1449 if (!sector_map[sec])
1450 {
1451 sector_map[sec] = ++type;
1452 sector_typ[type-1] = sec;
1453 for(int i=sec;i<64;++i)
1454 {
1455 if(Logic_sectors_are_equal(sec,i)) sector_map[i] = type;
1456 }
1457 }
1458 }
1459
1460 file << std::endl << "# Start Sector type enumeration" << std::endl;
1461 file << std::endl << "SECTOR TYPES ( 0 - 15 ): ";
1462 for (int i=0;i<16;++i) file << " " << std::setw(2) << sector_map[i];
1463 file << std::endl << "SECTOR TYPES ( 16 - 31 ): ";
1464 for (int i=16;i<32;++i) file << " " << std::setw(2) << sector_map[i];
1465 file << std::endl << std::endl << "SECTOR TYPES ( 32 - 47 ): ";
1466 for (int i=32;i<48;++i) file << " " << std::setw(2) << sector_map[i];
1467 file << std::endl << "SECTOR TYPES ( 48 - 63 ): ";
1468 for (int i=48;i<64;++i) file << " " << std::setw(2) << sector_map[i];
1469 file << std::endl << std::endl;
1470
1471 int dumped_type = 0;
1472 while(dumped_type != type)
1473 {
1474 int sec = sector_typ[dumped_type];
1475 ++dumped_type;
1476 file << std::endl << std::endl << "########################################";
1477 file << "#######################################" << std::endl;
1478 file << "#################################### TYPE " << dumped_type;
1479 file << " ###################################" << std::endl;
1480 file << "########################################";
1481 file << "#######################################" << std::endl;
1482 for (int RPC_station=First;RPC_station<=Extended;++RPC_station)
1483 {
1484 if(!LogicSector[sec][RPC_station][0]) break;
1485
1486 file << std::endl << "RPC GEOM " << std::setw(2) << dumped_type;
1487 file << " : station " << std::setw(1) << RPC_station + 1;
1488 int ncham = 0;
1489 if (sec > 32)
1490 {
1491 ncham = ((LogicSector[sec][RPC_station][0])->lvl1_z_index)?
1492 (LogicSector[sec][RPC_station][1])->lvl1_z_index :
1493 (LogicSector[sec][RPC_station][1])->lvl1_z_index + 1;
1494 } else
1495 {
1496 ncham = ((LogicSector[sec][RPC_station][1])->lvl1_z_index)?
1497 (LogicSector[sec][RPC_station][0])->lvl1_z_index :
1498 (LogicSector[sec][RPC_station][0])->lvl1_z_index + 1;
1499 }
1500 file << " made of " << std::setw(2) << ncham << " chamber. ";
1501 file << "Strips in connectors: eta ";
1502 int eta_conn = (RPC_station==Second)? 8 : 4;
1503 file << std::setw(1) << eta_conn << " phi 8" << std::endl;
1504
1505 file << "{" << std::endl;
1506 RPC_pack* pack = (sec>=32)? LogicSector[sec][RPC_station][0] :
1507 LogicSector[sec][RPC_station][1];
1508 int proj = (sec+1)%2;
1509
1510 while(pack)
1511 {
1512 file << " cham " << std::setw(2) << pack->lvl1_z_index;
1513 file << " type " << pack->cham->name << "-";
1514 file << " eta_side " << std::setw(2) << pack->nstrip_z;
1515 file << " " << std::setw(2) << pack->nstrip_z/eta_conn << " 2";
1516 file << " phi_side " << std::setw(2) << pack->nstrip_f;
1517 file << " " << std::setw(2) << pack->nstrip_f/8 << " 2";
1518 file << " #MDT chamber " << abs(pack->cham->amdb_z_index)
1519 << std::endl;
1520
1521 pack=(sec>=32)? pack->next_in_z[proj]:pack->prev_in_z[proj];
1522 }
1523
1524 file << "}" << std::endl;
1525 }
1526 }
1527
1528
1529 }
1530
1531
1532
1533
1534
1535 const encodeProj*
1536 RPCGeometry::find_proj(int Jobj) const
1537 {
1538 encodeProj* tmp = projMap;
1539 while(tmp)
1540 {
1541 if (tmp->Jobj == Jobj) return tmp;
1542 tmp = tmp->next;
1543 }
1544
1545 return 0;
1546 }
1547
1548 void
1549 RPCGeometry::add_proj(encodeProj* item)
1550 {
1551 if (!projMap)
1552 {
1553 projMap = item;
1554 return;
1555 }
1556
1557 encodeProj* tmp = projMap;
1558 while(tmp)
1559 {
1560 if(tmp->Jobj==item->Jobj && tmp->proj==tmp->proj) return;
1561 if(tmp->Jobj==item->Jobj && tmp->proj!=tmp->proj)
1562 {
1563 std::cout << "RPCGeometry ERROR in projMap:"
1564 << " try to insert Jobj=" << item->Jobj
1565 << " with proj=" << item->proj << " in" << std::endl;
1566 dump_projMap();
1567 exit(1);
1568 }
1569 if (tmp->next == 0)
1570 {
1571 tmp->next = item;
1572 return;
1573 }
1574 tmp = tmp->next;
1575 }
1576 }
1577
1578
1579
1580
1581
1582
1583
1584
1585
1586
1587
1588
1589
1590
1591
1592
1593
1594
1595
1596
1597
1598
1599
1600 void
1601 RPCGeometry::dump_projMap()
1602 {
1603 encodeProj* tmp = projMap;
1604 while(tmp)
1605 {
1606 std::cout << "RPCGeometry projMap: Jobj=" << tmp->Jobj
1607 << ", proj=" << tmp->proj << std::endl;
1608 tmp = tmp->next;
1609 }
1610 }
1611
1612
1613
1614
1615
1616 #ifndef LVL1_STANDALONE
1617 bool
1618 RPCGeometry::give_strip_index(unsigned long int code,int& Jtype,
1619 int& StationEta, int& StationPhi, int& DoubletR,
1620 int& DoubletZ, int& DoubletP, int& GasGap,
1621 int& MeasuresPhi,int& Strip) const
1622 {
1623 AMDC;
1624
1625 int Jff;
1626 int Jzz;
1627 int Job;
1628 int Jspli;
1629 int Jsl;
1630 int Jsz;
1631 int Jstri;
1632
1633 RPCdecoder decode(code);
1634
1635 if (decode)
1636 {
1637 int z_index = decode.rpc_z_index();
1638 int lvl1_station = decode.lvl1_station() - 1;
1639 int logic_sector = decode.logic_sector();
1640 int proj = (logic_sector+1)%2;
1641
1642 RPC_pack *pack = LogicSector[logic_sector][lvl1_station][0];
1643 while(pack->lvl1_z_index != z_index) pack = pack->next_in_z[proj];
1644
1645 Jff = (int)((pack->cham->physics_sector)/2) + 1;
1646 Jzz = pack->cham->amdb_z_index;
1647 Job = pack->Jobj;
1648 Jspli = (pack->Readout_f==1)? 1 : proj + 1;
1649 Jsl = decode.rpc_layer() + 1;
1650 Jsz = decode.strip_type();
1651 Jstri = decode.strip_number() + 1;
1652
1653 std::string StationName = pack->cham->name;
1654 Jtype = amdc->GetJtyp(StationName);
1655
1656 amdc->GetRPCindexFromAMDC(StationName,Jff,Jzz,Job,Jspli,Jsl,Jsz,Jstri,
1657 StationEta,StationPhi,DoubletR,DoubletZ,
1658 DoubletP,GasGap,MeasuresPhi,Strip);
1659
1660 }
1661 else return false;
1662
1663 return true;
1664 }
1665
1666
1667 unsigned long int
1668 RPCGeometry::give_strip_code( std::string StationName, int StationEta,
1669 int StationPhi, int DoubletR, int DoubletZ,
1670 int DoubletP, int GasGap, int MeasuresPhi,
1671 int Strip) const
1672 {
1673 AMDC;
1674
1675 int Jtype = amdc->GetJtyp(StationName);
1676 int Jff;
1677 int Jzz;
1678 int Job;
1679 int Jspli;
1680 int Jsl;
1681 int Jsz;
1682 int Jstri;
1683
1684
1685
1686
1687
1688
1689
1690
1691 amdc->GetAMDCindexFromRPC(StationName,StationEta,StationPhi,DoubletR,
1692 DoubletZ,DoubletP,GasGap,MeasuresPhi,Strip,
1693 Jff,Jzz,Job,Jspli,Jsl,Jsz,Jstri);
1694
1695 if ( Jtype==4 && Jff!=6 && Jff!=7 && Job>=s_JobMinForEncodeProj &&
1696 (Jzz==2 || Jzz==4 || Jzz==-2 || Jzz==-4) )
1697 {
1698 const encodeProj* tmp = find_proj(Job);
1699 if (!tmp)
1700 {
1701 std::cout << "RPCGeometry ERROR: can't find proj index"
1702 << std::endl;
1703 exit(1);
1704 }
1705 else Jspli = tmp->proj;
1706
1707 }
1708
1709 return this->give_strip_code(Jtype,Jff,Jzz,Job,Jsl,Jsz,Jstri,Jspli);
1710 }
1711 #endif
1712
1713 unsigned long int
1714 RPCGeometry::give_strip_code(int Jtype,int Jff,int Jzz,int Jobj,int Jlay,
1715 int Jsty,int Jstr,int proj) const
1716 {
1717
1718
1719
1720 MultiChamber *cham = T1muonGeo[Jtype-1][(Jff-1)*(2*Mzz+1)+Jzz+Mzz];
1721
1722
1723
1724
1725
1726
1727
1728 RPC_pack *pack = 0;
1729
1730 for(int stat=First;stat<=Extended;++stat)
1731 {
1732 RPC_pack *tmp = cham->firstRPC_z[proj-1][stat];
1733
1734 while(tmp && tmp->cham == cham)
1735 {
1736
1737 if(tmp->Jobj == Jobj)
1738 {
1739 if(pack)
1740 {
1741 printf("RPCGeometry, give_strip_code: ");
1742 printf("more than 1 pack found!\n");
1743 return 0;
1744 } else
1745 {
1746 pack = tmp;
1747 }
1748 }
1749 tmp = tmp->next_in_z[proj-1];
1750 }
1751 }
1752
1753 if(!pack)
1754 {
1755 proj++;
1756 for(int stat=First;stat<=Extended;++stat)
1757 {
1758 RPC_pack *tmp = cham->firstRPC_z[proj-1][stat];
1759 while(tmp && tmp->cham == cham)
1760 {
1761
1762 if(tmp->Jobj == Jobj)
1763 {
1764 if(pack)
1765 {
1766 printf("RPCGeometry, give_strip_code: ");
1767 printf("more than 1 pack found!\n");
1768 return 0;
1769 } else
1770 {
1771 pack = tmp;
1772 }
1773 }
1774 tmp = tmp->next_in_z[proj-1];
1775 }
1776 }
1777 }
1778
1779 int station = pack->station;
1780
1781 int z_index = pack->lvl1_z_index;
1782
1783 int logic_sector = cham->physics_sector*2 + proj - 2;
1784 if (logic_sector < 0 ) logic_sector = 31;
1785
1786 int nstrip = (Jsty-1)? pack->nstrip_z : pack->nstrip_f;
1787
1788 if (Jstr > nstrip)
1789 {
1790 printf("RPCGeometry, give_strip_code: ");
1791 printf("inconsistence on strip number!\n");
1792 printf(" RPC strips = %2d, digitized strip number = %2d\n",
1793 nstrip,Jstr);
1794 return 0;
1795 }
1796
1797
1798
1799 if (Jzz >= 0) logic_sector += 32;
1800
1801
1802
1803 RPC_pack *check = LogicSector[logic_sector][station][0];
1804
1805 while (check)
1806 {
1807 if(check->lvl1_z_index == z_index) break;
1808 check = check->next_in_z[proj-1];
1809 }
1810
1811 if (check != pack)
1812 {
1813 printf("RPCGeometry, give_strip_code: ");
1814 printf("RPC pack inconsistence in Sector Logic map!\n");
1815 exit(1);
1816 return 0;
1817 }
1818
1819 RPCdecoder decode(Jsty,logic_sector,station+1,Jlay-1,z_index,Jstr-1);
1820
1821
1822
1823
1824
1825
1826
1827 return decode.code();
1828 }
1829
1830
1831 bool
1832 RPCGeometry::give_strip_coordinates(unsigned long int code, float *xyz) const
1833 {
1834 RPCdecoder decode(code);
1835 if (!decode) return false;
1836
1837 int strip_number = decode.strip_number();
1838 int z_index = decode.rpc_z_index();
1839 int rpc_layer = decode.rpc_layer();
1840 int lvl1_station = decode.lvl1_station() - 1;
1841 int logic_sector = decode.logic_sector();
1842 int Jsty = decode.strip_type();
1843 int proj = (logic_sector+1)%2;
1844
1845
1846
1847
1848
1849
1850
1851
1852
1853
1854 RPC_pack *pack = LogicSector[logic_sector][lvl1_station][0];
1855 while(pack->lvl1_z_index != z_index) pack = pack->next_in_z[proj];
1856
1857 if(Jsty == 1)
1858 {
1859 xyz[0] = *(pack->f_proj[proj] + strip_number*4 + rpc_layer*2);
1860 xyz[1] = *(pack->f_proj[proj] + strip_number*4 + rpc_layer*2 + 1);
1861 xyz[2] = pack->pack_coo[2];
1862 } else
1863 {
1864
1865
1866
1867 xyz[2] = *(pack->z_proj + strip_number);
1868 float DeltaZ = pack->pack_ref_coo[2][rpc_layer][0] - xyz[2];
1869 xyz[0] = pack->pack_ref_coo[0][rpc_layer][0] -
1870 DeltaZ*pack->pack_ref_deltas[0][rpc_layer];
1871 xyz[1] = pack->pack_ref_coo[1][rpc_layer][0] -
1872 DeltaZ*pack->pack_ref_deltas[1][rpc_layer];
1873 }
1874
1875 return true;
1876 }
1877
1878
1879 bool
1880 RPCGeometry::give_strip_radius(unsigned long int code,float& radius) const
1881 {
1882 RPCdecoder decode(code);
1883 if(!decode) return false;
1884
1885 int z_index = decode.rpc_z_index();
1886 int rpc_layer = decode.rpc_layer();
1887 int lvl1_station = decode.lvl1_station() - 1;
1888 int logic_sector = decode.logic_sector();
1889 int proj = (logic_sector+1)%2;
1890
1891 RPC_pack *pack = LogicSector[logic_sector][lvl1_station][0];
1892 while(pack->lvl1_z_index != z_index) pack = pack->next_in_z[proj];
1893
1894 radius = pack->gas_layer_radius[rpc_layer];
1895 return true;
1896 }
1897
1898 bool
1899 RPCGeometry::give_eta_strip_size(unsigned int code,float& size) const
1900 {
1901 RPCdecoder decode(code);
1902 if(!decode) return false;
1903
1904 int z_index = decode.rpc_z_index();
1905 int lvl1_station = decode.lvl1_station() - 1;
1906 int logic_sector = decode.logic_sector();
1907 int proj = (logic_sector+1)%2;
1908
1909 RPC_pack *pack = LogicSector[logic_sector][lvl1_station][0];
1910 while(pack->lvl1_z_index != z_index) pack = pack->next_in_z[proj];
1911
1912 size = pack->stripsize_z;
1913 return true;
1914 }
1915
1916 bool
1917 RPCGeometry::give_phi_strip_size(unsigned int code,float& size) const
1918 {
1919 RPCdecoder decode(code);
1920 if(!decode) return false;
1921
1922 int z_index = decode.rpc_z_index();
1923 int lvl1_station = decode.lvl1_station() - 1;
1924 int logic_sector = decode.logic_sector();
1925 int proj = (logic_sector+1)%2;
1926
1927 RPC_pack *pack = LogicSector[logic_sector][lvl1_station][0];
1928 while(pack->lvl1_z_index != z_index) pack = pack->next_in_z[proj];
1929
1930 size = pack->stripsize_f;
1931 return true;
1932 }
1933
1934
1935 bool
1936 RPCGeometry::give_station_radius(unsigned int code,float& radius) const
1937 {
1938 RPCdecoder decode(code);
1939 if(!decode) return false;
1940
1941 int z_index = decode.rpc_z_index();
1942 int rpc_layer = decode.rpc_layer();
1943 int lvl1_station = decode.lvl1_station() - 1;
1944 int logic_sector = decode.logic_sector();
1945 int proj = (logic_sector+1)%2;
1946
1947 RPC_pack *pack = LogicSector[logic_sector][lvl1_station][0];
1948 while(pack->lvl1_z_index != z_index) pack = pack->next_in_z[proj];
1949
1950 radius = pack->gas_layer_radius[rpc_layer];
1951
1952 return true;
1953 }
1954
1955
1956 bool
1957 RPCGeometry::give_station_phi(unsigned int code,float& phi) const
1958 {
1959 RPCdecoder decode(code);
1960 if(!decode) return false;
1961
1962 int z_index = decode.rpc_z_index();
1963 int lvl1_station = decode.lvl1_station() - 1;
1964 int logic_sector = decode.logic_sector();
1965 int proj = (logic_sector+1)%2;
1966
1967 RPC_pack *pack = LogicSector[logic_sector][lvl1_station][0];
1968 while(pack->lvl1_z_index != z_index) pack = pack->next_in_z[proj];
1969
1970 phi = pack->rotation_angle;
1971
1972
1973
1974
1975
1976 return true;
1977 }
1978
1979 bool
1980 RPCGeometry::give_rpc_boundaries(int logic_sector,int station,int z_index,
1981 float& zmin, float& zmax) const
1982 {
1983 int proj = (logic_sector+1)%2;
1984
1985 RPC_pack* pack = LogicSector[logic_sector][station-1][0];
1986 while(pack->lvl1_z_index != z_index) pack = pack->next_in_z[proj];
1987
1988 zmin = fabsf(pack->pack_coo[2] - pack->half_dim[1]);
1989 zmax = fabsf(pack->pack_coo[2] + pack->half_dim[1]);
1990
1991 if (zmin > zmax)
1992 {
1993 float tmp = zmax;
1994 zmax = zmin;
1995 zmin = tmp;
1996 }
1997
1998 return true;
1999 }
2000
2001 bool
2002 RPCGeometry::isFull() const
2003 {
2004 return !s_testbeam;
2005 }
2006
2007
2008
2009
2010
2011
2012
2013
2014 void
2015 RPCGeometry::buildgeo()
2016 {
2017
2018
2019
2020
2021
2022
2023
2024
2025 AMDC;
2026
2027 Mtype = amdc->Mtyp();
2028 Mff = amdc->Mphi();
2029 Mzz = amdc->Mzz();
2030
2031 Mtype = amdc->StationTypeEnd();
2032 Mff = amdc->PositionPhiEnd();
2033 Mzz = amdc->PositionZEnd();
2034
2035
2036
2037 for (int Jtype=0;Jtype<MaxJtype;++Jtype)
2038 {
2039 for (int J=0;J<Mff*(2*Mzz+1);++J)
2040 {
2041 T1muonGeo[Jtype][J] = 0;
2042 }
2043 }
2044
2045
2046 for (int Jtype=1;Jtype<=MaxJtype;++Jtype)
2047 {
2048 std::string name = amdc->GetStationType(Jtype);
2049 int MJzz = amdc->PositionZEnd();
2050
2051
2052 int IsValid;
2053 int MJgeo = amdc->StationGeomEnd();
2054 double DimLocX,DimLocY,DimLocZ,CenterLocX,CenterLocY,CenterLocZ;
2055
2056 for(int Jgeo = amdc->StationGeomBegin();Jgeo < MJgeo;++Jgeo) {
2057
2058 int Jcut = 0;
2059 amdc->GetStationDimensions(name,Jgeo,Jcut,IsValid,
2060 DimLocX,DimLocY,DimLocZ,
2061 CenterLocX,CenterLocY,CenterLocZ );
2062
2063 if (IsValid > 0){
2064
2065 std::string GeomTechnoName;
2066 int GeomIsValid,GeomTechnoIndex;
2067 int GeomSplitX,GeomSplitY,GeomShapeIndex;
2068 double GeomDx,GeomDy,GeomDz,GeomWs,GeomWl;
2069 double GeomLe,GeomEx,GeomD1,GeomD2,GeomD3;
2070
2071 int NJobj = amdc->GetNumberOfStationObjects(name,Jgeo);
2072
2073
2074
2075 enum tech {SPA,MDT,RPC};
2076 int exist_tech[3] = {0,0,0};
2077
2078 for(int Jobj=1;Jobj<=NJobj;Jobj++) {
2079
2080 int Itec = amdc->ITEC(Jtype,Jgeo,Jobj);
2081 std::string tech_name = amdc->TechnologyName(Itec);
2082
2083 if (tech_name=="SPA") exist_tech[SPA] = Jobj;
2084 if (tech_name=="CHV") exist_tech[SPA] = Jobj;
2085 if (tech_name=="MDT") exist_tech[MDT] = Jobj;
2086 if (tech_name=="RPC") exist_tech[RPC] = Jobj;
2087 }
2088
2089 int JobjMid = 0;
2090 int IstaMid = 0;
2091 int ItecMid = 0;
2092 double widthSpa = 0;
2093
2094
2095 if (exist_tech[SPA]) {
2096 JobjMid = exist_tech[SPA];
2097 IstaMid = amdc->ISTA(Jtype,Jgeo,JobjMid);
2098 ItecMid = amdc->ITEC(Jtype,Jgeo,JobjMid);
2099
2100 std::string name = amdc->TechnologyName(ItecMid);
2101
2102 widthSpa = (name == "CHV")? (amdc->STAPP(ItecMid,IstaMid))/2. :
2103 (amdc->STAEE(ItecMid,IstaMid))/2.;
2104 } else if (exist_tech[MDT]) {
2105 JobjMid = exist_tech[MDT];
2106 IstaMid = amdc->ISTA(Jtype,Jgeo,JobjMid);
2107 ItecMid = amdc->ITEC(Jtype,Jgeo,JobjMid);
2108 widthSpa = (amdc->STAEE(ItecMid,IstaMid))/2.;
2109 } else if (exist_tech[RPC]) {
2110 JobjMid = exist_tech[RPC];
2111 IstaMid = amdc->ISTA(Jtype,Jgeo,JobjMid);
2112 ItecMid = amdc->ITEC(Jtype,Jgeo,JobjMid);
2113 widthSpa = (amdc->STATT(ItecMid,IstaMid,9))/2.;
2114 }
2115
2116
2117 amdc->GetStationObjectParam( name,Jgeo,JobjMid,GeomIsValid,
2118 GeomTechnoName,GeomTechnoIndex,
2119 GeomSplitX,GeomSplitY,
2120 GeomShapeIndex,GeomDx,GeomDy,
2121 GeomDz,GeomWs,GeomWl,GeomLe,
2122 GeomEx,GeomD1,GeomD2,GeomD3 );
2123
2124 double RmidOffset = GeomDz + widthSpa;
2125
2126
2127 if (!RmidOffset) {
2128 printf("RPCGeometry: error in computing RmidOffset!\n");
2129 exit(1);
2130 }
2131
2132
2133
2134 for (int Jff=amdc->PositionPhiBegin();Jff<=Mff;++Jff) {
2135 for (int Jzz=amdc->PositionZBegin();Jzz<=MJzz;++Jzz) {
2136
2137 int PosiIsValid,PosiJgeo,PosiJcut,PosiIsBarrel;
2138 double PosiPhi,PosiZ,PosiR,PosiS;
2139 double PosiAlfa, PosiBeta,PosiGamma;
2140
2141 amdc->GetStationPositionParam(name,Jff,Jzz,PosiIsValid,
2142 PosiJgeo,PosiJcut,PosiIsBarrel,
2143 PosiPhi,PosiZ,PosiR,PosiS,
2144 PosiAlfa,PosiBeta,PosiGamma);
2145
2146 if( PosiIsValid > 0 && PosiJgeo == Jgeo ) {
2147
2148
2149
2150 MultiChamber* cham = 0;
2151 RPC_pack* last_pack = 0;
2152
2153 double Rmid = nor(PosiR + RmidOffset);
2154
2155 for(int Jobj=1;Jobj<=NJobj;Jobj++) {
2156
2157 amdc->GetStationObjectParam( name,Jgeo,Jobj,GeomIsValid,
2158 GeomTechnoName,GeomTechnoIndex,
2159 GeomSplitX,GeomSplitY,
2160 GeomShapeIndex,GeomDx,GeomDy,
2161 GeomDz,GeomWs,GeomWl,GeomLe,
2162 GeomEx,GeomD1,GeomD2,GeomD3 );
2163
2164 int Ista = amdc->ISTA(Jtype,Jgeo,Jobj);
2165 int Itec = amdc->ITEC(Jtype,Jgeo,Jobj);
2166 std::string tech_name = amdc->TechnologyName(Itec);
2167
2168 if(tech_name == "RPC")
2169 {
2170
2171
2172
2173 if (!cham)
2174 {
2175 cham = new (MultiChamber);
2176 cham->Rmid = Rmid;
2177 cham->name = name;
2178 cham->type = Jgeo;
2179 cham->amdb_z_index = Jzz;
2180 cham->physics_sector =(name.c_str()[2] == 'L')?
2181 (Jff-1)*2 : (Jff-1)*2+1;
2182 cham->nRPC_packs = 0;
2183 cham->RPC_packs = 0;
2184 cham->next = 0;
2185 for(int i=0;i<2;++i)for(int j=First;j<=Extended;++j)
2186 cham->firstRPC_z[i][j] = 0;
2187 T1muonGeo[Jtype-1][(Jff-1)*(2*Mzz+1)+Jzz+Mzz] =
2188 cham;
2189 }
2190
2191 RPC_pack* pack = new (RPC_pack);
2192 pack->Jobj = Jobj;
2193 pack->Readout_f = 0;
2194 pack->cham = cham;
2195 pack->next = 0;
2196 pack->z_proj = 0;
2197 pack->lvl1_z_index = 0;
2198 for(int i=0;i<2;++i)
2199 {
2200 pack->f_proj[i] = 0;
2201 pack->next_in_z[i] = 0;
2202 pack->prev_in_z[i] = 0;
2203 for(int j=0;j<2;++j) for(int k=0;k<2;++k)
2204 pack->gas_layer_XY[i][j][k] = 0.;
2205 }
2206
2207
2208 pack->stripsize_f = nor(amdc->STATT(Itec,Ista,10));
2209 pack->stripsize_z = nor(amdc->STATT(Itec,Ista,11));
2210
2211 float Readout_f = amdc->STAOO(Itec,Ista,2);
2212
2213
2214 float Readout_z = 1.;
2215 float Length = nor(GeomLe);
2216 float Width = nor(GeomWs);
2217 float Dead = nor(GeomD1);
2218
2219 pack->Readout_f = (int) Readout_f;
2220
2221 pack->nstrip_z=(int)((Length/Readout_z - Dead)/
2222 pack->stripsize_z);
2223 pack->nstrip_f=(int)((Width/Readout_f - Dead)/
2224 pack->stripsize_f);
2225
2226
2227 pack->half_dim[0] = Width/2.;
2228 pack->half_dim[1] = Length/2.;
2229 pack->half_dim[2] = nor(amdc->STATT(Itec,Ista,9))/2.;
2230
2231
2232 float shift_rad = nor(GeomDz);
2233 float shift_z = nor(GeomDy);
2234 float shift_s = nor(GeomDx);
2235 float orto_rad = nor(PosiS);
2236 bool isOrtoRad = (orto_rad>15.) ? true : false;
2237 float x = nor(PosiR)+ shift_rad + pack->half_dim[2];
2238 float y = shift_s + orto_rad;
2239 float z = nor(PosiZ)+ shift_z + pack->half_dim[1];
2240 float phi = PosiPhi;
2241
2242
2243 float tmp[3] = {x,y,z};
2244
2245 rotate(nor(PosiZ),nor(PosiR),orto_rad,
2246 PosiAlfa,PosiBeta,PosiGamma,tmp);
2247
2248 pack->rotation_angle = phi;
2249
2250 x = tmp[0];
2251 y = tmp[1];
2252 z = tmp[2];
2253
2254 pack->pack_coo[0] = x*cos(phi) - y*sin(phi);
2255 pack->pack_coo[1] = x*sin(phi) + y*cos(phi);
2256 pack->pack_coo[2] = z;
2257
2258
2259
2260 if (name.c_str()[1] == 'O')
2261 {
2262 pack->station = Third;
2263 if (name.c_str()[2] == 'F' ||
2264 name.c_str()[2] == 'G' )
2265 {
2266 pack->station = (x<Rmid)? Third : Extended;
2267 }
2268 } else
2269 {
2270 pack->station = (x<Rmid)? First : Second;
2271 if( !exist_tech[SPA] && !exist_tech[MDT])
2272 pack->station = Second;
2273
2274
2275 }
2276
2277
2278
2279 float rpc_thickness= nor(amdc->STATT(Itec,Ista,4));
2280 float low_honeycomb= nor(amdc->STATT(Itec,Ista,7));
2281 float cen_honeycomb= nor(amdc->STATT(Itec,Ista,5));
2282 float gas_volume = nor(amdc->STATT(Itec,Ista,13));
2283 float bakelite = nor(amdc->STATT(Itec,Ista,14));
2284 float gas_gap_pet = nor(amdc->STAOO(Itec,Ista,7));
2285 float gas_gap_air = nor(amdc->STAOO(Itec,Ista,8));
2286 float str_thickness= nor(amdc->STATT(Itec,Ista,15));
2287 float str_support = nor(amdc->STATT(Itec,Ista,16));
2288 float str_pet_upper= nor(amdc->STAOO(Itec,Ista,9));
2289 float str_pet_lower= nor(amdc->STAOO(Itec,Ista,10));
2290
2291
2292
2293
2294
2295
2296
2297
2298
2299
2300
2301
2302
2303 float gas_center = str_thickness +
2304 str_support*2 +
2305 str_pet_upper +
2306 str_pet_lower +
2307 gas_gap_pet +
2308 gas_gap_air +
2309 bakelite +
2310 gas_volume/2. ;
2311
2312 float mid_layer_1 = low_honeycomb +
2313 gas_center;
2314
2315 float mid_layer_2 = low_honeycomb +
2316 rpc_thickness +
2317 cen_honeycomb +
2318 gas_center;
2319
2320 float radius_lay1 = nor(PosiR) + shift_rad +
2321 mid_layer_1;
2322 float radius_lay2 = nor(PosiR) + shift_rad +
2323 mid_layer_2;
2324
2325 if(amdc->ISHAPE(Jtype,Jgeo,Jobj) == -1)
2326 {
2327 radius_lay1 -= low_honeycomb;
2328 radius_lay2 -= low_honeycomb;
2329 }
2330
2331 pack->gas_layer_radius[0] = radius_lay1;
2332 pack->gas_layer_radius[1] = radius_lay2;
2333
2334
2335
2336 float Ref1[3];
2337 Ref1[0] = radius_lay1;
2338 Ref1[1] = shift_s + orto_rad;
2339 Ref1[2] = nor(PosiZ) + shift_z + pack->half_dim[1];
2340
2341 float Ref2[3];
2342 Ref2[0] = radius_lay1;
2343 Ref2[1] = shift_s + orto_rad;
2344 Ref2[2] = nor(PosiZ) + shift_z;
2345
2346 rotate(nor(PosiZ),nor(PosiR),orto_rad,
2347 PosiAlfa,PosiBeta,PosiGamma,Ref1);
2348
2349 rotate(nor(PosiZ),nor(PosiR),orto_rad,
2350 PosiAlfa,PosiBeta,PosiGamma,Ref2);
2351
2352 for (int nco=0;nco<3;++nco)
2353 {
2354 pack->pack_ref_coo[nco][0][0] = Ref2[nco];
2355 pack->pack_ref_coo[nco][0][1] = Ref1[nco];
2356 }
2357
2358 float DeltaX = pack->pack_ref_coo[0][0][0] -
2359 pack->pack_ref_coo[0][0][1];
2360 float DeltaY = pack->pack_ref_coo[1][0][0] -
2361 pack->pack_ref_coo[1][0][1];
2362 float DeltaZ = pack->pack_ref_coo[2][0][0] -
2363 pack->pack_ref_coo[2][0][1];
2364
2365 pack->pack_ref_deltas[0][0] = DeltaX/DeltaZ;
2366 pack->pack_ref_deltas[1][0] = DeltaY/DeltaZ;
2367
2368
2369
2370 Ref1[0] = radius_lay2;
2371 Ref1[1] = shift_s + orto_rad;
2372 Ref1[2] = nor(PosiZ) + shift_z + pack->half_dim[1];
2373
2374 Ref2[0] = radius_lay2;
2375 Ref2[1] = shift_s + orto_rad;
2376 Ref2[2] = nor(PosiZ) + shift_z;
2377
2378 rotate(nor(PosiZ),nor(PosiR),orto_rad,
2379 PosiAlfa,PosiBeta,PosiGamma,Ref1);
2380
2381 rotate(nor(PosiZ),nor(PosiR),orto_rad,
2382 PosiAlfa,PosiBeta,PosiGamma,Ref2);
2383
2384 for (int nco=0;nco<3;++nco)
2385 {
2386 pack->pack_ref_coo[nco][1][0] = Ref2[nco];
2387 pack->pack_ref_coo[nco][1][1] = Ref1[nco];
2388 }
2389
2390 DeltaX = pack->pack_ref_coo[0][1][0] -
2391 pack->pack_ref_coo[0][1][1];
2392 DeltaY = pack->pack_ref_coo[1][1][0] -
2393 pack->pack_ref_coo[1][1][1];
2394 DeltaZ = pack->pack_ref_coo[2][1][0] -
2395 pack->pack_ref_coo[2][1][1];
2396
2397 pack->pack_ref_deltas[0][1] = DeltaX/DeltaZ;
2398 pack->pack_ref_deltas[1][1] = DeltaY/DeltaZ;
2399
2400 float Zp = nor(PosiZ);
2401 float Rp = nor(PosiR);
2402 float Sp = orto_rad;
2403 float Aa = PosiAlfa;
2404 float Ba = PosiBeta;
2405 float Ga = PosiGamma;
2406
2407 if((int)Readout_f == 2)
2408 {
2409 for(int i=0;i<2;++i) for(int j=0;j<2;++j)
2410 {
2411 float sy = (i)? pack->half_dim[0]/2. :
2412 -pack->half_dim[0]/2.;
2413 float sx = pack->gas_layer_radius[j];
2414 float z = nor(PosiZ)+ shift_z +
2415 pack->half_dim[1];
2416
2417 float Ce[3] = {sx,sy,z};
2418 rotate(Zp,Rp,Sp,Aa,Ba,Ga,Ce);
2419 float X = Ce[0]*cos(phi) - Ce[1]*sin(phi);
2420 float Y = Ce[0]*sin(phi) + Ce[1]*cos(phi);
2421 pack->gas_layer_XY[i][j][0] = X;
2422 pack->gas_layer_XY[i][j][1] = Y;
2423 }
2424 } else
2425 {
2426
2427
2428
2429
2430 if (Jtype!=4 || Jff==6 || Jff==7 ||
2431 !(Jzz==-2 || Jzz==-4 || Jzz==2 || Jzz==4))
2432 {
2433 std::cout << "RPCGeometry WARNING:"
2434 << " ambiguities in data decoding"
2435 << std::endl;
2436 }
2437 if(Jobj <= s_JobMinForEncodeProj)
2438 s_JobMinForEncodeProj = Jobj;
2439 encodeProj* tmp = new (encodeProj);
2440 tmp->Jobj = Jobj;
2441 tmp->proj = (shift_s<=0)? 1 : 2;
2442 tmp->next = 0;
2443 add_proj(tmp);
2444
2445
2446 int i = (shift_s<0||(orto_rad<0 && isOrtoRad))? 0:1;
2447 for(int j=0;j<2;++j)
2448 {
2449 float sy = shift_s + orto_rad;
2450 float sx = pack->gas_layer_radius[j];
2451 float z = nor(PosiZ)+ shift_z +
2452 pack->half_dim[1];
2453
2454
2455 float Ce[3] = {sx,sy,z};
2456 rotate(Zp,Rp,Sp,Aa,Ba,Ga,Ce);
2457 float X = Ce[0]*cos(phi) - Ce[1]*sin(phi);
2458 float Y = Ce[0]*sin(phi) + Ce[1]*cos(phi);
2459 pack->gas_layer_XY[i][j][0] = X;
2460 pack->gas_layer_XY[i][j][1] = Y;
2461 }
2462 }
2463
2464 float *str_z = new float [pack->nstrip_z];
2465
2466 float sign_z = (Jzz)? (Jzz)/abs(Jzz) : +1;
2467
2468 float start_position = (sign_z > 0)?
2469 nor(PosiZ) + shift_z:
2470 nor(PosiZ) + shift_z+ Length;
2471
2472 float offset = (pack->half_dim[1]*2 -
2473 pack->stripsize_z * pack->nstrip_z)/2.;
2474 float start_z = start_position + (offset +
2475 pack->stripsize_z/2.)*sign_z;
2476
2477
2478
2479
2480 for (int i=0;i<pack->nstrip_z;++i)
2481 {
2482 float x = pack->gas_layer_radius[0];
2483 float y = shift_s + orto_rad;
2484 float z = start_z + i*pack->stripsize_z*sign_z;
2485 float Ce[3] = {pack->gas_layer_radius[0],0,z};
2486 rotate(Zp,Rp,Sp,Aa,Ba,Ga,Ce);
2487 *(str_z+i) = Ce[2];
2488 }
2489
2490 pack->z_proj = str_z;
2491
2492
2493
2494
2495
2496
2497 if ((int)Readout_f == 2)
2498 {
2499 float* str_f1 = new float [4*pack->nstrip_f];
2500 float* str_f2 = new float [4*pack->nstrip_f];
2501
2502 float offset = (pack->half_dim[0] -
2503 pack->stripsize_f * pack->nstrip_f)/2.;
2504
2505
2506
2507
2508
2509 float start_y = offset + shift_s + orto_rad +
2510 pack->stripsize_f/2.;
2511
2512 for (int i=0;i<pack->nstrip_f;++i)
2513 for(int r=0;r<2;++r)
2514 {
2515
2516
2517 float sy = start_y+i*pack->stripsize_f;
2518 float sx = pack->gas_layer_radius[r];
2519 float z = nor(PosiZ)+ shift_z +
2520 pack->half_dim[1];
2521
2522 float Ce[3] = {sx,sy,z};
2523 rotate(Zp,Rp,Sp,Aa,Ba,Ga,Ce);
2524 float X = Ce[0]*cos(phi) - Ce[1]*sin(phi);
2525 float Y = Ce[0]*sin(phi) + Ce[1]*cos(phi);
2526 *(str_f1+i*4+r*2) = X;
2527 *(str_f1+i*4+r*2+1) = Y;
2528 }
2529
2530
2531
2532
2533
2534
2535 start_y = - pack->half_dim[0] + offset
2536 + shift_s + orto_rad +
2537 pack->stripsize_f/2.;
2538
2539 for (int i=0;i<pack->nstrip_f;++i)
2540 for(int r=0;r<2;++r)
2541 {
2542
2543
2544 float sy = start_y+i*pack->stripsize_f;
2545 float sx = pack->gas_layer_radius[r];
2546 float z = nor(PosiZ)+ shift_z +
2547 pack->half_dim[1];
2548
2549 float Ce[3] = {sx,sy,z};
2550 rotate(Zp,Rp,Sp,Aa,Ba,Ga,Ce);
2551 float X = Ce[0]*cos(phi) - Ce[1]*sin(phi);
2552 float Y = Ce[0]*sin(phi) + Ce[1]*cos(phi);
2553 *(str_f2+i*4+r*2) = X;
2554 *(str_f2+i*4+r*2+1) = Y;
2555 }
2556
2557 pack->f_proj[1] = str_f1;
2558 pack->f_proj[0] = str_f2;
2559
2560 } else
2561 {
2562 float* str_f1 = new float [4*pack->nstrip_f];
2563
2564
2565
2566
2567
2568 float start_y = -pack->half_dim[0] + Dead +
2569 shift_s + orto_rad +
2570 pack->stripsize_f/2.;
2571
2572 for (int i=0;i<pack->nstrip_f;++i)
2573 for(int r=0;r<2;++r)
2574 {
2575
2576
2577 float sy = start_y+i*pack->stripsize_f;
2578 float sx = pack->gas_layer_radius[r];
2579 float z = nor(PosiZ)+ shift_z +
2580 pack->half_dim[1];
2581
2582 float Ce[3] = {sx,sy,z};
2583 rotate(Zp,Rp,Sp,Aa,Ba,Ga,Ce);
2584 float X = Ce[0]*cos(phi) - Ce[1]*sin(phi);
2585 float Y = Ce[0]*sin(phi) + Ce[1]*cos(phi);
2586 *(str_f1+i*4+r*2) = X;
2587 *(str_f1+i*4+r*2+1) = Y;
2588 }
2589
2590 if(shift_s!=0&&(orto_rad!=0 && isOrtoRad))
2591 {
2592 printf("RPCGeometry: inconsistence!\n");
2593 printf(" shift_s=%f, orto_rad=%f\n",
2594 shift_s,orto_rad);
2595 exit(1);
2596 }
2597
2598 if(shift_s<0||(orto_rad<0 && isOrtoRad))
2599 {
2600 pack->f_proj[1] = 0;
2601 pack->f_proj[0] = str_f1;
2602 } else
2603 {
2604 pack->f_proj[1] = str_f1;
2605 pack->f_proj[0] = 0;
2606 }
2607
2608
2609 }
2610
2611
2612
2613
2614
2615
2616
2617
2618
2619
2620
2621
2622
2623
2624
2625 ++(cham->nRPC_packs);
2626 if (!cham->RPC_packs)
2627 {
2628 cham->RPC_packs = pack;
2629 } else
2630 {
2631 last_pack->next = pack;
2632 }
2633 last_pack = pack;
2634
2635
2636
2637
2638
2639
2640
2641
2642
2643
2644
2645
2646
2647
2648
2649
2650
2651
2652
2653
2654
2655
2656
2657
2658
2659
2660
2661
2662
2663
2664
2665
2666
2667
2668
2669
2670 }
2671
2672 }
2673 }
2674
2675 }
2676 }
2677
2678
2679 }
2680 }
2681
2682
2683
2684
2685
2686
2687
2688
2689
2690
2691
2692
2693
2694
2695
2696
2697
2698
2699
2700
2701
2702
2703
2704
2705
2706
2707
2708
2709
2710
2711
2712
2713
2714
2715
2716
2717
2718
2719
2720
2721
2722
2723
2724
2725
2726
2727
2728
2729
2730
2731
2732
2733
2734
2735
2736
2737
2738
2739
2740
2741
2742
2743
2744
2745
2746
2747
2748
2749
2750
2751
2752
2753
2754
2755
2756
2757
2758
2759
2760
2761
2762
2763
2764
2765
2766
2767
2768
2769
2770
2771
2772
2773
2774
2775
2776
2777
2778
2779
2780
2781
2782
2783
2784
2785
2786
2787
2788
2789
2790
2791
2792
2793
2794
2795
2796
2797
2798
2799
2800
2801
2802
2803
2804
2805
2806
2807
2808
2809
2810
2811
2812
2813
2814
2815
2816
2817
2818
2819
2820
2821
2822
2823
2824
2825
2826
2827
2828
2829
2830
2831
2832
2833
2834
2835
2836
2837
2838
2839
2840
2841
2842
2843
2844
2845
2846
2847
2848
2849
2850
2851
2852
2853
2854
2855
2856
2857
2858
2859
2860
2861
2862
2863
2864
2865
2866
2867
2868
2869
2870
2871
2872
2873
2874
2875
2876
2877
2878
2879
2880
2881
2882
2883
2884
2885
2886
2887
2888
2889
2890
2891
2892
2893
2894
2895
2896
2897
2898
2899
2900
2901
2902
2903
2904
2905
2906
2907
2908
2909
2910
2911
2912
2913
2914
2915
2916
2917
2918
2919
2920
2921
2922
2923
2924
2925
2926
2927
2928
2929
2930
2931
2932
2933
2934
2935
2936
2937
2938
2939
2940
2941
2942
2943
2944
2945
2946
2947
2948
2949
2950
2951
2952
2953
2954
2955
2956
2957
2958
2959
2960
2961
2962
2963
2964
2965
2966
2967
2968
2969
2970
2971
2972
2973
2974
2975
2976
2977
2978
2979
2980
2981
2982
2983
2984
2985
2986
2987
2988
2989
2990
2991
2992
2993
2994
2995
2996
2997
2998
2999
3000
3001
3002
3003
3004
3005
3006
3007
3008
3009
3010
3011
3012
3013
3014
3015
3016
3017
3018
3019
3020
3021
3022
3023
3024
3025
3026
3027
3028
3029
3030
3031
3032
3033
3034
3035
3036
3037
3038
3039
3040
3041
3042
3043
3044
3045
3046
3047
3048
3049
3050
3051
3052
3053
3054
3055
3056
3057
3058
3059
3060
3061
3062
3063
3064
3065
3066
3067
3068
3069
3070
3071
3072
3073
3074
3075
3076
3077
3078
3079
3080
3081
3082
3083
3084
3085
3086
3087
3088
3089
3090
3091
3092
3093
3094
3095
3096
3097
3098
3099
3100
3101 }
3102 }
3103
3104
3105
3106 void
3107 RPCGeometry::rotate(float ZposCha, float RposCha, float OrtoRadp,
3108 float Alfa, float Beta, float Gamma,
3109 float xyz[3]) const
3110 {
3111
3112 Alfa = -Alfa*s_pi/180;
3113 Beta = -Beta*s_pi/180;
3114 Gamma = -Gamma*s_pi/180;
3115
3116 float t = xyz[0] - RposCha;
3117 float s = xyz[1] - OrtoRadp;
3118 float z = xyz[2] - ZposCha;
3119
3120
3121 s = s*cos(Alfa) + z*sin(Alfa);
3122 z = -s*sin(Alfa) + z*cos(Alfa);
3123
3124
3125 t = t*cos(Beta) + s*sin(Beta);
3126 s = -t*sin(Beta) + s*cos(Beta);
3127
3128
3129 z = z*cos(Gamma) + t*sin(Gamma);
3130 t = -z*sin(Gamma) + t*cos(Gamma);
3131
3132 xyz[0] = t + RposCha;
3133 xyz[1] = s + OrtoRadp;
3134 xyz[2] = z + ZposCha;
3135
3136 }
3137
3138 void
3139 RPCGeometry::readgeo() const
3140 {
3141
3142
3143
3144
3145
3146
3147
3148
3149
3150
3151
3152
3153
3154
3155
3156
3157
3158
3159
3160
3161
3162 }
3163
3164
3165
|
This page was automatically generated by the
LXR engine.
|
|