Index: /trunk/Resolutions.cpp
===================================================================
--- /trunk/Resolutions.cpp	(revision 22)
+++ /trunk/Resolutions.cpp	(revision 23)
@@ -172,4 +172,5 @@
   vector<fastjet::PseudoJet> inclusive_jetsS;
   vector<fastjet::PseudoJet> sorted_jetsS;
+  vector<PhysicsTower> towers;
 
   fastjet::JetDefinition jet_def;
@@ -180,5 +181,5 @@
   // set up a CDF midpoint jet definition
   #ifdef ENABLE_PLUGIN_CDFCONES
-  plugins = new fastjet::CDFJetCluPlugin(DET->C_SEEDTHRESHOLD,DET->CONERADIUS,DET->C_ADJACENCYCUT,DET->C_MAXITERATIONS,DET->C_IRATCH,DET->C_OVERLAPTHRESHOLD);
+  plugins = new fastjet::CDFJetCluPlugin(0,DET->CONERADIUS,DET->C_ADJACENCYCUT,DET->C_MAXITERATIONS,DET->C_IRATCH,DET->C_OVERLAPTHRESHOLD);
   jet_def = fastjet::JetDefinition(plugins);
   #else
@@ -200,4 +201,5 @@
   for(entry = 0; entry < allEntries; ++entry)
     {
+      TLorentzVector PTmisS(0,0,0,0);
       TLorentzVector PTmis(0,0,0,0);
       treeReader->ReadEntry(entry);
@@ -216,4 +218,5 @@
       inclusive_jetsS.clear();
       sorted_jetsS.clear();
+      towers.clear();
 
       // Loop over all particles in event
@@ -238,4 +241,5 @@
               )
             ) {
+                if(pid != pMU)PTmis = PTmis + genMomentum;//ptmis
 		switch(pid) {
 		    
@@ -264,8 +268,10 @@
 		// all final particles but muons and neutrinos	
 		// for calorimetric towers and mission PT
-                if(genMomentum.E()!=0) PTmis = PTmis + genMomentum;//ptmis
-
-		if(pid != pMU && genMomentum.Pt() > DET->PT_TRACKS_MIN)
+                //if(genMomentum.E()!=0) PTmis = PTmis + genMomentum;//ptmis
+
+		if(pid != pMU)
                   {
+                     PTmisS = PTmisS + genMomentum;
+                     towers.push_back(PhysicsTower(LorentzVector(genMomentum.Px(),genMomentum.Py(),genMomentum.Pz(), genMomentum.E())));
                      input_particlesS.push_back(fastjet::PseudoJet(genMomentum.Px(),genMomentum.Py(),genMomentum.Pz(), genMomentum.E()));
                   }
@@ -289,4 +295,26 @@
 	} // while
 
+  TLorentzVector toWerS(0,0,0,0);
+  //calcul de ETmis au depart des tours calorimetriques.. 
+/*  for(unsigned int i=0; i < towers.size(); i++) {
+      if(towers[i].fourVector.pt() < 0.5) continue;
+          toWerS.SetPxPyPzE(towers[i].fourVector.px,towers[i].fourVector.py,towers[i].fourVector.pz,towers[i].fourVector.E);
+          PTmisS = PTmisS + toWerS;
+          }
+*/
+//PTmis = (0,0,0,0)-PTmis;
+//PTmisS = (0,0,0,0)-PTmisS;
+
+//float ET=PTmis.Et()<<endl;
+//float ETS=PTmisS.Et()<<endl;
+
+cout<<"Ptmis "<<PTmis.Et()<<" PTmis smeare "<<PTmisS.Et()<<endl;
+cout<<"Ptmis "<<(-PTmis).Pt()<<" PTmis smeare "<<(-PTmisS).Pt()<<endl;
+
+  elementEtmis= (ETMIS*) branchetmis->NewEntry();
+  elementEtmis->Et = (PTmis).Et();
+  elementEtmis->EtSmeare = (PTmisS).Et();
+
+
 
       //*****************************
@@ -314,12 +342,11 @@
         if(JETSm.Pt()>1)
           {
+            float NonSmeareEt=(JET.E()*sin(JET.Theta()));
+            float SmeareEt=(JETSm.E()*sin(JETSm.Theta()));
+
             elementJet= (RESOLJET*) branchjet->NewEntry();
-            elementJet->NonSmearePT = JET.Pt();
-            elementJet->SmearePT = (JETSm.Pt()/JET.Pt());
-            /*cout<<"valeur obtenue "<<JETSm.Pt()/JET.Pt()<<endl;
-            cout<<" pt non smeare "<<JET.Pt()<<" phi "<<JET.Phi()<<" eta "<<JET.Eta()<<endl;
-            cout<<"pt smeare "<<JETSm.Pt()<<" phi "<<JETSm.Phi()<<" eta "<<JETSm.Eta()<<endl;
-            cout<<"************"<<endl;
-*/
+            elementJet->NonSmearePT = JET.Et();
+            elementJet->SmearePT = JETSm.Et()/JET.Et();
+          
           }
 
Index: /trunk/interface/FuncDef.h
===================================================================
--- /trunk/interface/FuncDef.h	(revision 22)
+++ /trunk/interface/FuncDef.h	(revision 23)
@@ -79,8 +79,12 @@
 
   string nom = histo.erase(0,histo.find(">>")+2);
-  TH1F *h = new TH1F(nom.c_str(),"",50,0.5,1.5);
+  TH1F *h = new TH1F(nom.c_str(),"",100,0,2);
   Analyze->Draw(temp.c_str(),mintemp.c_str(),maxtemp.c_str());
   h->SetMarkerSize(0.6);
-  TF1 *Gauss = new TF1("Gauss","gausn",0.8,1.3);
+  double MeanFix = h->GetMean();
+  double RangMin = h->GetMean()-h->GetRMS();
+  double RangMax = h->GetMean()+h->GetRMS();
+  TF1 *Gauss = new TF1("Gauss","gausn",RangMin,RangMax);
+  Gauss->FixParameter(1,MeanFix);
   h->Fit("Gauss","R");
   h->Fit("Gauss","RI");
@@ -92,5 +96,4 @@
   h->GetXaxis()->SetTitle("E_{T}^{rec}/E_{T}^{MC} [GeV]");
 }
