Calculate the sigma between results of the same mesurement

Hello, I measured some physics quantities values using 3 different methods then I have the three values, for example one is the ratio of two parameter and I got :

a. 0.692 +/-0.011
b. 0.6833 +/- 0.0053
c. 0.6712 +/- 0.0036

is there an automatic way to calculate the according between the three result (For example how many sigma)?

Thank you


Please read tips for efficient and successful posting and posting code

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


Try searching the forum for that. Here is one example result

Thank you @dastudillo I read @moneta 's code, but I didn’t understand so much…


std::vector<double> x(100);
gRandom->RndmArray(x.size(),x.data());
TMath::StdDev(x.begin(), x.end() );

for example, to calcucalte how many sigma is the according between the three results

 a. 0.692 +/-0.011
b. 0.6833 +/- 0.0053
c. 0.6712 +/- 0.0036

where do I have to write these number in @moneta 's code?
I tried to write

gRandom->RndmArray(x.size(3),x.data(0.692, 0.6833, 0.6712));

but I got this error

I also tried @eguiraud’s code in this way

using namespace ROOT::VecOps;
RVec<float> rv = {0.11, 0.0053, 0.0036};
auto srddev = StdDev(rv);

but
a. Does this methos use only the errors of the results and not the values (0.692, 0.6833, 0.6712)?
b. You see…ROOT closes automatically when I run the

auto srddev = StdDev(rv);
thank you

Yes, RVec::StdDev does a “dumb” calculation of the square root of the variance of the population, nothing fancy

That really should not happen. In fact I typed those exact same commands in a ROOT prompt on my laptop and things are fine:

~ root -l
root [0] using namespace ROOT::VecOps;
root [1] RVec<float> rv = {0.11, 0.0053, 0.0036};
root [2] auto srddev = StdDev(rv);
root [3] srddev
(double) 0.060945248
root [4]

Not sure what’s going on there :confused:

Thank you @eguiraud I tried the code using gentoo and it worked

I use Windows… @bellenot could it be a bug on Root for windows?

@eguiraud you said that this way is a dumb way… then, do you think that it’s better @moneta’s code

std::vector<double> x(100);
gRandom->RndmArray(x.size(),x.data());
TMath::StdDev(x.begin(), x.end() );

to get the according between the 3 measurements?

@moneta please, can you check the problem? I resume… I measured physics quantities using 3 instruments…for example for a ratio I got these values (one for each instrument)

a. 0.692 +/-0.011
b. 0.6833 +/- 0.0053
c. 0.6712 +/- 0.0036

I need a test to estimate the according between the 3 results.
I followed your code

std::vector<double> x(100);
gRandom->RndmArray(x.size(),x.data());
TMath::StdDev(x.begin(), x.end() );

then I wrote

std::vector<double> x(100);
gRandom->RndmArray(x.size(3),x.data(0.692, 0.6833, 0.6712));

but I got this error

Can you check please what I did wrong?
and, is there a better way for my goal how many sigma is the according?
Thank you

No no I meant that RVec::StdDev is dumb, i.e. it does not do anything smart like taking into account error values attached to the inputs :slight_smile:

As far as I understand TMath::StdDev does the same by default but it also has an overload that can take into account input weights, you might be able to do something with that.

This is not a valid syntax (just try on Linux)

And this code works for me on Windows:

C:\Users\bellenot\DAQ>root -l
root [0] using namespace ROOT::VecOps;
root [1] RVec<float> rv = {0.11, 0.0053, 0.0036};
root [2] auto srddev = StdDev(rv);
root [3] srddev
(double) 0.060945248
root [4]

@eguiraud yes, I understood that you meant the RVec::StdDev is dumb, that’s because I asked to @moneta if there is a better way for my goal!

@moneta: important… I said that I made 3 measurements using three instruments…but I got results using the measurements and making some calucations and fit …
then, for example

 a. 0.692 +/-0.011
b. 0.6833 +/- 0.0053
c. 0.6712 +/- 0.0036

are NOT the values read on the instruments, but are results of calculations using the values read by the instruments

1 Like

Thank you @bellenot

[quote=“bellenot, post:7, topic:43901, full:true”]

then how should I do?

I don’t know, try to follow what @eguiraud or @dastudillo proposed…

Thank you @bellenot I tried to follow the example showed by @dastudillo… but @eguiraud’s code doesn’t look like what I need… I was trying @moneta’s code but you see, I dont’ understand how to apply it…so I hope @moneta can help me, I know he is the most suitable to help regarding on the analysis way using ROOT…

or maybe, also @jblomer could know a way to calculcate the according between the results (because of he works in particle physics)

