[u/mrichter/AliRoot.git] / TPC / CalibMacros / CalibExB.C

void Init(){
  //
  // Initialize
  //
  AliMagFMaps* field = new AliMagFMaps("Maps","Maps", 2, 1., 10., AliMagFMaps::k5kG);
  AliTPCExB::RegisterField(0,field);
  AliMagF* fieldC0 = new AliMagWrapCheb("Maps","Maps", 2, 1, 10., AliMagWrapCheb::k5kG);
  AliTPCExB::RegisterField(1,fieldC0);
  AliMagF* fieldC1 = new AliMagWrapCheb("Maps","Maps", 2, 1, 10., AliMagWrapCheb::k5kG,kTRUE,"$(ALICE_ROOT)/data/maps/mfchebKGI_sym.root");
  AliTPCExB::RegisterField(2,fieldC1);

  gSystem->Load("libSTAT.so");
  AliTPCExBFirst *exbfirst1  = new  AliTPCExBFirst(fieldC1,0.88*2.6400e+04,50,50,50);
  AliTPCExB::SetInstance(exbfirst1);

}

void FitField(){
  //
  //
  //
  AliTPCExB * exb =  AliTPCExB::Instance();
  exb->TestExB("field.root");
  TFile f("field.root");
  TTree * tree = (TTree*)f.Get("positions");
  //
  TStatToolkit toolkit;
  Double_t chi2;
  TVectorD fitParam;
  TMatrixD covMatrix;
  Int_t npoints;
  //
  // SetAliases
  //
  tree->SetAlias("sa","sin(phi+0.0)");
  tree->SetAlias("ca","cos(phi+0.0)");
  tree->SetAlias("sa2","sin(phi*2+0.0)");
  tree->SetAlias("ca2","cos(phi*2+0.0)");
  tree->SetAlias("zn","(x2/250.)");
  tree->SetAlias("rn","(r/250.)");
    
  TString fstringSym="";
  //  
  fstringSym+="zn++";
  fstringSym+="rn++";
  fstringSym+="zn*rn++";
  fstringSym+="zn*zn++";
  fstringSym+="zn*zn*rn++";
  fstringSym+="zn*rn*rn++";
  //
  fstringSym+="sa++";
  fstringSym+="ca++";  
  fstringSym+="ca2++";
  fstringSym+="sa2++";
  fstringSym+="ca*zn++";
  fstringSym+="sa*zn++";
  fstringSym+="ca2*zn++";
  fstringSym+="sa2*zn++";
  fstringSym+="ca*zn*zn++";
  fstringSym+="sa*zn*zn++";
  fstringSym+="ca*zn*rn++";
  fstringSym+="sa*zn*rn++";

  // BrBz
  TString *strBrBz = toolkit.FitPlane(tree,"br/bz",fstringSym->Data(), "abs(r)<250&&abs(r)>90", chi2,npoints,fitParam,covMatrix);
  strBrBz->Tokenize("+")->Print();
  printf("Chi2/npoints = %f\n",TMath::Sqrt(chi2/npoints));
  tree->SetAlias("fitBrBz",strBrBz->Data());
  exb->fMatBrBz = new TVectorD(fitParam);
  // BrfiBz
  TString *strBrfiBz = toolkit.FitPlane(tree,"brfi/bz",fstringSym->Data(), "abs(r)<250&&abs(r)>90", chi2,npoints,fitParam,covMatrix);
  strBrfiBz->Tokenize("+")->Print();
  printf("Chi2/npoints = %f\n",TMath::Sqrt(chi2/npoints));
  tree->SetAlias("fitBrfiBz",strBrfiBz->Data());
  exb->fMatBrfiBz = new TVectorD(fitParam);

  
  //
  // BR integral parameterization
  //
  TString *strBRI0 = toolkit.FitPlane(tree,"bri",fstringSym->Data(), "abs(r)<250&&abs(r)>90&&x2>0&&x2<250", chi2,npoints,fitParam,covMatrix);
  strBRI0->Tokenize("+")->Print();
  printf("Chi2/npoints = %f\n",TMath::Sqrt(chi2/npoints));
  tree->SetAlias("fitBRI0",strBRI0->Data());
  exb->fMatBrBzI0 = new TVectorD(fitParam);

  TString *strBRI1 = toolkit.FitPlane(tree,"bri",fstringSym->Data(), "abs(r)<250&&abs(r)>90&&x2<-10&&x2>-250", chi2,npoints,fitParam,covMatrix);
  strBRI1->Tokenize("+")->Print();
  printf("Chi2/npoints = %f\n",TMath::Sqrt(chi2/npoints));
  tree->SetAlias("fitBRI1",strBRI1->Data());
  exb->fMatBrBzI1 = new TVectorD(fitParam);
  //

  TString *strBRFII0 = toolkit.FitPlane(tree,"brfii",fstringSym->Data(), "abs(r)<250&&abs(r)>90&&x2>0&&x2<250", chi2,npoints,fitParam,covMatrix);
  strBRFII0->Tokenize("+")->Print();
  printf("Chi2/npoints = %f\n",TMath::Sqrt(chi2/npoints));
  tree->SetAlias("fitBRFII0",strBRFII0->Data());
  exb->fMatBrfiBzI0 = new TVectorD(fitParam);
 
  TString *strBRFII1 = toolkit.FitPlane(tree,"brfii",fstringSym->Data(), "abs(r)<250&&abs(r)>90&&x2<-10&&x2>-240", chi2,npoints,fitParam,covMatrix);
  strBRFII1->Tokenize("+")->Print();
  printf("Chi2/npoints = %f\n",TMath::Sqrt(chi2/npoints));
  tree->SetAlias("fitBRFII1",strBRFII1->Data());
  exb->fMatBrfiBzI1 = new TVectorD(fitParam);

  TFile * fout = new TFile("bfit.root","recreate");
  exb->Write("bfit");
  fout->Close();

}