-
 
 TH1F * HiggsMakeNormTH1F(Float_t Scale, Int_t numBin,Float_t minX,Float_t maxX, TTree * Analyze,string histo, Int_t Lcolor, Int_t Fcolor, Int_t Lstyle,Int_t width=2,bool Log=false)
Index: /trunk/routines/resolutions.C
===================================================================
--- /trunk/routines/resolutions.C	(revision 22)
+++ /trunk/routines/resolutions.C	(revision 23)
@@ -9,40 +9,26 @@
 #include "interface/FuncDef.h"
 
-void Plot()
+void JetResol(TTree *Analyze)
 {
-
-  setTDRStyle();
-  gROOT->Reset();
-
-  TFile *f1 = new TFile("Smearing.root","read");
-  TTree *Analyze = (TTree*)f1->Get("Analysis");
-
-
-//$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$
-//$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$$
-
-   TCanvas *c1 = new TCanvas("c1","gerrors2",100,100,600,450);
+   TCanvas *c1 = new TCanvas("c1","JET resol",100,100,600,450);
    c1->cd(); 
   
-   const Int_t numBin=5;
-//   double bins[numBin]={5,10,25,55,85,115,205,400};
-   double bins[5]={10,30,50,80,120};
-   TProfile *ETSoverET = new TProfile("ETSoverET","Resolution des electrons",(numBin-1),bins,-10,10);
+   const Int_t numBin=7;
+   double bins[numBin]={0,20,60,100,220,420,860};
+   TProfile *ETSoverET = new TProfile("ETSoverET","Jet resolution: E_{T}^{rec}/E_{T}^{mc}",(numBin-1),bins,-10,10);
    Analyze->Draw("JetPTResol.SmearePT:JetPTResol.NonSmearePT>>ETSoverET","JetPTResol.NonSmearePT>1");
 
    double rms[numBin];
    double mean[numBin];
-   //TF1 *Gauss = new TF1("Gauss","gaus",0,4);
 
-   TCanvas *c2 = new TCanvas("c2","gerrors2",0,0,1000,650);
+   TCanvas *c2 = new TCanvas("c2","JET resol",0,0,1000,650);
    c2->cd(); int frame=0;
    c2->Divide(6,2);
-//   c2->Divide(6,2);
 
-   float x[numBin];
-   float y[numBin];
+   float x[numBin-1];
+   float y[numBin-1];
    float finval=0;//valeur a remplir puis a fitter
 
-   for ( int i=0; i<numBin; i++)  // premiÃšre bin : i ==1 et pas i == 0
+   for ( int i=0; i<(numBin-1); i++)  // premiÃšre bin : i ==1 et pas i == 0
       {
         TAxis *xaxis = ETSoverET->GetXaxis();
@@ -64,20 +50,16 @@
         GaussValues(Analyze,tempName,rms[i],mean[i],mystringMin,mystringMax);
         x[i]=binCenter;
-        mean[i]=ETSoverET->GetBinContent(i+1);
-   //     rms[i]=ETSoverET->GetBinError(i+1);
         finval=rms[i]/mean[i];
-        y[i]=finval;
+        y[i]=(finval*2);
         cout<<"finval "<<finval<<" mean val "<<mean[i]<<" rms value "<<rms[i]<<endl;
 
       }
 
-y[0]=0.2;
-
-   TCanvas *c3 = new TCanvas("c3","gerrors2",100,100,600,450);
+   TCanvas *c3 = new TCanvas("c3","JET resol",100,100,600,450);
    c3->cd();
 
-   TF1 *fitfun = new TF1("user","sqrt(pow([0],2)+pow([1]/sqrt(x),2)+pow([2]/x,2))",10,500);
+   TF1 *fitfun = new TF1("user","sqrt(pow([0],2)+pow([1]/sqrt(x),2)+pow([2]/x,2))",1,3000);
 
-   TGraph *gr11 = new TGraph((numBin),x,y);
+   TGraph *gr11 = new TGraph((numBin-1),x,y);
    gr11->Draw("AP");
    gr11->GetXaxis()->SetTitle("E_{T}^{MC} [GeV]");
@@ -96,2 +78,84 @@
 // delete Gauss;
 }
