Exit[]

MPI = 0.1349766; META = 0.547862; MP = 0.93827;
Mass = {MPI,MP,META,MPI,MP};

Dalitz and Van Hove boundraries

The border of the Dalitz Plots

uminmax[plab_,W_,mass_]:=Block[{ma,mb,mc,md,qi,qf,s,umin,umax},
ma = mass[[1]];mb = mass[[2]]; mc =W; md = mass[[5]];
s=ma^2+mb^2+2mb Sqrt[plab^2+ma^2];
qi = Sqrt[(s-(ma+mb)^2)(s-(ma-mb)^2)/(4.s)]; 
qf = Sqrt[(s-(mc+md)^2)(s-(mc-md)^2)/(4.s)]; 
umin = (ma^2-mb^2-mc^2+md^2)^2/4./s-(qi-qf)^2;
umax = (ma^2-mb^2-mc^2+md^2)^2/4./s-(qi+qf)^2;
{umin,umax}
]

snorm[pv_]:=pv[[4]]^2 - pv[[1]]^2- pv[[2]]^2- pv[[3]]^2; 
lam[a_,b_,c_]:=a^2+b^2+c^2-2(a b + b c + c a);

ToDalitz[event_,Nmax_]:=Block[{xtab,s12, s23, s31,p1,p2,p3},
xtab = ConstantArray[0,{Nmax,3}];
Do[
{p1,p2,p3}=event[[ie]];
s12 = snorm[p1+p2];s23 = snorm[p2+p3];s31 = snorm[p3+p1];
xtab[[ie,All]]={s12,s23,s31};
,{ie,Nmax}];
xtab
]

sijbounds[s_,MT_,Mi_,Mj_,Mk_]:=Block[{Ei,Ej,pi,pj},
Ej = (s-Mk^2+Mj^2)/(2Sqrt[s]);
Ei = (MT^2-s-Mi^2)/(2Sqrt[s]);
pj = Sqrt[Ej^2-Mj^2];
pi = Sqrt[Ei^2-Mi^2];
{(Ej+Ei)^2-(pj+pi)^2,(Ej+Ei)^2-(pj-pi)^2}];

M12 = Mass[[3]]^2; M22 = Mass[[4]]^2;M32 = Mass[[5]]^2;
s12min =(Sqrt[M12]+Sqrt[M22])^2 ;

plab = 3; S = Mass[[1]]^2+Mass[[2]]^2 + 2Mass[[2]]Sqrt[plab^2+Mass[[1]]^2];
MT = Sqrt[S];
s12max = (MT-Sqrt[M32])^2;
DalitzBorder[plab]=Plot[sijbounds[s,MT,MP,MPI,META],{s,s12min,s12max}, Frame->True ,PlotStyle->{Thick,Black},AxesOrigin->{0,0},FrameLabel->{Style["(Subscript[M, ηπ]^2)",20],Style["(Subscript[M, πp]^2)",20]}];

plab = 6; S = Mass[[1]]^2+Mass[[2]]^2 + 2Mass[[2]]Sqrt[plab^2+Mass[[1]]^2];
MT = Sqrt[S];
s12max = (MT-Sqrt[M32])^2;
DalitzBorder[plab]=Plot[sijbounds[s,MT,MP,MPI,META],{s,s12min,s12max}, Frame->True ,PlotStyle->{Thick,Black},AxesOrigin->{0,0},FrameLabel->{Style["(Subscript[M, ηπ]^2)",20],Style["(Subscript[M, πp]^2)",20]}];

plab = 9; S = Mass[[1]]^2+Mass[[2]]^2 + 2Mass[[2]]Sqrt[plab^2+Mass[[1]]^2];
MT = Sqrt[S];
s12max = (MT-Sqrt[M32])^2;
DalitzBorder[plab]=Plot[sijbounds[s,MT,MP,MPI,META],{s,s12min,s12max}, Frame->True ,PlotStyle->{Thick,Black},AxesOrigin->{0,0},FrameLabel->{Style["(Subscript[M, ηπ]^2)",20],Style["(Subscript[M, πp]^2)",20]}];

plab = 12; S = Mass[[1]]^2+Mass[[2]]^2 + 2Mass[[2]]Sqrt[plab^2+Mass[[1]]^2];
MT = Sqrt[S];
s12max = (MT-Sqrt[M32])^2;
DalitzBorder[plab]=Plot[sijbounds[s,MT,MP,MPI,META],{s,s12min,s12max}, Frame->True ,PlotStyle->{Thick,Black},AxesOrigin->{0,0},FrameLabel->{Style["(Subscript[M, ηπ]^2)",20],Style["(Subscript[M, πp]^2)",20]}];

GraphicsGrid[{{DalitzBorder[3],DalitzBorder[6]},{DalitzBorder[9],DalitzBorder[12]}}]

The border of the Van Hove Plots

S3 = MPI^2+MP^2+2MP Sqrt[3^2+MPI^2];
S6 = MPI^2+MP^2+2MP Sqrt[6^2+MPI^2];
S9 = MPI^2+MP^2+2MP Sqrt[9^2+MPI^2];
S12 = MPI^2+MP^2+2MP Sqrt[12^2+MPI^2];