void CalibExB(){
  // 
  // Macro to test ExB correction versus laser B filed scan data
  //
  // Before running it be sure the file laserScan.root exist
  //
  // Create this file from scan - See AliTPCCalibLaser.C
  //
  TFile fscan("laserScan.root");
  TTree * treeT = (TTree*)fscan.Get("Mean");
  TStatToolkit toolkit;
  Double_t chi2;
  TVectorD fitParam;
  TMatrixD covMatrix;
  Int_t npoints;
  
  TCut cutF0("abs(gphi1-(pphi0+pphi1*bz))<0.05");
  TCut cutF1("abs(gphiP1-(pphiP0+pphiP1*bz))<0.0005");
  TCut cutN("entries>2");
  
  TCut cutA = cutF0+cutF1+cutN;
  
  treeT->SetAlias("side","(-1+(LTr.fP[1]>0)*2)");       // side
  treeT->SetAlias("dr","(abs(LTr.fP[1]/250.))");
  treeT->SetAlias("sa","sin(atan2(lx1+0.0,lx0+0.0))");
  treeT->SetAlias("ca","cos(atan2(lx1+0.0,lx0+0.0))");
  treeT->SetAlias("ta","tan(asin(LTr.fP[2]+0.0))");

  TString fstring="";
  //
  fstring+="((dr)^3-1)*bz++";             //1
  fstring+="((dr)^3-1)*bz*sa++";        //3
  fstring+="((dr)^3-1)*bz*ca++";        //5
  //
  fstring+="((dr)^1-1)*bz++";           //1
  fstring+="((dr)^1-1)*bz*sa++";        //3
  fstring+="((dr)^1-1)*bz*ca++";        //5
  //
  fstring+="side*((dr)^3-1)*bz++";             //1
  fstring+="side*((dr)^3-1)*bz*sa++";        //3
  fstring+="side*((dr)^3-1)*bz*ca++";        //5
  //
  fstring+="side*((dr)^1-1)*bz++";           //7
  fstring+="side*((dr)^1-1)*bz*sa++";        //9
  fstring+="side*((dr)^1-1)*bz*ca++";        //11
  //
  fstring+="((dr)^3-1)*bz^2*ta++";        //2
  fstring+="((dr)^3-1)*bz^2*sa*ta++";     //4
  fstring+="((dr)^3-1)*bz^2*ca*ta++";     //6

  fstring+="((dr)^1-1)*bz^2*ta++";        //2
  fstring+="((dr)^1-1)*bz^2*sa*ta++";     //4
  fstring+="((dr)^1-1)*bz^2*ca*ta++";     //6
  //  
  fstring+="side*((dr)^1-1)*bz^2*ta++";        //8 
  fstring+="side*((dr)^1-1)*bz^2*sa*ta++";     //10
  fstring+="side*((dr)^1-1)*bz^2*ca*ta++";     //12


  TString *strq0 = toolkit.FitPlane(treeT,"gphi1-pphi0",fstring->Data(), cutA, chi2,npoints,fitParam,covMatrix);
  strq0->Tokenize("+")->Print();
  treeT->SetAlias("fit",strq0->Data());
  printf("Chi2/npoints = %f\n",TMath::Sqrt(chi2/npoints));


  treeT->SetAlias("bir1",  "-AliTPCExB::GetBrI(254,atan2(lx1,lx0+0),LTr.fP[1],2)");
  treeT->SetAlias("birfi1","-AliTPCExB::GetBrfiI(254,atan2(lx1,lx0+0),LTr.fP[1],2)");

  treeT->SetAlias("bir0",  "AliTPCExB::GetBrI(254,atan2(lx1,lx0+0),LTr.fP[1],2)");
  treeT->SetAlias("birfi0","AliTPCExB::GetBrfiI(254,atan2(lx1,lx0+0),LTr.fP[1],2)");

  treeT->SetAlias("fbz00", "(bir0+birfi0*ta)");
  treeT->SetAlias("fbz02", "(birfi0+bir0*ta)");
  treeT->SetAlias("fbz10", "(bir1+birfi1*ta)");
  treeT->SetAlias("fbz12", "(birfi1+bir1*ta)");
  //
  treeT->SetAlias("fbz0", "((fSide==0)*fbz00+(fSide==1)*fbz10)");
  treeT->SetAlias("fbz2", "((fSide==0)*fbz02+(fSide==1)*fbz12)");

  //
  TString fstringeb="";
  //
  fstringeb+="bz*fbz0++";                //1
  fstringeb+="sign(bz)*bz^2*fbz2++";     //2
  fstringeb+="((dr)^3-1)*bz*sa++";        //9
  fstringeb+="side*((dr)^3-1)*bz*sa++";        //9
  fstringeb+="side*((dr)^3-1)*bz++";           //7
  //fstringeb+="((dr)^3-1)*bz++";           //7
  //fstringeb+="((dr)^3-1)*bz*ca++";        //11
  //fstringeb+="side*((dr)^3-1)*bz*sa++";        //9
  //fstringeb+="side*((dr)^3-1)*bz*ca++";        //11
  
  //  
  TString *strExB = toolkit.FitPlane(treeT,"gphi1-pphi0",fstringeb->Data(), "abs(gphi1-pphi0-fit)<0.06&&abs(bz)>0.1"+cutA, chi2,npoints,fitParam,covMatrix);
  strExB->Tokenize("+")->Print();
  printf("Chi2/npoints = %f\n",TMath::Sqrt(chi2/npoints));
  treeT->SetAlias("fitEB",strExB->Data());


 TString fstringeb="";
 treeT->SetAlias("rr","sqrt(lx0^2+(lx1+0)^2)/250.");
 //
 fstringeb+="((dr)^3-1)*bz*rr++";        //9
 fstringeb+="((dr)^1-1)*bz*rr++";        //9
 fstringeb+="side*((dr)^3-1)*bz*rr++";        //9
 fstringeb+="side*((dr)^1-1)*bz*rr++";        //9
 //
 fstringeb+="((dr)^3-1)*bz*sa++";        //9
 fstringeb+="((dr)^1-1)*bz*sa++";        //9
 fstringeb+="side*((dr)^3-1)*bz*sa++";        //9
 fstringeb+="side*((dr)^1-1)*bz*sa++";        //9


  {
    Float_t otcor0    =0.7;
    Float_t omegatau0 =0.35;
    Float_t rms0      =10;
    for (Float_t otcor=0.6; otcor<1.1;otcor+=0.1)
      for (Float_t ometau=0.4; ometau<0.45;ometau+=0.01){
	char hname[100];
	sprintf(hname,"hist_ometau_%f(50,-0.10,0.10)",ometau);
	char expr[100];
	sprintf(expr,"gphi1-pphi0-bz*fbz0*%f-sign(bz)*(%f*%f*bz)^2*fbz2>>%s",ometau,ometau,otcor,hname);
	treeT->Draw(expr,"abs(gphi1-pphi0-fit)<0.06&&abs(bz)>0.1"+cutA);    
	printf("Ometau=%f\tCor=%f\tRMS=%f\n",ometau,otcor,treeT->GetHistogram()->GetRMS());
	if (rms0>treeT->GetHistogram()->GetRMS()){
	  otcor0=otcor;
	  omegatau0=ometau;
	  rms0=treeT->GetHistogram()->GetRMS();
	}
      }
  }
 
}


