# Control on a sum

Hello, in this macro I calculate the emittance
emit.cpp (4.9 KB)

``````ts->SetBranchAddress("SecondaryParticleVert.x", &xVv);
for(jentry=0; jentry<nentries;jentry++) {
ts->GetEntry(jentry);
x+=k*xv-k*xVv;
x2+=(k*xv-k*xVv)*(k*xv-k*xVv);
xp+=xpv;
xp2+=xpv*xpv;
xxp+=(k*xv-k*xVv)*xpv;
}
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/=1e6;
``````

I must you see, I get the values from the TTree, I save them in an array the I calculate the sums

``````
x+=k*xv-k*xVv;
x2+=(k*xv-k*xVv)*(k*xv-k*xVv);
xp+=xpv;
xp2+=xpv*xpv;
xxp+=(k*xv-k*xVv)*xpv;
``````

Currently I use all the values in the ROOT file, but now I need a check, i.e. I must sum the values only if the new value is lower than N times the previous one…
for example, lets set N=3
I start from x=0…then if using the first entry the new x value is 0.5, I use the entry, if the new value is 3.2 I don’t use that entry…
Moreover, I must count how many entries I used because currently I calulcate

``````   x2*=1./nentries;
xp*=1./nentries;
xp2*=1./nentries;
xxp*=1./nentries;
``````

but if I have 1000 entries and I use only 990, I must divide for 990 and not for 1000…
Thank you

``````int n = 1;
double previous_value = value;
double sum = previous_value;
for (int i = 1; i <nentries; i++) {
if (value[i] < previous_value/N) {
previous_value = value[i] ;
sum = sum + previous_value ;
n++;
}
}
``````

That’s `n` in the previous code.

Thank you @couet .
I maybe didn’t understand something, because I wrote this code but it doesn’t use any event

emit.cpp (5.5 KB)

Very likely be the test:

