Problem calculating the mean values

Hello, I wrote this code to calculate the emittance

Double_t emittance=0.;
			Double_t x=0;
			Double_t x2=0;	
			Double_t xp=0;
			Double_t xp2=0;
			Double_t xxp=0;
			Double_t sigmax2=0;
			Double_t sigmaxp2=0;
			Double_t sigmaxsigmaxp=0;
			int jentry;
			int nentries = ts->GetEntries();
			ts->SetBranchAddress("PositionDirection.x", &x);
			ts->SetBranchAddress("SecondaryParticleAng.x", &xp);
		    for(jentry=0; jentry<nentries;jentry++) {
				ts->GetEntry(jentry);
				x+=x;
				x2+=x*x;
				xp+=xp;
				xp2+=xp*xp;
				xxp+=x*xp;
			}
			x*=1./nentries;
			x2*=1./nentries;
			xp*=1./nentries;
			xp2*=1./nentries;
			xxp*=1./nentries;
			sigmax2=x2-x*x;
			sigmaxp2=xp2-xp*xp;
			sigmaxsigmaxp=xxp-x*xp;
			std::cout << "N = " << nentries << "\t x = " << x << "\t x2 = " << 2 << " \t xp = " << xp <<  " \t xp2= " << xp2 << " \t xxp= " <<  xxp << std::endl;
			emittance=sigmax2*sigmaxp2-sigmaxsigmaxp;
			emittance=TMath::Sqrt(emittance);
			emittance*=1000000;

You see, that I’ve to calculate some mean values (as , etc) but these values are very different than the values calculated by Root. For example, you see here my =-3,765328*10^-7

Instead Root calculated = -0.000227

and surely the Root value is fine because using my value, I get a very big emittance (impossible)…but I can’t understand what is wrong in my code.

this is the macro emit.cpp (3.0 KB)
and this is the Root file https://we.tl/t-TMHdHaqL5D

Thank you


Please read tips for efficient and successful posting and posting code

ROOT Version: Not Provided
Platform: Not Provided
Compiler: Not Provided


Hi,
The difference could be caused by the fact that in the ROOT mean value maybe the underflow and overflow;s you are having in y are not used. You should call before filling the histogram TH1::StatOverflows(true) and you should get the same value

Lorenzo

Hello @moneta thank you. Unfortunately I tried the code, and the problem was due to the code and not to the underflow and overflow!

Anyway I solvede myself
I had to use the array!

Double_t emittance=0.;
		Double_t x=0;
		Double_t x2=0;	
		Double_t xp=0;
		Double_t xp2=0;
		Double_t xxp=0;
		Double_t sigmax2=0;
		Double_t sigmaxp2=0;
		Double_t sigmaxsigmaxp=0;
		int jentry;
		int nentries = ts->GetEntries();
		double xv[10];
		double xpv[10];
		ts->SetBranchAddress("PositionDirection.x", &xv[0]);
		ts->SetBranchAddress("SecondaryParticleAng.x", &xpv[0]);
	    for(jentry=0; jentry<nentries;jentry++) {
			ts->GetEntry(jentry);
			x+=xv[0];
			x2+=xv[0]*xv[0];
			xp+=xpv[0];
			xp2+=xpv[0]*xpv[0];
			xxp+=xv[0]*xpv[0];
		}
		x*=1./nentries;
		x2*=1./nentries;
		xp*=1./nentries;
		xp2*=1./nentries;
		xxp*=1./nentries;
		sigmax2=x2-x*x;
		sigmaxp2=xp2-xp*xp;
		sigmaxsigmaxp=xxp-x*xp;
		emittance=sigmax2*sigmaxp2-sigmaxsigmaxp;
		emittance=TMath::Sqrt(emittance);
		emittance*=1000000;
		std::cout << "N = " << nentries << "\t x = " << x << "\t x2 = " << x2 << " \t xp = " << xp <<  " \t xp2= " << xp2 << " \t xxp= " <<  xxp << std::endl;

Now it works fine! Here the results and they are equal to the Root results!