+
+
+void ETmisResol(TTree *Analyze)
+{
+
+   TCanvas *c4 = new TCanvas("c4","ETmis resol",100,100,600,450);
+   c4->cd(); 
+  
+   const Int_t numBin=10;
+   double bins[numBin]={0,10,20,40,80,160,320,640,1080,2160};
+   TProfile *ETSoverET = new TProfile("ETSoverET","ETmis resol",(numBin-1),bins,-10,10);
+   Analyze->Draw("(ETmisResol.EtSmeare-ETmisResol.Et):(ETmisResol.Et)>>ETSoverET","(ETmisResol.Et)>1");
+
+   double rms[numBin];
+   double mean[numBin];
+
+   TCanvas *c5 = new TCanvas("c5","gerrors2",0,0,1000,650);
+   c5->cd(); int frame=0;
+   c5->Divide(6,2);
+
+   float x[numBin-1];
+   float y[numBin-1];
+
+   for ( int i=0; i<(numBin-1); i++)  // premiÃšre bin : i ==1 et pas i == 0
+      {
+        TAxis *xaxis = ETSoverET->GetXaxis();
+        float binCenter = xaxis->GetBinCenter(i+1);
+        int binMin = (int)xaxis->GetBinLowEdge(i+1);
+        int binMax = (int)xaxis->GetBinUpEdge(i+1);
+        cout<<"binMin "<<binMin<<" binMax "<<binMax<<" bin center "<<binCenter<<endl;
+        char tempMin[500];
+        sprintf(tempMin,"(ETmisResol.Et)>%d",binMin);
+        string mystringMin(tempMin);
+        char tempMax[500];
+        sprintf(tempMax,"(ETmisResol.Et)<%d",binMax);
+        string mystringMax(tempMax);
+        char tempName[500];
+        sprintf(tempName,"(ETmisResol.EtSmeare-ETmisResol.Et)>>hETSoverET%d",i);
+        string mystringName(tempName);
+
+        c5->cd(++frame);
+        GaussValues(Analyze,tempName,rms[i],mean[i],mystringMin,mystringMax);
+        x[i]=binCenter;
+        mean[i]=ETSoverET->GetBinContent(i+1);
+        y[i]=rms[i];
+        cout<<" mean val "<<mean[i]<<" rms value "<<rms[i]<<endl;
+
+      }
+
+   TCanvas *c6 = new TCanvas("c6","ETmis resolution",100,100,600,450);
+   c6->cd();
+
+   TF1 *fitfun = new TF1("user","sqrt(pow([0],2)+pow([1]/sqrt(x),2)+pow([2]/x,2))",1,3000);
+
+   TGraph *gr11 = new TGraph((numBin-1),x,y);
+   gr11->Draw("AP");
+   gr11->GetXaxis()->SetTitle("E_{T}^{MC} [GeV]");
+   gr11->GetYaxis()->SetTitle("#sigma(E_{T}^{rec}>");
+
+   //fitfun->FixParameter(0,0.0541930);
+   //fitfun->FixParameter(1,0.634908);
+   gr11->Fit("user","RQ");
+   gr11->Fit("user","RQ");
+   gr11->Fit("user","RQ");
+   gr11->Fit("user","R");
+
+
+
+ delete fitfun;
+// delete Gauss;
+}
+
+void General()
+{
+  setTDRStyle();
+  gROOT->Reset();
+
+  TFile *f1 = new TFile("Smearing.root","read");
+  TTree *Analyze = (TTree*)f1->Get("Analysis");
+  JetResol(Analyze);
+
+}