Hi,
I am not sure I have understood you completly.
You have 3 measurements, and what you would like to know, if they agree to each other ?
If this is the case, it is very simple, you can compute the mean, and compute their chisquared value, from which you can obtained a p-value. To do this, is easier by performing a simple fit with an horizontal line:

{
std::vector<double> x = {1,2,3};
std::vector<double> y = { 0.692,  0.6833, 0.6712 };
std::vector<double> ey = { 0.011,  0.0053,  0.0036};

auto g = new TGraphErrors(3,x.data(),y.data(), nullptr, ey.data()); 
auto result = g->Fit("pol0","S"); 
double chi2 = result->Chi2(); 
double ndf = result->Ndf(); 
double pvalue = ROOT::Math::chisquared_cdf_c(chi2, ndf); 
std::cout << "p value of the fit  : " << pvalue << std::endl;
 // number of standard sigma. Divide p-value by 2 since deviations can be in both sides
std::cout << "number of sigma is " << ROOT::Math::normal_quantile_c(0.5*pvalue,1) << std::endl;
}

Lorenzo

2 Likes

Thank you @moneta it was what I needed!

@bellenot I tried the code on my root (windows machine) but when I lunch the command

double ndf = result->Ndf();
Root crashes

I also tried the macro

sigmavalue.cpp (2.6 KB)

but I get this error


instead, both (the code on the terminal and the macro) work fine using gentoo by ssh connection.

Is it a problem on root for windows?
Thank you

It works for me with ROOT master (6.25/01):

C:\Users\bellenot\rootdev>root
   ------------------------------------------------------------------
  | Welcome to ROOT 6.25/01                        https://root.cern |
  | (c) 1995-2021, The ROOT Team; conception: R. Brun, F. Rademakers |
  | Built for win32 on Mar 09 2021, 17:28:03                         |
  | From heads/master@v6-25-01-74-g202a84bb2e                        |
  | With MSVC 19.24.28319.0                                          |
  | Try '.help', '.demo', '.license', '.credits', '.quit'/'.q'       |
   ------------------------------------------------------------------

root [0] std::vector<double> x = {1,2,3};
root [1] std::vector<double> y = { 0.692,  0.6833, 0.6712 };
root [2] std::vector<double> ey = { 0.011,  0.0053,  0.0036};
root [3] auto g = new TGraphErrors(3,x.data(),y.data(), nullptr, ey.data());
root [4] auto result = g->Fit("pol0","S");

****************************************
Minimizer is Linear / Migrad
Chi2                      =      5.78669
NDf                       =            2
p0                        =      0.67618   +/-   0.0028745
root [5] double chi2 = result->Chi2();
root [6] double ndf = result->Ndf();
root [7] double pvalue = ROOT::Math::chisquared_cdf_c(chi2, ndf);
root [8] std::cout << "p value of the fit  : " << pvalue << std::endl;
p value of the fit  : 0.0553906
root [9] std::cout << "number of sigma is " << ROOT::Math::normal_quantile_c(0.5*pvalue,1) << std::endl;
number of sigma is 1.9158
root [10]

And if you look at the (many) compilation errors with your code, you should be able to fix them. If you want me to try, please provide also the input file.

You are at least missing several headers:

#include "TF1.h"
#include "Math/ProbFuncMathCore.h"
#include "Math/QuantFuncMathCore.h"
#include <cstdlib>
#include <iostream>
#include <fstream>

Then it compiles fine:

C:\Users\bellenot\rootdev>root -l
root [0] .L sigmavalue.cpp+
Info in <TWinNTSystem::ACLiC>: creating shared library C:/Users/bellenot/rootdev/sigmavalue_cpp.dll
sigmavalue_cpp_ACLiC_dict.cxx
   Creating library C:/Users/bellenot/rootdev\sigmavalue_cpp.lib and object C:/Users/bellenot/rootdev\sigmavalue_cpp.exp
sigmavalue_cpp_ACLiC_dict.cxx
   Creating library C:/Users/bellenot/rootdev\sigmavalue_cpp.lib and object C:/Users/bellenot/rootdev\sigmavalue_cpp.exp
root [1] .q

Thank you @bellenot adding the headers the macro worked on windows too!

@moneta sorry…one more question please…
to quantify the agree between an experimental measurement and the theoretical value, is there some function on ROOT?
For example, the Planck’s constant value is h=6.62607004 *10^-34Js…
I got h=(6.424+/- 0.052)*10^-34 Js
I want to quantify the goodness of the experimental result… is there a function in ROOT?
Thank you

Hi,
In this case with a single measurement, you just compute the difference divide by the error.
This will give you number of sigma’s.

Lorenzo

Thank you @moneta

This topic was automatically closed 14 days after the last reply. New replies are no longer allowed.