hexa[3] = ContourPlot[(Sqrt[((2/3)(y)^2)])+(Sqrt[((2/3)((((-1)/2)y-((Sqrt[3])/2)x))^2)])+(Sqrt[((2/3)((((-1)/2)y+((Sqrt[3])/2)x))^2)])==Sqrt[S3],{x,-3,3},{y,-3,3},ContourStyle->{Thick, Blue}];
hexa[6] = ContourPlot[(Sqrt[((2/3)(y)^2)])+(Sqrt[((2/3)((((-1)/2)y-((Sqrt[3])/2)x))^2)])+(Sqrt[((2/3)((((-1)/2)y+((Sqrt[3])/2)x))^2)])==Sqrt[S6],{x,-3,3},{y,-3,3},ContourStyle->{Thick, Blue}];
hexa[9] = ContourPlot[(Sqrt[((2/3)(y)^2)])+(Sqrt[((2/3)((((-1)/2)y-((Sqrt[3])/2)x))^2)])+(Sqrt[((2/3)((((-1)/2)y+((Sqrt[3])/2)x))^2)])==Sqrt[S9],{x,-4,4},{y,-4,4},ContourStyle->{Thick, Blue}];
hexa[12] = ContourPlot[(Sqrt[((2/3)(y)^2)])+(Sqrt[((2/3)((((-1)/2)y-((Sqrt[3])/2)x))^2)])+(Sqrt[((2/3)((((-1)/2)y+((Sqrt[3])/2)x))^2)])==Sqrt[S12],{x,-4,4},{y,-4,4},ContourStyle->{Thick, Blue}];

(* the limit countour with masses *)
m1 = MPI; m2 = META; m3 = MP;
hexar[3] = ContourPlot[(Sqrt[((2/3)(y)^2+m1^2)])+(Sqrt[((2/3)((((-1)/2)y-((Sqrt[3])/2)x))^2+m2^2)])+(Sqrt[((2/3)((((-1)/2)y+((Sqrt[3])/2)x))^2+m3^2)])==Sqrt[S3],{x,-2,2},{y,-2,2},ContourStyle->{Thick, Darker[Green]}];
hexar[6] = ContourPlot[(Sqrt[((2/3)(y)^2+m1^2)])+(Sqrt[((2/3)((((-1)/2)y-((Sqrt[3])/2)x))^2+m2^2)])+(Sqrt[((2/3)((((-1)/2)y+((Sqrt[3])/2)x))^2+m3^2)])==Sqrt[S6],{x,-3,3},{y,-3,3},ContourStyle->{Thick, Darker[Green]}];
hexar[9] = ContourPlot[(Sqrt[((2/3)(y)^2+m1^2)])+(Sqrt[((2/3)((((-1)/2)y-((Sqrt[3])/2)x))^2+m2^2)])+(Sqrt[((2/3)((((-1)/2)y+((Sqrt[3])/2)x))^2+m3^2)])==Sqrt[S9],{x,-4,4},{y,-4,4},ContourStyle->{Thick, Darker[Green]}];
hexar[12] = ContourPlot[(Sqrt[((2/3)(y)^2+m1^2)])+(Sqrt[((2/3)((((-1)/2)y-((Sqrt[3])/2)x))^2+m2^2)])+(Sqrt[((2/3)((((-1)/2)y+((Sqrt[3])/2)x))^2+m3^2)])==Sqrt[S12],{x,-4,4},{y,-4,4},ContourStyle->{Thick, Darker[Green]}];

axe12[3] = Plot[{(Sqrt[3])x,-(Sqrt[3])x},{x,-0.5Sqrt[S3],0.5Sqrt[S3]},Axes->False,PlotStyle->{{Black,Thick},{Black,Thick}}];
axe3[3] = ParametricPlot[{y,0},{y,-Sqrt[S3],Sqrt[S3]},PlotStyle->{Black,Thick}];
axe12[6] = Plot[{(Sqrt[3])x,-(Sqrt[3])x},{x,-0.5Sqrt[S6],0.5Sqrt[S6]},Axes->False,PlotStyle->{{Black,Thick},{Black,Thick}}];
axe3[6] = ParametricPlot[{y,0},{y,-Sqrt[S6],Sqrt[S6]},PlotStyle->{Black,Thick}];
axe12[9] = Plot[{(Sqrt[3])x,-(Sqrt[3])x},{x,-0.5Sqrt[S9],0.5Sqrt[S9]},Axes->False,PlotStyle->{{Black,Thick},{Black,Thick}}];
axe3[9] = ParametricPlot[{y,0},{y,-Sqrt[S9],Sqrt[S9]},PlotStyle->{Black,Thick}];
axe12[12] = Plot[{(Sqrt[3])x,-(Sqrt[3])x},{x,-0.5Sqrt[S12],0.5Sqrt[S12]},Axes->False,PlotStyle->{{Black,Thick},{Black,Thick}}];
axe3[12] = ParametricPlot[{y,0},{y,-Sqrt[S12],Sqrt[S12]},PlotStyle->{Black,Thick}];

