Introduction of the Copyright and cvs Log
[u/mrichter/AliRoot.git] / MUON / AliMUONClusterFinderv0.cxx
CommitLineData
4c039060 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15
16/*
17$Log$
18*/
19
a897a37a 20#include "AliMUONClusterFinderv0.h"
21#include "AliMUONSegResV1.h"
22//#include "TTree.h"
23//#include "AliRun.h"
24//#include <TCanvas.h>
25//#include <TH1.h>
26//#include <TPad.h>
27//#include <TGraph.h>
28
29//----------------------------------------------------------
30ClassImp(AliMUONClusterFinderv0)
31
32 AliMUONClusterFinderv0::AliMUONClusterFinderv0
33(AliMUONsegmentation *segmentation, AliMUONresponse *response,
34 TClonesArray *digits, Int_t chamber) : AliMUONClusterFinder(segmentation,response,digits,chamber)
35{;}
36
37 AliMUONClusterFinderv0::AliMUONClusterFinderv0():AliMUONClusterFinder()
38{;}
39
40/*
41void AliMUONClusterFinder::AddRawCluster(const AliMUONRawCluster c)
42{
43 //
44 // Add a raw cluster copy to the list
45 //
46 AliMUON *MUON=(AliMUON*)gAlice->GetModule("MUON");
47 MUON->AddRawCluster(fChamber,c);
48 fNRawClusters++;
49}
50*/
51
52
53
54void AliMUONClusterFinderv0::Decluster(AliMUONRawCluster *cluster)
55{
56// AliMUONdigit *dig;
57// Int_t q;
58 static int done=0;
59 if (!done) {
60 printf("Calling decluster\n");
61 done=1;
62 }
63
64
65
66 Int_t mul = cluster->fMultiplicity;
67// printf("Decluster - multiplicity %d \n",mul);
68
69 if (mul == 1) {
70// printf("\n Nothing special for 1-clusters \n");
71//
72// Nothing special for 1-clusters
73//
74 AddRawCluster(*cluster);
75 } else if (mul ==2) {
76//
77// 2-cluster, compute offset
78//
79 SetOffset(cluster);
80 FillCluster(cluster);
81 AddRawCluster(*cluster);
82 } else if (mul ==3) {
83//
84// 3-cluster, check topology
85// printf("\n 3-cluster, check topology \n");
86//
87 if (Centered(cluster)) {
88//
89// ok, cluster is centered
90// printf("\n ok, cluster is centered \n");
91 } else {
92//
93// cluster is not centered, split into 2+1
94// printf("\n cluster is not centered, split into 2+1 \n");
95 }
96
97 } else {
98 if (mul >(50-5)) printf("Decluster - multiplicity %d approaching 50\n",mul);
99//
100// 4-and more-pad clusters
101//
102 SplitByLocalMaxima(cluster);
103 } // multiplicity
104}
105
106Int_t AliMUONClusterFinderv0::PeakOffsetAndCoordinates(Int_t DigitIndex, Float_t *X, Float_t *Y)
107//
108// Computes for which allowed offsets the digit has the highest neighbouring charge
109// Returns the value of the offset, and sets the pyisical coordinates of that pad
110// Loop on physical neighbours is specific to AliMUONsegmentationV1
111{
112Int_t nPara, offset, returnOffset=0 ;
113AliMUONdigit* dig= (AliMUONdigit*)fDigits->UncheckedAt(DigitIndex);
114AliMUONsegmentationV1* seg = (AliMUONsegmentationV1*) fSegmentation;
115seg->GetNParallelAndOffset(dig->fPadX,dig->fPadY,&nPara,&offset);
116if (nPara>1)
117 {
118 Float_t qMax=0;
119 for (Int_t i=0;i<nPara; i++)
120 {
121 // Compute the charge on the 9 neighbouring pads
122 // We assume that there are no pads connected in parallel in the neighbourhood
123 //
124 Float_t q=0;
125 for (Int_t dx=-1;dx<2;dx++)
126 for (Int_t dy=-1;dy<2;dy++)
127 {
128 if (dx==dy && dy==0)
129 continue;
130 Int_t padY=dig->fPadY+dy;
131 Int_t padX=seg->Ix((Int_t) (dig->fPadX+dx+i*offset) , padY);
132 if (fHitMap->TestHit(padX, padY)==empty)
133 continue;
134 AliMUONdigit* digt = (AliMUONdigit*) fHitMap->GetHit(padX,padY);
135 q += digt->fSignal;
136 }
137 if (q>qMax)
138 {
139 returnOffset=i*offset;
140 qMax=q;
141 }
142 }
143 }
144fSegmentation->GetPadCxy(dig->fPadX+returnOffset,dig->fPadY,*X,*Y);
145return returnOffset;
146}
147
148
149void AliMUONClusterFinderv0::SetOffset(AliMUONRawCluster *cluster)
150// compute the offsets assuming that there is only one peak !
151{
152//DumpCluster(cluster);
153Float_t X,Y;
154cluster->fOffsetMap[0]=PeakOffsetAndCoordinates(cluster->fIndexMap[0],&X,&Y);
155for (Int_t i=1;i<cluster->fMultiplicity;i++) {
156 AliMUONdigit* dig= (AliMUONdigit*)fDigits->UncheckedAt(cluster->fIndexMap[i]);
157 fSegmentation->Distance2AndOffset(dig->fPadX,dig->fPadY,X,Y,&(cluster->fOffsetMap[i]));
158 }
159}
160
161void AliMUONClusterFinderv0::DumpCluster(AliMUONRawCluster *cluster)
162{
163printf ("other cluster\n");
164for (Int_t i=0; i<cluster->fMultiplicity; i++)
165 {
166 AliMUONdigit* dig= (AliMUONdigit*)fDigits->UncheckedAt(cluster->fIndexMap[i]);
167 Int_t nPara, offset;
168 fSegmentation->GetNParallelAndOffset(dig->fPadX,dig->fPadY,&nPara,&offset);
169
170 printf("X %d Y %d Q %d NPara %d \n",dig->fPadX, dig->fPadY,dig->fSignal, nPara);
171 }
172}
173
174Bool_t AliMUONClusterFinderv0::Centered(AliMUONRawCluster *cluster)
175{
176 AliMUONdigit* dig;
177 dig= (AliMUONdigit*)fDigits->UncheckedAt(cluster->fIndexMap[0]);
178 Int_t ix=dig->fPadX;
179 Int_t iy=dig->fPadY;
180 Int_t nn;
181 Int_t X[kMaxNeighbours], Y[kMaxNeighbours], XN[kMaxNeighbours], YN[kMaxNeighbours];
182 fSegmentation->Neighbours(ix,iy,&nn,X,Y);
183
184 Int_t nd=0;
185 for (Int_t i=0; i<nn; i++) {
186 if (fHitMap->TestHit(X[i],Y[i]) == used) {
187 XN[nd]=X[i];
188 YN[nd]=Y[i];
189 nd++;
190 }
191 }
192 if (nd==2) {
193//
194// cluster is centered !
195 SetOffset(cluster);
196 FillCluster(cluster);
197 AddRawCluster(*cluster);
198 return kTRUE;
199 } else if (nd ==1) {
200//
201// Highest signal on an edge, split cluster into 2+1
202//
203// who is the neighbour ?
204 Int_t nind=fHitMap->GetHitIndex(XN[0], YN[0]);
205 Int_t i1= (nind==cluster->fIndexMap[1]) ? 1:2;
206 Int_t i2= (nind==cluster->fIndexMap[1]) ? 2:1;
207//
208// 2-cluster
209 AliMUONRawCluster cnew;
210 cnew.fMultiplicity=2;
211 cnew.fIndexMap[0]=cluster->fIndexMap[0];
212 cnew.fIndexMap[1]=cluster->fIndexMap[i1];
213 SetOffset(&cnew);
214 FillCluster(&cnew);
215 AddRawCluster(cnew);
216//
217// 1-cluster
218 cluster->fMultiplicity=1;
219 cluster->fIndexMap[0]=cluster->fIndexMap[i2];
220 cluster->fIndexMap[1]=0;
221 cluster->fIndexMap[2]=0;
222 FillCluster(cluster);
223 AddRawCluster(*cluster);
224 return kFALSE;
225 } else {
226 printf("\n Completely screwed up %d !! \n",nd);
227
228 }
229
230 return kFALSE;
231}
232
233
234void AliMUONClusterFinderv0::SplitByLocalMaxima(AliMUONRawCluster *c)
235{
236 AliMUONdigit* dig[50], *digt;
237 Int_t ix[50], iy[50], q[50];
238 Float_t x[50], y[50];
239 Int_t i; // loops over digits
240 Int_t j; // loops over local maxima
241
242 Int_t mul=c->fMultiplicity;
243//
244// dump digit information into arrays
245//
246 for (i=0; i<mul; i++)
247 {
248 dig[i]= (AliMUONdigit*)fDigits->UncheckedAt(c->fIndexMap[i]);
249 ix[i]= dig[i]->fPadX;
250 iy[i]= dig[i]->fPadY;
251 q[i] = dig[i]->fSignal;
252 fSegmentation->GetPadCxy(ix[i], iy[i], x[i], y[i]);
253 }
254//
255// Find local maxima
256//
257 Bool_t IsLocal[50];
258 Int_t NLocal=0;
259 Int_t AssocPeak[50];
260 Int_t IndLocal[50];
261 Int_t nn;
262 Int_t X[kMaxNeighbours], Y[kMaxNeighbours];
263 for (i=0; i<mul; i++) {
264 fSegmentation->Neighbours(ix[i], iy[i], &nn, X, Y);
265 IsLocal[i]=kTRUE;
266 for (j=0; j<nn; j++) {
267 if (fHitMap->TestHit(X[j], Y[j])==empty) continue;
268 digt=(AliMUONdigit*) fHitMap->GetHit(X[j], Y[j]);
269 if (digt->fSignal > q[i]) {
270 IsLocal[i]=kFALSE;
271 break;
272//
273// handle special case of neighbouring pads with equal signal
274 } else if (digt->fSignal == q[i]) {
275 if (NLocal >0) {
276 for (Int_t k=0; k<NLocal; k++) {
277 if (X[j]==ix[IndLocal[k]] && Y[j]==iy[IndLocal[k]]){
278 IsLocal[i]=kFALSE;
279 }
280 }
281 }
282 }
283 } // loop over next neighbours
284 if (IsLocal[i]) {
285 IndLocal[NLocal]=i;
286 // New for LYON : we guess which is the actual position of the pad hit
287 // But this would run like that for normal chamber !
288 c->fOffsetMap[i]=PeakOffsetAndCoordinates(c->fIndexMap[i], &(x[i]), &(y[i]));
289 NLocal++;
290 }
291 } // loop over all digits
292// printf("Found %d local Maxima",NLocal);
293//
294// Check if enough local clusters have been found,
295// if not add global maxima to the list
296//
297// But what the hell is that ? (Manu)
298//
299// Int_t nPerMax=mul/NLocal;
300// if (nPerMax > 5) {
301// Int_t nGlob=mul/5-NLocal+1;
302// if (nGlob > 0) {
303// Int_t nnew=0;
304// for (i=0; i<mul; i++) {
305// if (!IsLocal[i]) {
306// IndLocal[NLocal]=i;
307// IsLocal[i]=kTRUE;
308// NLocal++;
309// nnew++;
310// }
311// if (nnew==nGlob) break;
312// }
313// }
314// }
315//
316//
317// Associate hits to peaks
318//
319 for (i=0; i<mul; i++) {
320 //
321 // loop on digits
322 //
323 Float_t dmin=1.E10;
324 Float_t qmax=0;
325 Int_t offset;
326 if (IsLocal[i]) continue;
327 for (j=0; j<NLocal; j++) {
328 //
329 // Loop on peaks
330 //
331 Int_t il=IndLocal[j];
332// Float_t d=TMath::Sqrt((x[i]-x[il])*(x[i]-x[il])
333// +(y[i]-y[il])*(y[i]-y[il]));
334 // Can run like that for non-Lyon chambers
335 Float_t d = fSegmentation->Distance2AndOffset(ix[i],iy[i],x[il],y[il], &offset);
336 Float_t ql=q[il];
337//
338// Select nearest peak
339//
340 if (d<dmin) {
341 dmin=d;
342 qmax=ql;
343 AssocPeak[i]=j;
344 c->fOffsetMap[i]=offset;
345 } else if (d==dmin) {
346//
347// If more than one take highest peak
348//
349 if (ql>qmax) {
350 dmin=d;
351 qmax=ql;
352 AssocPeak[i]=j;
353 c->fOffsetMap[i]=offset;
354 }
355 } // end if
356 } // End loop on peaks
357 } // end loop on digits
358//
359// One cluster for each maximum
360//
361 for (j=0; j<NLocal; j++) {
362 AliMUONRawCluster cnew;
363 cnew.fIndexMap[0]=c->fIndexMap[IndLocal[j]];
364 cnew.fOffsetMap[0]=c->fOffsetMap[IndLocal[j]];
365 cnew.fMultiplicity=1;
366 for (i=0; i<mul; i++) {
367 if (IsLocal[i]) continue;
368 if (AssocPeak[i]==j) {
369 cnew.fIndexMap[cnew.fMultiplicity]=c->fIndexMap[i];
370 cnew.fOffsetMap[cnew.fMultiplicity]=c->fOffsetMap[i];
371 cnew.fMultiplicity++;
372 }
373 }
374 FillCluster(&cnew);
375 AddRawCluster(cnew);
376 }
377}
378
379/*
380void AliMUONClusterFinderv0::FillCluster(AliMUONRawCluster* c)
381{
382//
383// Completes cluster information starting from list of digits
384//
385 AliMUONdigit* dig;
386 Float_t x, y;
387 Int_t ix, iy;
388
389 c->fPeakSignal=0;
390 c->fX=0;
391 c->fY=0;
392 c->fQ=0;
393 for (Int_t i=0; i<c->fMultiplicity; i++)
394 {
395 dig= (AliMUONdigit*)fDigits->UncheckedAt(c->fIndexMap[i]);
396 ix=dig->fPadX + c.fOffsetMap[i]; // should be 0 for non-LYON
397 iy=dig->fPadY;
398 Int_t q=dig->fSignal;
399//
400//
401// peak signal and track list
402 if (q>c->fPeakSignal) {
403 c->fPeakSignal=0;
404 c->fTracks[0]=dig->fTracks[0];
405 c->fTracks[1]=dig->fTracks[1];
406 c->fTracks[2]=dig->fTracks[2];
407 }
408//
409// centre of gravity
410 fSegmentation->GetPadCxy(ix, iy, x, y);
411 c->fX += q*x;
412 c->fY += q*y;
413 c->fQ += q;
414 }
415 c->fX/=c->fQ;
416// Not valid for inclined tracks in X !!! (Manu)
417// c->fX=fSegmentation->GetAnod(c->fX);
418 c->fY/=c->fQ;
419//
420// apply correction to the coordinate along the anode wire
421//
422 if (fCogCorr) {
423 x=c->fX;
424 y=c->fY;
425 fSegmentation->GetPadIxy(x, y, ix, iy);
426 fSegmentation->GetPadCxy(ix, iy, x, y);
427 Float_t YonPad=(c->fY-y)/fSegmentation->Dpy();
428 c->fY=y-fCogCorr->Eval(YonPad, 0, 0);
429 }
430
431}
432*/
433/*
434void AliMUONClusterFinder::FindCluster(Int_t i, Int_t j, AliMUONRawCluster &c){
435//
436// Find clusters
437//
438//
439// Add i,j as element of the cluster
440//
441
442 Int_t idx = fHitMap->GetHitIndex(i,j);
443 AliMUONdigit* dig = (AliMUONdigit*) fHitMap->GetHit(i,j);
444 Int_t q=dig->fSignal;
445 if (q > TMath::Abs(c.fPeakSignal)) {
446 c.fPeakSignal=q;
447 c.fTracks[0]=dig->fTracks[0];
448 c.fTracks[1]=dig->fTracks[1];
449 c.fTracks[2]=dig->fTracks[2];
450 }
451//
452// Make sure that list of digits is ordered
453//
454 Int_t mu=c.fMultiplicity;
455 c.fIndexMap[mu]=idx;
456
457 if (mu > 0) {
458 for (Int_t ind=mu-1; ind>=0; ind--) {
459 Int_t ist=(c.fIndexMap)[ind];
460 Int_t ql=((AliMUONdigit*)fDigits
461 ->UncheckedAt(ist))->fSignal;
462 if (q>ql) {
463 c.fIndexMap[ind]=idx;
464 c.fIndexMap[ind+1]=ist;
465 } else {
466 break;
467 }
468 }
469 }
470
471 c.fMultiplicity++;
472
473 if (c.fMultiplicity >= 50 ) {
474 printf("FindCluster - multiplicity >50 %d \n",c.fMultiplicity);
475 c.fMultiplicity=50-1;
476 }
477
478// Prepare center of gravity calculation
479 Float_t x, y;
480 fSegmentation->GetPadCxy(i, j, x, y);
481 c.fX += q*x;
482 c.fY += q*y;
483 c.fQ += q;
484// Flag hit as taken
485 fHitMap->FlagHit(i,j);
486//
487// Now look recursively for all neighbours
488//
489 Int_t nn;
490 Int_t Xlist[kMaxNeighbours], Ylist[kMaxNeighbours];
491 fSegmentation->Neighbours(i,j,&nn,Xlist,Ylist);
492 for (Int_t in=0; in<nn; in++) {
493 Int_t ix=Xlist[in];
494 Int_t iy=Ylist[in];
495 if (fHitMap->TestHit(ix,iy)==unused) FindCluster(ix, iy, c);
496 }
497}
498*/
499
500//_____________________________________________________________________________
501
502void AliMUONClusterFinderv0::FindRawClusters()
503{
504 //
505 // simple MUON cluster finder from digits -- finds neighbours and
506 // fill the tree with raw clusters
507 //
508 if (!fNdigits) return;
509
510 fHitMap = new AliMUONHitMapA1(fSegmentation, fDigits);
511
512 AliMUONdigit *dig;
513
514 int ndig;
515 int nskip=0;
516
517 fHitMap->FillHits();
518 for (ndig=0; ndig<fNdigits; ndig++) {
519 dig = (AliMUONdigit*)fDigits->UncheckedAt(ndig);
520 Int_t i=dig->fPadX;
521 Int_t j=dig->fPadY;
522 if (fHitMap->TestHit(i,j)==used ||fHitMap->TestHit(i,j)==empty) {
523 nskip++;
524 continue;
525 }
526 AliMUONRawCluster c;
527 c.fMultiplicity=0;
528// c.fPeakSignal=dig->fSignal;
529// c.fTracks[0]=dig->fTracks[0];
530// c.fTracks[1]=dig->fTracks[1];
531// c.fTracks[2]=dig->fTracks[2];
532 c.fPeakSignal=0;
533 FindCluster(i,j, c);
534 // center of gravity
535 c.fX /= c.fQ;
536 c.fX=fSegmentation->GetAnod(c.fX);
537 c.fY /= c.fQ;
538//
539// apply correction to the coordinate along the anode wire
540//
541
542
543 if (fCogCorr) {
544 Int_t ix,iy;
545 Float_t x=c.fX;
546 Float_t y=c.fY;
547 fSegmentation->GetPadIxy(x, y, ix, iy);
548 fSegmentation->GetPadCxy(ix, iy, x, y);
549 Float_t YonPad=(c.fY-y)/fSegmentation->Dpy();
550 c.fY=y-fCogCorr->Eval(YonPad,0,0);
551 }
552//
553// Analyse cluster and decluster if necessary
554//
555 Decluster(&c);
556//
557//
558//
559// reset Cluster object
560 for (int k=0;k<c.fMultiplicity;k++) {
561 c.fIndexMap[k]=0;
562 c.fOffsetMap[k]=0;
563 }
564 c.fMultiplicity=0;
565 } // end loop ndig
566 delete fHitMap;
567}
568
569/*
570
571void AliMUONClusterFinder::
572CalibrateCOG()
573{
574 Float_t x[5];
575 Float_t y[5];
576 Int_t n, i;
577 TF1 func;
578
579 if (fSegmentation) {
580 fSegmentation->GiveTestPoints(n, x, y);
581 for (i=0; i<n; i++) {
582 Float_t xtest=x[i];
583 Float_t ytest=y[i];
584 SinoidalFit(xtest, ytest, func);
585 }
586 fCogCorr = new TF1(func);
587 }
588}
589*/
590/*
591
592void AliMUONClusterFinder::
593SinoidalFit(Float_t x, Float_t y, TF1 &func)
594{
595//
596 static Int_t count=0;
597 char canvasname[3];
598 count++;
599 sprintf(canvasname,"c%d",count);
600
601// MANU : without const, error on HP
602 const Int_t ns=101;
603 Float_t xg[ns], yg[ns], xrg[ns], yrg[ns];
604 Float_t xsig[ns], ysig[ns];
605
606 AliMUONsegmentation *segmentation=fSegmentation;
607
608 Int_t ix,iy;
609 segmentation->GetPadIxy(x,y,ix,iy);
610 segmentation->GetPadCxy(ix,iy,x,y);
611 Int_t isec=segmentation->Sector(ix,iy);
612// Pad Limits
613 Float_t xmin = x-segmentation->GetRealDpx(isec)/2;
614 Float_t ymin = y-segmentation->Dpy()/2;
615//
616// Integration Limits
617 Float_t dxI=fResponse->Nsigma()*fResponse->ChwX();
618 Float_t dyI=fResponse->Nsigma()*fResponse->ChwY();
619
620//
621// Scanning
622//
623 Int_t i;
624 Float_t qp;
625//
626// y-position
627 Float_t yscan=ymin;
628 Float_t dy=segmentation->Dpy()/(ns-1);
629
630 for (i=0; i<ns; i++) {
631//
632// Pad Loop
633//
634 Float_t sum=0;
635 Float_t qcheck=0;
636 segmentation->SigGenInit(x, yscan, 0);
637
638 for (segmentation->FirstPad(x, yscan, dxI, dyI);
639 segmentation->MorePads();
640 segmentation->NextPad())
641 {
642 qp=fResponse->IntXY(segmentation);
643 qp=TMath::Abs(qp);
644//
645//
646 if (qp > 1.e-4) {
647 qcheck+=qp;
648 Int_t ixs=segmentation->Ix();
649 Int_t iys=segmentation->Iy();
650 Float_t xs,ys;
651 segmentation->GetPadCxy(ixs,iys,xs,ys);
652 sum+=qp*ys;
653 }
654 } // Pad loop
655 Float_t ycog=sum/qcheck;
656 yg[i]=(yscan-y)/segmentation->Dpy();
657 yrg[i]=(ycog-y)/segmentation->Dpy();
658 ysig[i]=ycog-yscan;
659 yscan+=dy;
660 } // scan loop
661//
662// x-position
663 Float_t xscan=xmin;
664 Float_t dx=segmentation->GetRealDpx(isec)/(ns-1);
665
666 for (i=0; i<ns; i++) {
667//
668// Pad Loop
669//
670 Float_t sum=0;
671 Float_t qcheck=0;
672 segmentation->SigGenInit(xscan, y, 0);
673
674 for (segmentation->FirstPad(xscan, y, dxI, dyI);
675 segmentation->MorePads();
676 segmentation->NextPad())
677 {
678 qp=fResponse->IntXY(segmentation);
679 qp=TMath::Abs(qp);
680//
681//
682 if (qp > 1.e-2) {
683 qcheck+=qp;
684 Int_t ixs=segmentation->Ix();
685 Int_t iys=segmentation->Iy();
686 Float_t xs,ys;
687 segmentation->GetPadCxy(ixs,iys,xs,ys);
688 sum+=qp*xs;
689 }
690 } // Pad loop
691 Float_t xcog=sum/qcheck;
692 xcog=segmentation->GetAnod(xcog);
693
694 xg[i]=(xscan-x)/segmentation->GetRealDpx(isec);
695 xrg[i]=(xcog-x)/segmentation->GetRealDpx(isec);
696 xsig[i]=xcog-xscan;
697 xscan+=dx;
698 }
699
700 TCanvas *c1=new TCanvas(canvasname,canvasname,400,10,600,700);
701 TPad* pad11 = new TPad("pad11"," ",0.01,0.51,0.49,0.99);
702 TPad* pad12 = new TPad("pad12"," ",0.51,0.51,0.99,0.99);
703 TPad* pad13 = new TPad("pad13"," ",0.01,0.01,0.49,0.49);
704 TPad* pad14 = new TPad("pad14"," ",0.51,0.01,0.99,0.49);
705 pad11->SetFillColor(11);
706 pad12->SetFillColor(11);
707 pad13->SetFillColor(11);
708 pad14->SetFillColor(11);
709 pad11->Draw();
710 pad12->Draw();
711 pad13->Draw();
712 pad14->Draw();
713 TGraph *graphx = new TGraph(ns,xg ,xsig);
714 TGraph *graphxr= new TGraph(ns,xrg,xsig);
715 TGraph *graphy = new TGraph(ns,yg ,ysig);
716 TGraph *graphyr= new TGraph(ns,yrg,ysig);
717//
718// Creates a Root function based on function sinoid above
719// and perform the fit
720//
721 Double_t sinoid(Double_t *x, Double_t *par);
722 TF1 *sinoidf = new TF1("sinoidf",sinoid,0.5,0.5,5);
723 graphyr->Fit("sinoidf","V");
724 sinoidf->Copy(func);
725 func.Eval(0,0,0);
726//
727 pad11->cd();
728 graphx->SetFillColor(42);
729 graphx->SetMarkerColor(4);
730 graphx->SetMarkerStyle(21);
731 graphx->Draw("AC");
732 graphx->GetHistogram()->SetXTitle("x on pad");
733 graphx->GetHistogram()->SetYTitle("xcog-x");
734
735
736 pad12->cd();
737 graphxr->SetFillColor(42);
738 graphxr->SetMarkerColor(4);
739 graphxr->SetMarkerStyle(21);
740 graphxr->Draw("AP");
741 graphxr->GetHistogram()->SetXTitle("xcog on pad");
742 graphxr->GetHistogram()->SetYTitle("xcog-x");
743
744
745 pad13->cd();
746 graphy->SetFillColor(42);
747 graphy->SetMarkerColor(4);
748 graphy->SetMarkerStyle(21);
749 graphy->Draw("AF");
750 graphy->GetHistogram()->SetXTitle("y on pad");
751 graphy->GetHistogram()->SetYTitle("ycog-y");
752
753
754
755 pad14->cd();
756 graphyr->SetFillColor(42);
757 graphyr->SetMarkerColor(4);
758 graphyr->SetMarkerStyle(21);
759 graphyr->Draw("AF");
760 graphyr->GetHistogram()->SetXTitle("ycog on pad");
761 graphyr->GetHistogram()->SetYTitle("ycog-y");
762
763 c1->Update();
764
765}
766*/
767/*
768Double_t sinoid(Double_t *x, Double_t *par)
769{
770 Double_t arg = -2*TMath::Pi()*x[0];
771 Double_t fitval= par[0]*TMath::Sin(arg)+
772 par[1]*TMath::Sin(2*arg)+
773 par[2]*TMath::Sin(3*arg)+
774 par[3]*TMath::Sin(4*arg)+
775 par[4]*TMath::Sin(5*arg);
776 return fitval;
777 }
778*/
779
780
781
782
783
784
785