``````            if (( k*xv[i]-k*xVv[i])<previous_x/N  && ((k*xv[i]-k*xVv[i])*(k*xv[i]-k*xVv[i]))<previous_x2/N && xpv[i]<previous_xp/N && (xpv[i]*xpv[i])<previous_xp2/N && ((k*xv[i]-k*xVv[i])*xpv[i])<previous_xxp/N) {

``````

filters all events.

Maybe in my first message I didn’t explain what I need.

I start from the first value of the entries…and I assign the x, x2,xp, xp2 and xxp values.
Then I check the second entry if the new x,x2, xp,xp2 and xxp2 (ie. the sum of the first and secon entries) are lower than N times the old x,x2, xp, xp2 I use the second entry and I sum i.e.

``````x+=...
x2+=...
etc
``````

otherwise I don’t use that entry. Same for the third entry… if the the sum of previous value and third entry is lower that N times the previous value I use the third entry

so maybe
` if (value[i] < previous_value/N) {`
isn’t the right check…or didn’t I understand the code?

I write an example

**First enrtry: **
x=0.5
xp=0.03

Second enrty
x=0.3
xp=0.04

I set the control N=3

the new x is 0.5+0.3=0.8
the new xp is 0.03+0.04=0.07
…then because `0.8<3*0.5 and 0.07< 3 * 0.03`
I can use this entry

Second example
**First enrtry: **
x=0.5
xp=0.03

Second enrty
x=0.3
xp=0.4

the new x is `0.5+0.3=0.8` but the new xp is` 0.03+0.4=0.43>3*0.03` then I don’t use the second entrry

That’s what I understood from your first message but if your algorithm is a bit different simply program it . You have described it. You should be able to transcribe it into C++.

Thank you @couet I think that I did it!

Hello @couet sorry to write again regarding this problem.
I modified the macro to have this algorithm:

a. I’ve a vector xpv
b. I get the first entry of xpv and I set xp=xpv
c. Regarding of netxt entries, If xp+xpv is in the range [-cxp;cxp] I use the entry otherwise I dont’ use it. For example, I set c=1
First entry is 2 then I set xp=2
The second entry is -1, because 2-1=1 and 1 is in the range [-2;2] I use this value and my new xp is 1
The third entry is 1 because 1+1=2 and 2 isn’t in the range [-1;1], I don’t use this value and my xp is still 1.

I wrote this macro
emitcont.cpp (5.1 KB)

i.e. you see, I wrote the line
` if(xp+xpv <= abs(c*xp) ){`

but I have a problem…you see here

for entry 3881 I’ve xp=-0.00411267 , the entry 3882 is -1.4384 and the macro get it, but it should reject this entry, because -0.00411267-1.4384=-1.44251267 and -1.44251267 is out of the range [-0.00411267;0.00411267]

Can you help to fix it please?

``````root 
Processing emiter.cpp...
In file included from input_line_9:1:
n++
^
;
``````

Not really good either:

``````root emitcut.cpp
------------------------------------------------------------------
| Welcome to ROOT 6.25/01                        https://root.cern |
| (c) 1995-2021, The ROOT Team; conception: R. Brun, F. Rademakers |
| Built for macosx64 on Mar 09 2021, 09:09:23                      |
| With Apple clang version 12.0.0 (clang-1200.0.32.29)             |
| Try '.help', '.demo', '.license', '.credits', '.quit'/'.q'       |
------------------------------------------------------------------

root 
Processing emitcut.cpp...
[/Users/couet/git/roottrunk-bin/lib/libCling.so] cling_runtime_internal_throwIfInvalidPointer (no debug info)
[<unknown binary>] (no debug info)
[<unknown binary>] (no debug info)
[/Users/couet/git/roottrunk-bin/lib/libCling.so] cling::IncrementalExecutor::executeWrapper(llvm::StringRef, cling::Value*) const (no debug info)
[/Users/couet/git/roottrunk-bin/lib/libCling.so] cling::Interpreter::RunFunction(clang::FunctionDecl const*, cling::Value*) (no debug info)
[/Users/couet/git/roottrunk-bin/lib/libCling.so] cling::Interpreter::EvaluateInternal(std::__1::basic_string<char, std::__1::char_traits<char>, std::__1::allocator<char> > const&, cling::CompilationOptions, cling::Value*, cling::Transaction**, unsigned long) (no debug info)
[/Users/couet/git/roottrunk-bin/lib/libCling.so] cling::MetaSema::actOnxCommand(llvm::StringRef, llvm::StringRef, cling::Value*) (no debug info)
[/Users/couet/git/roottrunk-bin/lib/libCling.so] cling::MetaParser::isXCommand(cling::MetaSema::ActionResult&, cling::Value*) (no debug info)
[/Users/couet/git/roottrunk-bin/lib/libCling.so] cling::MetaParser::isCommand(cling::MetaSema::ActionResult&, cling::Value*) (no debug info)
[/Users/couet/git/roottrunk-bin/lib/libCling.so] cling::MetaProcessor::process(llvm::StringRef, cling::Interpreter::CompilationResult&, cling::Value*, bool) (no debug info)
[/Users/couet/git/roottrunk-bin/lib/libCling.so] HandleInterpreterException(cling::MetaProcessor*, char const*, cling::Interpreter::CompilationResult&, cling::Value*) (no debug info)
[/Users/couet/git/roottrunk-bin/lib/libCling.so] TCling::ProcessLine(char const*, TInterpreter::EErrorCode*) (no debug info)
[/Users/couet/git/roottrunk-bin/lib/libCling.so] TCling::ProcessLineSynch(char const*, TInterpreter::EErrorCode*) (no debug info)
[/Users/couet/git/roottrunk-bin/lib/libCore.so] TApplication::ExecuteFile(char const*, int*, bool) (no debug info)
[/Users/couet/git/roottrunk-bin/lib/libRint.so] TRint::ProcessLineNr(char const*, char const*, int*) (no debug info)
[/Users/couet/git/roottrunk-bin/lib/libRint.so] TRint::Run(bool) (no debug info)
[/Users/couet/git/roottrunk-bin/bin/root.exe] main (no debug info)
[/usr/lib/system/libdyld.dylib] start (no debug info)
Error in <HandleInterpreterException>: Trying to dereference null pointer or trying to call routine taking non-null arguments.
Execution of your code was aborted.
In file included from input_line_9:1:
/Users/couet/Downloads/emitcut.cpp:23:4: warning: null passed to a callee that requires a non-null argument [-Wnonnull]
fin->ls ();
^~~
root 

``````

@couet I see you miss the root file… I’m uploading it…but it’s 800MB…so it will need time! I will give you the link as soon as it will finish

`emitcut.cpp` does not have any `cout << "jentry = "` so I don’t think that’s the code that produces the output you are showing?

Yes sorry… @Axel you are right

this is the right one

emitcont.cpp (5.1 KB)

``````if (xp+xpv <= abs(c*xp) ){
``````

Here, `c` is 1, `xp` is `-1.44251` and `xpv` is `-1.4384`. So the `if` evaluates

``````if (-1.44251 + (-1.4384) <= abs(1.*(-1.44251)))
``````

i.e.

``````if (-2.88091 <= 1.44251)
``````

and the answer is: yes, totally.

Look - it’s better if you use `printf` to trace things you don’t understand in your code, or you use ACLiC (`.L emitcont.cpp+`) and use `gdb` to debug: it allows you to see all the values while things happen! It’s not really possible for us to help all our many many users get over their non-ROOT-related programming issues…

Don’t worry @axel and @couet
I think I solved my self
writing

``````ts->GetEntry(0);
x+=k*xv-k*xVv;
x2+=(k*xv-k*xVv)*(k*xv-k*xVv);
xp+=xpv;
xp2+=xpv*xpv;
xxp+=(k*xv-k*xVv)*xpv;
for(jentry=1; jentry<nentries;jentry++) {
ts->GetEntry(jentry);
if(abs(xp+xpv )<= abs(c*xp)  && abs(x+(k*xv-k*xVv)) <= abs(c*x) ){
x+=k*xv-k*xVv;
x2+=(k*xv-k*xVv)*(k*xv-k*xVv);
xp+=xpv;
xp2+=xpv*xpv;
xxp+=(k*xv-k*xVv)*xpv;
``````
1 Like

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