vh[3] = Show[hexa[3],hexar[3],axe12[3],axe3[3],Epilog->{Text[Style["p",25],{-0.8,-2}],Text[Style["π",20],{2.2,0.2}],Text[Style["η",20],{-1.4,2}],Text[Style["Plab = 3 GeV",18],{0,2.5}]}];
vh[6] = Show[hexa[6],hexar[6],axe12[6],axe3[6],Epilog->{Text[Style["p",25],{-1.1,-2.5}],Text[Style["π",20],{2.7,0.2}],Text[Style["η",20],{-1.6,2.2}],Text[Style["Plab = 6 GeV",18],{0,2.8}]}];
vh[9] = Show[hexa[9],hexar[9],axe12[9],axe3[9],Epilog->{Text[Style["p",25],{-1.3,-3}],Text[Style["π",20],{3.5,0.4}],Text[Style["η",20],{-2,2.7}],Text[Style["Plab = 9 GeV",18],{0,3.5}]}];vh[12] = Show[hexa[12],hexar[12],axe12[12],axe3[12],Epilog->{Text[Style["p",25],{-1.6,-3.5}],Text[Style["π",20],{3.7,0.4}],Text[Style["η",20],{-2.1,3}],Text[Style["Plab = 12 GeV",18],{0,3.7}]}];


GraphicsGrid[{{vh[3],vh[6]},{vh[9],vh[12]}}]

Events from ROOT

Module

Boost from Lab frame to the center-of-mass frame

Boost[dat_,plab_]:=Block[{
πE,πX,πY,πZ,πEt,πZt,
pE,pX,pY,pZ,pEt,pZt,
ηE,ηX,ηY,ηZ,ηEt,ηZt,bZ,bZt,bE,bEt,tZ,tZt,tEt,tE,
dim,S,Ep,Pp,Elab,Ea,Eb,pa,pb,p1,p2,p3,pt},
dim =Dimensions[dat][[1]];
Elab = Sqrt[plab^2+MPI^2];
S = MP^2+MPI^2+2MP Elab; (* Mandelstam S variable *)
Pp = Sqrt[(S-(MPI+MP)^2)(S-(MPI-MP)^2)/S]/2; (* c.o.m. initial breakup mom. betweem π and target *)
Ea = Sqrt[Pp^2+MPI^2]; (* energy of beam in the c.o.m. frame *)
Eb = Sqrt[Pp^2+MP^2]; (* energy of target in the c.o.m. frame *)
pa = {0,0,+Sqrt[Ea^2-MPI^2], Ea}; (* four vector of a in c.o.m. frame *)
pb = {0,0,-Sqrt[Eb^2-MP^2  ], Eb}; (* four vector of a in c.o.m. frame *)

ηZ = ConstantArray[0,dim];
πZ = ConstantArray[0,dim];
pZ = ConstantArray[0,dim];

Do[
{p1,p2,p3}=dat[[i]];
pt=pb/Eb-pa/Ea;
ηZ[[i]]=p1[[4]]pt[[4]]-(p1[[1]]pt[[1]]+p1[[2]]pt[[2]]+p1[[3]]pt[[3]]);
πZ[[i]]=p2[[4]]pt[[4]]-(p2[[1]]pt[[1]]+p2[[2]]pt[[2]]+p2[[3]]pt[[3]]);
pZ[[i]]=p3[[4]]pt[[4]]-(p3[[1]]pt[[1]]+p3[[2]]pt[[2]]+p3[[3]]pt[[3]]);
,{i,dim}];
ηZ=ηZ (Ea Eb)/Pp/(Ea+Eb);
πZ=πZ (Ea Eb)/Pp/(Ea+Eb);
pZ=pZ (Ea Eb)/Pp/(Ea+Eb);
{ηZ,πZ,pZ}]

Load data

(* CHANGE according to your working directory *)
SetDirectory["/Users/vincent/Dropbox/mes articles/VanHovePlot/Model2"];
dat03=Import["events03.dat", "Table"];Dimensions[dat03]
dat06=Import["events06.dat", "Table"];Dimensions[dat06]
dat09=Import["events09.dat", "Table"];Dimensions[dat09]
dat12=Import["events12.dat", "Table"];Dimensions[dat12]

from (px,py,pz,E) to (E,px,py,pz)

dat = dat12; dim12 = Dimensions[dat][[1]];
Allevents12 = Table[{dat[[i,{2,3,4,1}]],dat[[i,{6,7,8,5}]],dat[[i,{10,11,12,9}]]},{i,dim12}];
dat = dat09; dim09 = Dimensions[dat][[1]];
Allevents09 = Table[{dat[[i,{2,3,4,1}]],dat[[i,{6,7,8,5}]],dat[[i,{10,11,12,9}]]},{i,dim09}];
dat = dat06; dim06 = Dimensions[dat][[1]];
Allevents06 = Table[{dat[[i,{2,3,4,1}]],dat[[i,{6,7,8,5}]],dat[[i,{10,11,12,9}]]},{i,dim06}];
dat = dat03; dim03 = Dimensions[dat][[1]];
Allevents03 = Table[{dat[[i,{2,3,4,1}]],dat[[i,{6,7,8,5}]],dat[[i,{10,11,12,9}]]},{i,dim03}];