void MakePic(){
  //
  //
  //

}
Commit	Line	Data
ec5bd527	1	void Init(){
	2	//
	3	// Initialize
	4	//
ec5bd527	5	AliMagFMaps* field = new AliMagFMaps("Maps","Maps", 2, 1., 10., AliMagFMaps::k5kG);
ec5bd527	6	AliTPCExB::RegisterField(0,field);
58388598	7	AliMagF* fieldC0 = new AliMagWrapCheb("Maps","Maps", 2, 1, 10., AliMagWrapCheb::k5kG);
	8	AliTPCExB::RegisterField(1,fieldC0);
	9	AliMagF* fieldC1 = new AliMagWrapCheb("Maps","Maps", 2, 1, 10., AliMagWrapCheb::k5kG,kTRUE,"$(ALICE_ROOT)/data/maps/mfchebKGI_sym.root");
	10	AliTPCExB::RegisterField(2,fieldC1);
	11
ec5bd527	12	gSystem->Load("libSTAT.so");
58388598	13	AliTPCExBFirst exbfirst1 = new AliTPCExBFirst(fieldC1,0.882.6400e+04,50,50,50);
ec5bd527	14	AliTPCExB::SetInstance(exbfirst1);
	15
	16	}
	17
	18	void FitField(){
	19	//
	20	//
	21	//
	22	AliTPCExB * exb = AliTPCExB::Instance();
	23	exb->TestExB("field.root");
	24	TFile f("field.root");
	25	TTree * tree = (TTree*)f.Get("positions");
	26	//
	27	TStatToolkit toolkit;
	28	Double_t chi2;
	29	TVectorD fitParam;
	30	TMatrixD covMatrix;
	31	Int_t npoints;
	32	//
	33	// SetAliases
	34	//
	35	tree->SetAlias("sa","sin(phi+0.0)");
	36	tree->SetAlias("ca","cos(phi+0.0)");
	37	tree->SetAlias("sa2","sin(phi*2+0.0)");
	38	tree->SetAlias("ca2","cos(phi*2+0.0)");
	39	tree->SetAlias("zn","(x2/250.)");
	40	tree->SetAlias("rn","(r/250.)");
	41
	42	TString fstringSym="";
	43	//
	44	fstringSym+="zn++";
	45	fstringSym+="rn++";
	46	fstringSym+="zn*rn++";
	47	fstringSym+="zn*zn++";
	48	fstringSym+="znznrn++";
	49	fstringSym+="znrnrn++";
	50	//
	51	fstringSym+="sa++";
	52	fstringSym+="ca++";
	53	fstringSym+="ca2++";
	54	fstringSym+="sa2++";
	55	fstringSym+="ca*zn++";
	56	fstringSym+="sa*zn++";
	57	fstringSym+="ca2*zn++";
	58	fstringSym+="sa2*zn++";
	59	fstringSym+="caznzn++";
	60	fstringSym+="saznzn++";
	61	fstringSym+="caznrn++";
	62	fstringSym+="saznrn++";
	63
	64	// BrBz
	65	TString *strBrBz = toolkit.FitPlane(tree,"br/bz",fstringSym->Data(), "abs(r)<250&&abs(r)>90", chi2,npoints,fitParam,covMatrix);
	66	strBrBz->Tokenize("+")->Print();
	67	printf("Chi2/npoints = %f\n",TMath::Sqrt(chi2/npoints));
	68	tree->SetAlias("fitBrBz",strBrBz->Data());
	69	exb->fMatBrBz = new TVectorD(fitParam);
	70	// BrfiBz
	71	TString *strBrfiBz = toolkit.FitPlane(tree,"brfi/bz",fstringSym->Data(), "abs(r)<250&&abs(r)>90", chi2,npoints,fitParam,covMatrix);
	72	strBrfiBz->Tokenize("+")->Print();
	73	printf("Chi2/npoints = %f\n",TMath::Sqrt(chi2/npoints));
	74	tree->SetAlias("fitBrfiBz",strBrfiBz->Data());
	75	exb->fMatBrfiBz = new TVectorD(fitParam);
	76
	77
78	//
79	// BR integral parameterization
80	//
81	TString *strBRI0 = toolkit.FitPlane(tree,"bri",fstringSym->Data(), "abs(r)<250&&abs(r)>90&&x2>0&&x2<250", chi2,npoints,fitParam,covMatrix);
82	strBRI0->Tokenize("+")->Print();
83	printf("Chi2/npoints = %f\n",TMath::Sqrt(chi2/npoints));
84	tree->SetAlias("fitBRI0",strBRI0->Data());
85	exb->fMatBrBzI0 = new TVectorD(fitParam);
86
87	TString *strBRI1 = toolkit.FitPlane(tree,"bri",fstringSym->Data(), "abs(r)<250&&abs(r)>90&&x2<-10&&x2>-250", chi2,npoints,fitParam,covMatrix);
88	strBRI1->Tokenize("+")->Print();
89	printf("Chi2/npoints = %f\n",TMath::Sqrt(chi2/npoints));
90	tree->SetAlias("fitBRI1",strBRI1->Data());
91	exb->fMatBrBzI1 = new TVectorD(fitParam);
92	//
93
94	TString *strBRFII0 = toolkit.FitPlane(tree,"brfii",fstringSym->Data(), "abs(r)<250&&abs(r)>90&&x2>0&&x2<250", chi2,npoints,fitParam,covMatrix);
95	strBRFII0->Tokenize("+")->Print();
96	printf("Chi2/npoints = %f\n",TMath::Sqrt(chi2/npoints));
97	tree->SetAlias("fitBRFII0",strBRFII0->Data());
98	exb->fMatBrfiBzI0 = new TVectorD(fitParam);
99
100	TString *strBRFII1 = toolkit.FitPlane(tree,"brfii",fstringSym->Data(), "abs(r)<250&&abs(r)>90&&x2<-10&&x2>-240", chi2,npoints,fitParam,covMatrix);
101	strBRFII1->Tokenize("+")->Print();
102	printf("Chi2/npoints = %f\n",TMath::Sqrt(chi2/npoints));
103	tree->SetAlias("fitBRFII1",strBRFII1->Data());
104	exb->fMatBrfiBzI1 = new TVectorD(fitParam);
105
106	TFile * fout = new TFile("bfit.root","recreate");
107	exb->Write("bfit");
108	fout->Close();
109
110	}
111
112
113
114
115
116
27f5d48d	117	void CalibExB(){
	118	//
	119	// Macro to test ExB correction versus laser B filed scan data
	120	//
	121	// Before running it be sure the file laserScan.root exist
27f5d48d	122	//
	123	// Create this file from scan - See AliTPCCalibLaser.C
	124	//
	125	TFile fscan("laserScan.root");
	126	TTree * treeT = (TTree*)fscan.Get("Mean");
27f5d48d	127	TStatToolkit toolkit;
	128	Double_t chi2;
	129	TVectorD fitParam;
	130	TMatrixD covMatrix;
	131	Int_t npoints;
	132
	133	TCut cutF0("abs(gphi1-(pphi0+pphi1*bz))<0.05");
	134	TCut cutF1("abs(gphiP1-(pphiP0+pphiP1*bz))<0.0005");
	135	TCut cutN("entries>2");
	136
	137	TCut cutA = cutF0+cutF1+cutN;
	138
ec5bd527	139	treeT->SetAlias("side","(-1+(LTr.fP[1]>0)*2)"); // side
27f5d48d	140	treeT->SetAlias("dr","(abs(LTr.fP[1]/250.))");
	141	treeT->SetAlias("sa","sin(atan2(lx1+0.0,lx0+0.0))");
	142	treeT->SetAlias("ca","cos(atan2(lx1+0.0,lx0+0.0))");
	143	treeT->SetAlias("ta","tan(asin(LTr.fP[2]+0.0))");
	144
	145	TString fstring="";
	146	//
5fd9157a	147	fstring+="((dr)^3-1)*bz++"; //1
27f5d48d	148	fstring+="((dr)^3-1)bzsa++"; //3
27f5d48d	149	fstring+="((dr)^3-1)bzca++"; //5
5fd9157a	150	//
	151	fstring+="((dr)^1-1)*bz++"; //1
	152	fstring+="((dr)^1-1)bzsa++"; //3
	153	fstring+="((dr)^1-1)bzca++"; //5
	154	//
ec5bd527	155	fstring+="side((dr)^3-1)bz++"; //1
	156	fstring+="side((dr)^3-1)bz*sa++"; //3
	157	fstring+="side((dr)^3-1)bz*ca++"; //5
	158	//
27f5d48d	159	fstring+="side((dr)^1-1)bz++"; //7
27f5d48d	160	fstring+="side((dr)^1-1)bz*sa++"; //9
27f5d48d	161	fstring+="side((dr)^1-1)bz*ca++"; //11
5fd9157a	162	//
	163	fstring+="((dr)^3-1)bz^2ta++"; //2
	164	fstring+="((dr)^3-1)bz^2sa*ta++"; //4
	165	fstring+="((dr)^3-1)bz^2ca*ta++"; //6
	166
	167	fstring+="((dr)^1-1)bz^2ta++"; //2
	168	fstring+="((dr)^1-1)bz^2sa*ta++"; //4
	169	fstring+="((dr)^1-1)bz^2ca*ta++"; //6
	170	//
	171	fstring+="side((dr)^1-1)bz^2*ta++"; //8
	172	fstring+="side((dr)^1-1)bz^2sata++"; //10
	173	fstring+="side((dr)^1-1)bz^2cata++"; //12
27f5d48d	174
	175
	176
	177
	178	TString *strq0 = toolkit.FitPlane(treeT,"gphi1-pphi0",fstring->Data(), cutA, chi2,npoints,fitParam,covMatrix);
	179	strq0->Tokenize("+")->Print();
	180	treeT->SetAlias("fit",strq0->Data());
ec5bd527	181	printf("Chi2/npoints = %f\n",TMath::Sqrt(chi2/npoints));
	182
	183
27f5d48d	184
58388598	185	treeT->SetAlias("bir1", "-AliTPCExB::GetBrI(254,atan2(lx1,lx0+0),LTr.fP[1],2)");
58388598	186	treeT->SetAlias("birfi1","-AliTPCExB::GetBrfiI(254,atan2(lx1,lx0+0),LTr.fP[1],2)");
ec5bd527	187
58388598	188	treeT->SetAlias("bir0", "AliTPCExB::GetBrI(254,atan2(lx1,lx0+0),LTr.fP[1],2)");
58388598	189	treeT->SetAlias("birfi0","AliTPCExB::GetBrfiI(254,atan2(lx1,lx0+0),LTr.fP[1],2)");
ec5bd527	190
	191	treeT->SetAlias("fbz00", "(bir0+birfi0*ta)");
	192	treeT->SetAlias("fbz02", "(birfi0+bir0*ta)");
	193	treeT->SetAlias("fbz10", "(bir1+birfi1*ta)");
	194	treeT->SetAlias("fbz12", "(birfi1+bir1*ta)");
	195	//
	196	treeT->SetAlias("fbz0", "((fSide==0)fbz00+(fSide==1)fbz10)");
	197	treeT->SetAlias("fbz2", "((fSide==0)fbz02+(fSide==1)fbz12)");
	198
	199	//
27f5d48d	200	TString fstringeb="";
27f5d48d	201	//
ec5bd527	202	fstringeb+="bz*fbz0++"; //1
	203	fstringeb+="sign(bz)bz^2fbz2++"; //2
	204	fstringeb+="((dr)^3-1)bzsa++"; //9
	205	fstringeb+="side((dr)^3-1)bz*sa++"; //9
	206	fstringeb+="side((dr)^3-1)bz++"; //7
	207	//fstringeb+="((dr)^3-1)*bz++"; //7
	208	//fstringeb+="((dr)^3-1)bzca++"; //11
	209	//fstringeb+="side((dr)^3-1)bz*sa++"; //9
	210	//fstringeb+="side((dr)^3-1)bz*ca++"; //11
27f5d48d	211
ec5bd527	212	//
58388598	213	TString *strExB = toolkit.FitPlane(treeT,"gphi1-pphi0",fstringeb->Data(), "abs(gphi1-pphi0-fit)<0.06&&abs(bz)>0.1"+cutA, chi2,npoints,fitParam,covMatrix);
27f5d48d	214	strExB->Tokenize("+")->Print();
ec5bd527	215	printf("Chi2/npoints = %f\n",TMath::Sqrt(chi2/npoints));
27f5d48d	216	treeT->SetAlias("fitEB",strExB->Data());
ec5bd527	217
	218
	219
	220	TString fstringeb="";
	221	treeT->SetAlias("rr","sqrt(lx0^2+(lx1+0)^2)/250.");
	222	//
	223	fstringeb+="((dr)^3-1)bzrr++"; //9
	224	fstringeb+="((dr)^1-1)bzrr++"; //9
	225	fstringeb+="side((dr)^3-1)bz*rr++"; //9
	226	fstringeb+="side((dr)^1-1)bz*rr++"; //9
	227	//
	228	fstringeb+="((dr)^3-1)bzsa++"; //9
	229	fstringeb+="((dr)^1-1)bzsa++"; //9
	230	fstringeb+="side((dr)^3-1)bz*sa++"; //9
	231	fstringeb+="side((dr)^1-1)bz*sa++"; //9
	232
	233
	234
	235
	236
	237
	238	{
	239	Float_t otcor0 =0.7;
	240	Float_t omegatau0 =0.35;
	241	Float_t rms0 =10;
	242	for (Float_t otcor=0.6; otcor<1.1;otcor+=0.1)
	243	for (Float_t ometau=0.4; ometau<0.45;ometau+=0.01){
	244	char hname[100];
	245	sprintf(hname,"hist_ometau_%f(50,-0.10,0.10)",ometau);
	246	char expr[100];
	247	sprintf(expr,"gphi1-pphi0-bzfbz0%f-sign(bz)(%f%fbz)^2fbz2>>%s",ometau,ometau,otcor,hname);
	248	treeT->Draw(expr,"abs(gphi1-pphi0-fit)<0.06&&abs(bz)>0.1"+cutA);
	249	printf("Ometau=%f\tCor=%f\tRMS=%f\n",ometau,otcor,treeT->GetHistogram()->GetRMS());
	250	if (rms0>treeT->GetHistogram()->GetRMS()){
	251	otcor0=otcor;
	252	omegatau0=ometau;
	253	rms0=treeT->GetHistogram()->GetRMS();
	254	}
	255	}
	256	}
	257
27f5d48d	258	}
ec5bd527	259
	260
	261
	262
58388598	263	void MakePic(){
	264	//
	265	//
	266	//
	267
	268	}
	269