diff --git a/Examples/Example1/Analysis.cxx b/Examples/Example1/Analysis.cxx
index fa8554aebb8d4a91f48286322f34a9f005aa5fbf..0a60fe37dc998bf3977f6574983f3a109b90a94f 100644
--- a/Examples/Example1/Analysis.cxx
+++ b/Examples/Example1/Analysis.cxx
@@ -75,6 +75,9 @@ void Analysis::TreatEvent(){
   OriginalELab = ReactionConditions->GetKineticEnergy(0);
   OriginalThetaLab = ReactionConditions->GetTheta(0);
   OriginalBeamEnergy = ReactionConditions->GetBeamEnergy();
+  ReactionVertexX = ReactionConditions->GetVertexPositionX(); 
+  ReactionVertexY = ReactionConditions->GetVertexPositionY(); 
+  ReactionVertexZ = ReactionConditions->GetVertexPositionZ(); 
   // Get the Init information on beam position and energy
   // and apply by hand the experimental resolution
   // This is because the beam diagnosis are not simulated
@@ -207,6 +210,9 @@ void Analysis::InitOutputBranch() {
   RootOutput::getInstance()->GetTree()->Branch("OriginalELab",&OriginalELab,"OriginalELab/D");
   RootOutput::getInstance()->GetTree()->Branch("OriginalThetaLab",&OriginalThetaLab,"OriginalThetaLab/D");
   RootOutput::getInstance()->GetTree()->Branch("OriginalBeamEnergy",&OriginalBeamEnergy,"OriginalBeamEnergy/D");
+  RootOutput::getInstance()->GetTree()->Branch("ReactionVertexX",&ReactionVertexX,"ReactionVertexX/D");
+  RootOutput::getInstance()->GetTree()->Branch("ReactionVertexY",&ReactionVertexY,"ReactionVertexY/D");
+  RootOutput::getInstance()->GetTree()->Branch("ReactionVertexZ",&ReactionVertexZ,"ReactionVertexZ/D");
 }
 
 ////////////////////////////////////////////////////////////////////////////////
diff --git a/Examples/Example1/Analysis.h b/Examples/Example1/Analysis.h
index 908a9f782c1812208524ef2cab6eae4dde6b4de7..7c2408efc4df64191df2b28f36793f17c080017b 100644
--- a/Examples/Example1/Analysis.h
+++ b/Examples/Example1/Analysis.h
@@ -54,6 +54,9 @@ class Analysis: public NPL::VAnalysis{
     double OriginalELab;
     double OriginalThetaLab;
     double OriginalBeamEnergy;
+    double ReactionVertexX;
+    double ReactionVertexY;
+    double ReactionVertexZ;
     NPL::Reaction* He10Reaction;
 
     // intermediate variable
diff --git a/Examples/Example1/ShowResults.C b/Examples/Example1/ShowResults.C
index 8ea68944e989e3f7845c980a10f7837f31c3488e..3c8f9dfa9bdf91acea8e3538dee8f460449ab6a3 100644
--- a/Examples/Example1/ShowResults.C
+++ b/Examples/Example1/ShowResults.C
@@ -86,7 +86,7 @@ void ShowResults(){
     f->SetNpx(1000);
     
     TCanvas* c2 = new TCanvas("Simulated","Simulated",600,0,600,600);
-    c2->Divide(2,2);
+    c2->Divide(2,3);
     
     c2->cd(1);
     chain->Draw("OriginalELab:OriginalThetaLab>>hS(1000,0,90,1000,0,30)","","col");
@@ -110,6 +110,7 @@ void ShowResults(){
     
     TLine* lT = new TLine(0,0,90,90);
     lT->Draw();
+
     c2->cd(4);
     chain->Draw("OriginalBeamEnergy:BeamEnergy>>hS4(1000,500,600,1000,500,600)","BeamEnergy>0","col");
     TH1F* hS4 = (TH1F*) gDirectory->FindObjectAny("hS4"); 
@@ -120,5 +121,16 @@ void ShowResults(){
     lB->Draw();
     
 
+    c2->cd(5);
+    chain->Draw("ReactionVertexY:ReactionVertexX>>hVXY(1000,-20,20,1000,-20,20)","","col");
+    TH2F* hVXY = (TH2F*) gDirectory->FindObjectAny("hVXY"); 
+    hVXY->GetYaxis()->SetTitle("Reaction vertex Y (mm)");
+    hVXY->GetXaxis()->SetTitle("Reaction vertex X (mm)");  
+
+    c2->cd(6);
+    chain->Draw("ReactionVertexX:ReactionVertexZ*1000.>>hVXZ(1000,-15,15,1000,-20,20)","","col");
+    TH2F* hVXZ = (TH2F*) gDirectory->FindObjectAny("hVXZ"); 
+    hVXZ->GetZaxis()->SetTitle("Reaction vertex X (mm)");
+    hVXZ->GetXaxis()->SetTitle("Reaction vertex Z (um)");  
 
 }
diff --git a/Examples/Example1/sim.sh b/Examples/Example1/sim.sh
new file mode 100755
index 0000000000000000000000000000000000000000..561a386a8ee24fc95ae23f21b4c9e7881b5c3a76
--- /dev/null
+++ b/Examples/Example1/sim.sh
@@ -0,0 +1,3 @@
+npsimulation -D Example1.detector -E Example1.reaction -O Example1 -B run.mac
+npanalysis --last-sim  -O Example1
+root ShowResults.C