Dalitz Plots

tabDalitzAll03=ToDalitz[Allevents03,dim03];
FlatDalitz[03] =Show[DalitzBorder[03],DensityHistogram[tabDalitzAll03[[All,1;;2]],ColorFunction->"Rainbow"],DalitzBorder[03]];
tabDalitzAll06=ToDalitz[Allevents06,dim06];
FlatDalitz[06] =Show[DalitzBorder[06],DensityHistogram[tabDalitzAll06[[All,1;;2]],ColorFunction->"Rainbow"],DalitzBorder[06]];
tabDalitzAll09=ToDalitz[Allevents09,dim09];
FlatDalitz[09] =Show[DalitzBorder[09],DensityHistogram[tabDalitzAll09[[All,1;;2]],ColorFunction->"Rainbow"],DalitzBorder[09]];
tabDalitzAll12=ToDalitz[Allevents12,dim12];
FlatDalitz[12] =Show[DalitzBorder[12],DensityHistogram[tabDalitzAll12[[All,1;;2]],ColorFunction->"Rainbow"],DalitzBorder[12]];

DalitzAll=GraphicsGrid[{{FlatDalitz[03],FlatDalitz[06]},{FlatDalitz[09],FlatDalitz[12]}}]
(*Export["DalitzAll.pdf",DalitzAll,ImageSize->Scaled[2.0]];*)

Van Hove Plots

funω[q_,ηZ_,πZ_,pZ_]:=ArcTan[(pZ-ηZ)/Sqrt[3],πZ];

 {ηZ,πZ,pZ} = Boost[Allevents03,3.0];
q03 = Sqrt[πZ^2 + ηZ^2+pZ^2];
ω03 = Table[funω[q03[[i]],ηZ[[i]],πZ[[i]],pZ[[i]]],{i,dim03}];
vhflatA03 = Table[{ω03[[i]],q03[[i]]},{i,dim03}];
vhflatB03= Table[{vhflatA03[[i,2]] Cos[vhflatA03[[i,1]]],vhflatA03[[i,2]] Sin[vhflatA03[[i,1]]]},{i,dim03}];

 {ηZ,πZ,pZ} = Boost[Allevents06,6.0];
q06 = Sqrt[πZ^2 + ηZ^2+pZ^2];
ω06 = Table[funω[q06[[i]],ηZ[[i]],πZ[[i]],pZ[[i]]],{i,dim06}];
vhflatA06 = Table[{ω06[[i]],q06[[i]]},{i,dim06}];
vhflatB06= Table[{vhflatA06[[i,2]] Cos[vhflatA06[[i,1]]],vhflatA06[[i,2]] Sin[vhflatA06[[i,1]]]},{i,dim06}];

 {ηZ,πZ,pZ} = Boost[Allevents09,9.0];
q09 = Sqrt[πZ^2 + ηZ^2+pZ^2];
ω09 = Table[funω[q09[[i]],ηZ[[i]],πZ[[i]],pZ[[i]]],{i,dim09}];
vhflatA09 = Table[{ω09[[i]],q09[[i]]},{i,dim09}];
vhflatB09= Table[{vhflatA09[[i,2]] Cos[vhflatA09[[i,1]]],vhflatA09[[i,2]] Sin[vhflatA09[[i,1]]]},{i,dim09}];


 {ηZ,πZ,pZ} = Boost[Allevents12,12.0];
q12 = Sqrt[πZ^2 + ηZ^2+pZ^2];
ω12 = Table[funω[q12[[i]],ηZ[[i]],πZ[[i]],pZ[[i]]],{i,dim12}];
vhflatA12 = Table[{ω12[[i]],q12[[i]]},{i,dim12}];
vhflatB12= Table[{vhflatA12[[i,2]] Cos[vhflatA12[[i,1]]],vhflatA12[[i,2]] Sin[vhflatA12[[i,1]]]},{i,dim12}];

lplotB03 =DensityHistogram[vhflatB03,ColorFunction->"SolarColors"];
plotB[3]=Show[vh[3],lplotB03 ];

lplotB06 =DensityHistogram[vhflatB06,ColorFunction->"SolarColors"];
plotB[6]=Show[vh[6],lplotB06 ];

lplotB09 =DensityHistogram[vhflatB09,ColorFunction->"SolarColors"];
plotB[9]=Show[vh[9],lplotB09 ];

lplotB12 =DensityHistogram[vhflatB12,ColorFunction->"SolarColors"];
plotB[12]=Show[vh[12],lplotB12 ];

VanHoveAll=GraphicsGrid[{{ plotB[3], plotB[6]},{ plotB[9], plotB[12]}}]
(*Export["VanHoveAll.pdf",VanHoveAll,ImageSize->Scaled[2.0]];*)

Projections

framestyle = {Frame->True,FrameStyle->15};(* ChartStyle->EdgeForm[None] *)
dbin=0.05;
bspec03a ={Range[(META+MPI)^2,(Sqrt[S3]-MP)^2,dbin]}; 
bspec03b ={Range[(MP+MPI)^2,(Sqrt[S3]-META)^2,dbin]}; 
bspec03c ={Range[(META+MP)^2,(Sqrt[S3]-MPI)^2,dbin]}; 

projA[03,1]=Histogram[tabDalitzAll03[[All,1]],bspec03a,ChartStyle->Blue,framestyle,FrameLabel->Style["(Subscript[M, ηπ]^2)",20]]
projA[03,2]=Histogram[tabDalitzAll03[[All,2]],bspec03b,ChartStyle->Blue,framestyle,FrameLabel->Style["(Subscript[M, πp]^2)",20]]
projA[03,3]=Histogram[tabDalitzAll03[[All,3]],bspec03c,ChartStyle->Blue,framestyle,FrameLabel->Style["(Subscript[M, ηp]^2)",20]]

framestylea = {Frame->True,FrameStyle->15, PlotRange->{{0,6},All},AxesOrigin->{0,0}};
framestyleb = {Frame->True,FrameStyle->15, PlotRange->{{1,6},All},AxesOrigin->{0,0}};framestylec = {Frame->True,FrameStyle->15, PlotRange->{{2,6},All},AxesOrigin->{0,0}};

dbin=0.050; (* bin size: 50 MeV *)
bspec03a ={Range[(META+MPI)^2,(Sqrt[S3]-MP)^2,dbin]}; 
bspec03b ={Range[(MP+MPI)^2,(Sqrt[S3]-META)^2,dbin]}; 
bspec03c ={Range[(META+MP)^2,(Sqrt[S3]-MPI)^2,dbin]}; 

bspec06a ={Range[(META+MPI)^2,(Sqrt[S6]-MP)^2,dbin]}; 
bspec06b ={Range[(MP+MPI)^2,(Sqrt[S6]-META)^2,dbin]}; 
bspec06c ={Range[(META+MP)^2,(Sqrt[S6]-MPI)^2,dbin]}; 

bspec09a ={Range[(META+MPI)^2,(Sqrt[S9]-MP)^2,dbin]}; 
bspec09b ={Range[(MP+MPI)^2,(Sqrt[S9]-META)^2,dbin]}; 
bspec09c ={Range[(META+MP)^2,(Sqrt[S9]-MPI)^2,dbin]}; 

bspec12a ={Range[(META+MPI)^2,(Sqrt[S12]-MP)^2,dbin]}; 
bspec12b ={Range[(MP+MPI)^2,(Sqrt[S12]-META)^2,dbin]}; 
bspec12c ={Range[(META+MP)^2,(Sqrt[S12]-MPI)^2,dbin]}; 

proj0[03,1]=Histogram[tabDalitzAll03[[All,1]],bspec03a,ChartStyle->Blue,framestylea,FrameLabel->Style["(Subscript[M, ηπ]^2)",20]];
proj0[03,2]=Histogram[tabDalitzAll03[[All,2]],bspec03b,ChartStyle->Blue,framestyleb,FrameLabel->Style["(Subscript[M, πp]^2)",20]];
proj0[03,3]=Histogram[tabDalitzAll03[[All,3]],bspec03c,ChartStyle->Blue,framestylec,FrameLabel->Style["(Subscript[M, ηp]^2)",20]];

proj0[06,1]=Histogram[tabDalitzAll06[[All,1]],bspec06a,ChartStyle->Blue,framestylea,FrameLabel->Style["(Subscript[M, ηπ]^2)",20]];
proj0[06,2]=Histogram[tabDalitzAll06[[All,2]],bspec06b,ChartStyle->Blue,framestyleb,FrameLabel->Style["(Subscript[M, πp]^2)",20]];
proj0[06,3]=Histogram[tabDalitzAll06[[All,3]],bspec06c,ChartStyle->Blue,framestylec,FrameLabel->Style["(Subscript[M, ηp]^2)",20]];

proj0[09,1]=Histogram[tabDalitzAll09[[All,1]],bspec09a,ChartStyle->Blue,framestylea,FrameLabel->Style["(Subscript[M, ηπ]^2)",20]];
proj0[09,2]=Histogram[tabDalitzAll09[[All,2]],bspec09b,ChartStyle->Blue,framestyleb,FrameLabel->Style["(Subscript[M, πp]^2)",20]];
proj0[09,3]=Histogram[tabDalitzAll09[[All,3]],bspec09c,ChartStyle->Blue,framestylec,FrameLabel->Style["(Subscript[M, ηp]^2)",20]];

proj0[12,1]=Histogram[tabDalitzAll12[[All,1]],bspec12a,ChartStyle->Blue,framestylea,FrameLabel->Style["(Subscript[M, ηπ]^2)",20]];
proj0[12,2]=Histogram[tabDalitzAll12[[All,2]],bspec12b,ChartStyle->Blue,framestyleb,FrameLabel->Style["(Subscript[M, πp]^2)",20]];
proj0[12,3]=Histogram[tabDalitzAll12[[All,3]],bspec12c,ChartStyle->Blue,framestylec,FrameLabel->Style["(Subscript[M, ηp]^2)",20]];

proj0All=GraphicsGrid[{{proj0[03,1],proj0[03,2],proj0[03,3]},{proj0[06,1],proj0[06,2],proj0[06,3]},
{proj0[09,1],proj0[09,2],proj0[09,3]},{proj0[12,1],proj0[12,2],proj0[12,3]}}]
Export["proj0All.pdf",proj0All,ImageSize->Scaled[2.0]];

Cut in Van Hove Angle

cut A the triangle  (η , π) > 0

indexA03 = Position[ω03,_?(π>#>2π/3&)]//Flatten; 
indexA06 = Position[ω06,_?(π>#>2π/3&)]//Flatten; 
indexA09 = Position[ω09,_?(π>#>2π/3&)]//Flatten; 
indexA12 = Position[ω12,_?(π>#>2π/3&)]//Flatten; 

cutA03 = Allevents03[[indexA03,All,All]];
cutA06 = Allevents06[[indexA06,All,All]];
cutA09 = Allevents09[[indexA09,All,All]];
cutA12 = Allevents12[[indexA12,All,All]];

{ηZ,πZ,pZ} = Boost[cutA12,12.0];
xq12 = Sqrt[πZ^2 + ηZ^2+pZ^2];
xω12 = Table[funω[xq12[[i]],ηZ[[i]],πZ[[i]],pZ[[i]]],{i,Length[indexA12]}];
xvhflatA12 = Table[{xω12[[i]],xq12[[i]]},{i,Length[indexA12]}];
xvhflatB12= Table[{xvhflatA12[[i,2]] Cos[xvhflatA12[[i,1]]],xvhflatA12[[i,2]] Sin[xvhflatA12[[i,1]]]},{i,Length[indexA12]}];
xlplotB12 =DensityHistogram[xvhflatB12,ColorFunction->"SolarColors"];
Show[vh[12],xlplotB12 ]

tabDalitzCutA03=ToDalitz[cutA03,Length[indexA03]];
CutADalitz[03] =Show[DalitzBorder[03],DensityHistogram[tabDalitzCutA03[[All,1;;2]],ColorFunction->"Rainbow"],DalitzBorder[03]];

tabDalitzCutA06=ToDalitz[cutA06,Length[indexA06]];
CutADalitz[06] =Show[DalitzBorder[06],DensityHistogram[tabDalitzCutA06[[All,1;;2]],ColorFunction->"Rainbow"],DalitzBorder[06]];

tabDalitzCutA09=ToDalitz[cutA09,Length[indexA09]];
CutADalitz[09] =Show[DalitzBorder[09],DensityHistogram[tabDalitzCutA09[[All,1;;2]],ColorFunction->"Rainbow"],DalitzBorder[09]];

tabDalitzCutA12=ToDalitz[cutA12,Length[indexA12]];
CutADalitz[12] =Show[DalitzBorder[12],DensityHistogram[tabDalitzCutA12[[All,1;;2]],ColorFunction->"Rainbow"],DalitzBorder[12]];

DalitzCutA=GraphicsGrid[{{CutADalitz[03],CutADalitz[06]},{CutADalitz[09],CutADalitz[12]}}]
(*Export["DalitzCutA.pdf",DalitzCutA,ImageSize->Scaled[2.0]];*)

projA[03,1]=Histogram[tabDalitzCutA03[[All,1]],bspec03a,ChartStyle->Red];
projA[03,2]=Histogram[tabDalitzCutA03[[All,2]],bspec03b,ChartStyle->Red];
projA[03,3]=Histogram[tabDalitzCutA03[[All,3]],bspec03c,ChartStyle->Red];

projA[06,1]=Histogram[tabDalitzCutA06[[All,1]],bspec06a,ChartStyle->Red];
projA[06,2]=Histogram[tabDalitzCutA06[[All,2]],bspec06b,ChartStyle->Red];
projA[06,3]=Histogram[tabDalitzCutA06[[All,3]],bspec06c,ChartStyle->Red];

projA[09,1]=Histogram[tabDalitzCutA09[[All,1]],bspec09a,ChartStyle->Red];
projA[09,2]=Histogram[tabDalitzCutA09[[All,2]],bspec09b,ChartStyle->Red];
projA[09,3]=Histogram[tabDalitzCutA09[[All,3]],bspec09c,ChartStyle->Red];

projA[12,1]=Histogram[tabDalitzCutA12[[All,1]],bspec12a,ChartStyle->Red];
projA[12,2]=Histogram[tabDalitzCutA12[[All,2]],bspec12b,ChartStyle->Red];
projA[12,3]=Histogram[tabDalitzCutA12[[All,3]],bspec12c,ChartStyle->Red];

projCutA=GraphicsGrid[{
{Show[proj0[03,1],projA[03,1]],Show[proj0[03,2],projA[03,2]],Show[proj0[03,3],projA[03,3]]},{Show[proj0[06,1],projA[06,1]],Show[proj0[06,2],projA[06,2]],Show[proj0[06,3],projA[06,3]]},
{Show[proj0[09,1],projA[09,1]],Show[proj0[09,2],projA[09,2]],Show[proj0[09,3],projA[09,3]]},
{Show[proj0[12,1],projA[12,1]],Show[proj0[12,2],projA[12,2]],Show[proj0[12,3],projA[12,3]]}}]
(*Export["projCutA.pdf",projCutA,ImageSize->Scaled[2.0]];*)

cut B the triangle  (η , p) < 0

indexB03 = Position[ω03,_?(2π/3>#>π/3&)]//Flatten; 
indexB06 = Position[ω06,_?(2π/3>#>π/3&)]//Flatten; 
indexB09 = Position[ω09,_?(2π/3>#>π/3&)]//Flatten; 
indexB12 = Position[ω12,_?(2π/3>#>π/3&)]//Flatten; 

cutB03 = Allevents03[[indexB03,All,All]];
cutB06 = Allevents06[[indexB06,All,All]];
cutB09 = Allevents09[[indexB09,All,All]];
cutB12 = Allevents12[[indexB12,All,All]];

{ηZ,πZ,pZ} = Boost[cutB12,12.0];
xq12 = Sqrt[πZ^2 + ηZ^2+pZ^2];
xω12 = Table[funω[xq12[[i]],ηZ[[i]],πZ[[i]],pZ[[i]]],{i,Length[indexB12]}];
xvhflatA12 = Table[{xω12[[i]],xq12[[i]]},{i,Length[indexB12]}];
xvhflatB12= Table[{xvhflatA12[[i,2]] Cos[xvhflatA12[[i,1]]],xvhflatA12[[i,2]] Sin[xvhflatA12[[i,1]]]},{i,Length[indexB12]}];
xlplotB12 =DensityHistogram[xvhflatB12,ColorFunction->"SolarColors"];
Show[vh[12],xlplotB12 ]

tabDalitzCutB03=ToDalitz[cutB03,Length[indexB03]];
CutBDalitz[03] =Show[DalitzBorder[03],DensityHistogram[tabDalitzCutB03[[All,1;;2]],ColorFunction->"Rainbow"],DalitzBorder[03]];

tabDalitzCutB06=ToDalitz[cutB06,Length[indexB06]];
CutBDalitz[06] =Show[DalitzBorder[06],DensityHistogram[tabDalitzCutB06[[All,1;;2]],ColorFunction->"Rainbow"],DalitzBorder[06]];

tabDalitzCutB09=ToDalitz[cutB09,Length[indexB09]];
CutBDalitz[09] =Show[DalitzBorder[09],DensityHistogram[tabDalitzCutB09[[All,1;;2]],ColorFunction->"Rainbow"],DalitzBorder[09]];

tabDalitzCutB12=ToDalitz[cutB12,Length[indexB12]];
CutBDalitz[12] =Show[DalitzBorder[12],DensityHistogram[tabDalitzCutB12[[All,1;;2]],ColorFunction->"Rainbow"],DalitzBorder[12]];

DalitzCutB=GraphicsGrid[{{CutBDalitz[03],CutBDalitz[06]},{CutBDalitz[09],CutBDalitz[12]}}]
(*Export["DalitzCutB.pdf",DalitzCutB,ImageSize->Scaled[2.0]];*)

projB[03,1]=Histogram[tabDalitzCutB03[[All,1]],bspec03a,ChartStyle->Darker[Green]];
projB[03,2]=Histogram[tabDalitzCutB03[[All,2]],bspec03b,ChartStyle->Darker[Green]];
projB[03,3]=Histogram[tabDalitzCutB03[[All,3]],bspec03c,ChartStyle->Darker[Green]];

projB[06,1]=Histogram[tabDalitzCutB06[[All,1]],bspec06a,ChartStyle->Darker[Green]];
projB[06,2]=Histogram[tabDalitzCutB06[[All,2]],bspec06b,ChartStyle->Darker[Green]];
projB[06,3]=Histogram[tabDalitzCutB06[[All,3]],bspec06c,ChartStyle->Darker[Green]];

projB[09,1]=Histogram[tabDalitzCutB09[[All,1]],bspec09a,ChartStyle->Darker[Green]];
projB[09,2]=Histogram[tabDalitzCutB09[[All,2]],bspec09b,ChartStyle->Darker[Green]];
projB[09,3]=Histogram[tabDalitzCutB09[[All,3]],bspec09c,ChartStyle->Darker[Green]];

projB[12,1]=Histogram[tabDalitzCutB12[[All,1]],bspec12a,ChartStyle->Darker[Green]];
projB[12,2]=Histogram[tabDalitzCutB12[[All,2]],bspec12b,ChartStyle->Darker[Green]];
projB[12,3]=Histogram[tabDalitzCutB12[[All,3]],bspec12c,ChartStyle->Darker[Green]];

projCutB=GraphicsGrid[{
{Show[proj0[03,1],projB[03,1]],Show[proj0[03,2],projB[03,2]],Show[proj0[03,3],projB[03,3]]},{Show[proj0[06,1],projB[06,1]],Show[proj0[06,2],projB[06,2]],Show[proj0[06,3],projB[06,3]]},
{Show[proj0[09,1],projB[09,1]],Show[proj0[09,2],projB[09,2]],Show[proj0[09,3],projB[09,3]]},
{Show[proj0[12,1],projB[12,1]],Show[proj0[12,2],projB[12,2]],Show[proj0[12,3],projB[12,3]]}}]
(*Export["projCutB.pdf",projCutB,ImageSize->Scaled[2.0]];*)

cut C the triangle  (π , p) < 0

indexC03 = Position[ω03,_?(-2π/3>#>-π&)]//Flatten; 
indexC06 = Position[ω06,_?(-2π/3>#>-π&)]//Flatten; 
indexC09 = Position[ω09,_?(-2π/3>#>-π&)]//Flatten; 
indexC12 = Position[ω12,_?(-2π/3>#>-π&)]//Flatten; 

cutC03 = Allevents03[[indexC03,All,All]];
cutC06 = Allevents06[[indexC06,All,All]];
cutC09 = Allevents09[[indexC09,All,All]];
cutC12 = Allevents12[[indexC12,All,All]];

{ηZ,πZ,pZ} = Boost[cutC12,12.0];
xq12 = Sqrt[πZ^2 + ηZ^2+pZ^2];
xω12 = Table[funω[xq12[[i]],ηZ[[i]],πZ[[i]],pZ[[i]]],{i,Length[indexC12]}];
xvhflatA12 = Table[{xω12[[i]],xq12[[i]]},{i,Length[indexC12]}];
xvhflatB12= Table[{xvhflatA12[[i,2]] Cos[xvhflatA12[[i,1]]],xvhflatA12[[i,2]] Sin[xvhflatA12[[i,1]]]},{i,Length[indexC12]}];
xlplotB12 =DensityHistogram[xvhflatB12,ColorFunction->"SolarColors"];
Show[vh[12],xlplotB12 ]

tabDalitzCutC03=ToDalitz[cutC03,Length[indexC03]];
CutCDalitz[03] =Show[DalitzBorder[03],DensityHistogram[tabDalitzCutC03[[All,1;;2]],ColorFunction->"Rainbow"],DalitzBorder[03]];

tabDalitzCutC06=ToDalitz[cutC06,Length[indexC06]];
CutCDalitz[06] =Show[DalitzBorder[06],DensityHistogram[tabDalitzCutC06[[All,1;;2]],ColorFunction->"Rainbow"],DalitzBorder[06]];

tabDalitzCutC09=ToDalitz[cutC09,Length[indexC09]];
CutCDalitz[09] =Show[DalitzBorder[09],DensityHistogram[tabDalitzCutC09[[All,1;;2]],ColorFunction->"Rainbow"],DalitzBorder[09]];

tabDalitzCutC12=ToDalitz[cutC12,Length[indexC12]];
CutCDalitz[12] =Show[DalitzBorder[12],DensityHistogram[tabDalitzCutC12[[All,1;;2]],ColorFunction->"Rainbow"],DalitzBorder[12]];

DalitzCutC=GraphicsGrid[{{CutCDalitz[03],CutCDalitz[06]},{CutCDalitz[09],CutCDalitz[12]}}]
(*Export["DalitzCutC.pdf",DalitzCutC,ImageSize->Scaled[2.0]];*)

projC[03,1]=Histogram[tabDalitzCutC03[[All,1]],bspec03a,ChartStyle->Yellow];
projC[03,2]=Histogram[tabDalitzCutC03[[All,2]],bspec03b,ChartStyle->Yellow];
projC[03,3]=Histogram[tabDalitzCutC03[[All,3]],bspec03c,ChartStyle->Yellow];

projC[06,1]=Histogram[tabDalitzCutC06[[All,1]],bspec06a,ChartStyle->Yellow];
projC[06,2]=Histogram[tabDalitzCutC06[[All,2]],bspec06b,ChartStyle->Yellow];
projC[06,3]=Histogram[tabDalitzCutC06[[All,3]],bspec06c,ChartStyle->Yellow];

projC[09,1]=Histogram[tabDalitzCutC09[[All,1]],bspec09a,ChartStyle->Yellow];
projC[09,2]=Histogram[tabDalitzCutC09[[All,2]],bspec09b,ChartStyle->Yellow];
projC[09,3]=Histogram[tabDalitzCutC09[[All,3]],bspec09c,ChartStyle->Yellow];

projC[12,1]=Histogram[tabDalitzCutC12[[All,1]],bspec12a,ChartStyle->Yellow];
projC[12,2]=Histogram[tabDalitzCutC12[[All,2]],bspec12b,ChartStyle->Yellow];
projC[12,3]=Histogram[tabDalitzCutC12[[All,3]],bspec12c,ChartStyle->Yellow];

projCutC=GraphicsGrid[{
{Show[proj0[03,1],projC[03,1]],Show[proj0[03,2],projC[03,2]],Show[proj0[03,3],projC[03,3]]},{Show[proj0[06,1],projC[06,1]],Show[proj0[06,2],projC[06,2]],Show[proj0[06,3],projC[06,3]]},
{Show[proj0[09,1],projC[09,1]],Show[proj0[09,2],projC[09,2]],Show[proj0[09,3],projC[09,3]]},
{Show[proj0[12,1],projC[12,1]],Show[proj0[12,2],projC[12,2]],Show[proj0[12,3],projC[12,3]]}}]
(*Export["projCutC.pdf",projCutC,ImageSize->Scaled[2.0]];*)