#define data_cxx
#include "data.hpp"
#include "TF1.h"
#include "TH1.h"
#include "TH2.h"
#include "TStyle.h"
#include "TCanvas.h"
#include "TPostScript.h"
#include "TRandom.h"
#include "TMath.h"
#include "TLorentzVector.h"
#include <iostream>
#include <map>
#include <vector>
#include <algorithm>

#include "lm_access/BadLBNList.hpp"
//#include "BadLBNList.hpp"
#include "TauNNoutput.h"
#include "AnglesUtil.hpp"


float data::MC_TrkMatchEffCorr (float eta)
{
  float tmp=0.;

  TF1 *CC_corr = new TF1("CC_corr","pol1",-1.2,1.2);
  CC_corr->SetParameters(8.94579e-01,1.13404e-02);

  TF1 *ECpos_corr = new TF1("ECpos_corr","pol2",1.4,2.6);
  ECpos_corr->SetParameters(3.20992e+00,-2.18638e+00,4.56033e-01);
  
  TF1 *ECneg_corr = new TF1("ECneg_corr","pol5",-1.2,1.2);
  ECneg_corr->SetParameters(7.88176e+01,2.34414e+02,2.75933e+02,1.58811e+02,4.47067e+01,4.93190e+00);

  if( fabs(eta)<1.1 ) 
    tmp=0.927; //CC_corr->Eval(eta);
  if( eta>1.5 && eta<2.5 )
    tmp=ECpos_corr->Eval(eta);
  if( eta>-2.5 && eta<-1.5 )
    tmp=ECneg_corr->Eval(eta);

  return tmp;
}

void data::SetQCD (const int algo_int)
{
  is_qcd=algo_int;
  cout << "is_qcd: " << algo_int << endl;
}

void data::SetEventName (const char *fname)
{
  strcpy (event_name, fname);
  printf("event_name = %s \n",event_name);
}

/*
void data::Smear(TLorentzVector &v1, bool is_fiducial)
{
  float fiducial_width=0.03;
  float nonfiducial_width=0.06;
  float smear,scale;
  TRandom *rand1 = new TRandom(0);

  if (is_fiducial) {
    smear=rand1->Gaus(0.,fiducial_width);
    scale=1.;    //scaling: 1.00 in fiducial region
  } else {
    smear=rand1->Gaus(0.,nonfiducial_width);
    scale=1.;   //scaling: 0.94 in non-fiducial region
  }

  v1*=scale;
  v1*=(1+smear);
}


bool data::SortPT (const TLorentzVector& p1, const TLorentzVector& p2) 
{
  return p1.Px()<p2.Px();
};
*/

void data::Loop()
{
  //   In a ROOT session, you can do:
  //      Root > .L Tau.C
  //      Root > Tau t
  //      Root > t.GetEntry(12); // Fill t data members with entry number 12
  //      Root > t.Show();       // Show values of entry 12
  //      Root > t.Show(16);     // Read and show values of entry 16
  //      Root > t.Loop();       // Loop on all entries
  //

  //     This is the loop skeleton
  //       To read only selected branches, Insert statements like:
  // METHOD1:
  //    fChain->SetBranchStatus("*",0);  // disable all branches
  //    fChain->SetBranchStatus("branchname",1);  // activate branchname
  // METHOD2: replace line
  //    fChain->GetEntry(i);  // read all branches
  //by  b_branchname->GetEntry(i); //read only this branch
  if (fChain == 0) return;

  Int_t nentries = Int_t(fChain->GetEntriesFast());
  cout<<"Number of entries:  "<<nentries<<endl;
  cout<<"Number of entries:  "<<Int_t(fChain->GetEntries())<<endl;
  int doubleev = 0;
  int badlumi11 = 0;
  //int badlumi12 = 0;

  map<TString, int> StoreEvents;

  lm_access::BadLBNList* pbadlbnlist_v11 = new lm_access::BadLBNList();
  pbadlbnlist_v11->loadBadLBNs("./recorded_filelist_cs_em1trk_EM_HI_0-0.badlbns");

  float pi=TMath::Pi();

  int event_cnt=0;
  
  // random numbers for smearing
  TRandom *rand1 = new TRandom(0);

  // inital neural network stuff
  //TauNNoutput *NN0 = new TauNNoutput("NN_type0.dat");
  TauNNoutput *NN1 = new TauNNoutput("NN_type1_p13.08_2.dat");
  TauNNoutput *NN2 = new TauNNoutput("NN_type2_p13.08_2.dat");
  //TauNNoutput *NN3 = new TauNNoutput("NN_type3_p13.08_2.dat");

  //vector<TLorentzVector> em;

  //TFile *f = new TFile("outfile.root","RECREATE");

  // root tree
  TTree *_ttreeGood;
  Int_t _el_n;  
  Float_t _el_pt[3];  
  Float_t _el_px[3];
  Float_t _el_py[3];
  Float_t _el_pz[3];
  Float_t _el_E[3];
  Float_t _el_opt[3];  
  Float_t _el_opx[3];
  Float_t _el_opy[3];
  Float_t _el_opz[3];
  Float_t _el_oE[3];
  Float_t _el_emfr[3];
  Float_t _el_iso[3];
  Float_t _el_hmx8[3];
  Float_t _el_eta[3];
  Float_t _el_deta[3];
  Float_t _el_phi[3];
  Float_t _el_dphi[3];
  //Float_t _el_dphmet[3];
  Float_t _el_theta[3];  
  Float_t _el_fiducial[3];  
  Float_t _el_fiducial_eta[3];  
  Float_t _el_fiducial_phi[3];
  Float_t _el_calEM1_E[3];  
  Float_t _el_calEM2_E[3];  
  Float_t _el_calEM3_E[3];  
  Float_t _el_calEM4_E[3];  
  Float_t _el_calPS_E[3];  
  Float_t _el_calFH1_E[3];  
  Float_t _el_cellE[3];  
  Float_t _el_cellEt[3];
  Int_t _el_hasTM[3];
  Float_t _el_TMchi2[3];
  Float_t _el_TMspatialchi2[3];
  Float_t _el_recow[3];
  Float_t _el_tmw[3];
  Float_t _el_trkqpt[3];
  Float_t _el_trkphi[3];
  Float_t _el_trktheta[3];
  Float_t _el_trkz0[3];
  Float_t _el_trkr0[3];
  Float_t _inv_mass;
  Int_t _el_l3em[3];
  Int_t _el_emhi[3];
  Int_t _el_emhish[3];
  Int_t _el_emhishtr[3];
  //---------------------
  Float_t _ev_vtx[3];
  Int_t   _ev_run;
  Int_t   _ev_evt;
  Int_t   _ev_goodrun;
  Int_t   _ev_doubleev;
  Int_t   _ev_goodLBN;
  Int_t   _ev_dataMC;
  //-----------------------
  Int_t   _trg_emhi;
  Int_t   _trg_emhish;
  Int_t   _trg_emhitr;
  Int_t   _trg_emhishtr;
  Int_t   _trg_emmx;
  Int_t   _trg_emmxsh;
  Int_t   _trg_emmxtr;
  Int_t   _trg_emmxshtr;
  //-----------------------
  Int_t _l3e_n;
  Float_t         _l3e_et[40];  
  Float_t         _l3e_deteta[40]; 
  Float_t         _l3e_phi[40];  
  Float_t         _l3e_emfrac[40]; 
  Float_t         _l3e_chi2[40];  
  Float_t         _l3e_isolation[40]; 
  Float_t         _l3e_eoverp[40];   
  Float_t         _l3e_EM1width[40]; 
  Float_t         _l3e_EM2width[40]; 
  Float_t         _l3e_EM3width[40];  
  Float_t         _l3e_EM1Rwidth[40]; 
  Float_t         _l3e_EM2Rwidth[40]; 
  Float_t         _l3e_EM3Rwidth[40]; 
  Int_t           _l3e_L3Tool[40];   
  Int_t           _l3e_PSmatch[40];  
  Int_t           _l3e_TRmatch[40];  
  Int_t           _l3e_CALTRmatch[40];
  Int_t           _l3e_CPSTRmatch[40];
  //-----------------------
  Int_t           _l3trk_n;
  Float_t         _l3trk_eta[40];  
  Float_t         _l3trk_tanl[40]; 
  Float_t         _l3trk_phi[40];  
  Float_t         _l3trk_pt[40];   
  //-----------------------

  
  _ttreeGood  = new TTree("treeGood","preselected objects"); 
  _ttreeGood->Branch("el_n",      &_el_n,      "el_n/I");
  _ttreeGood->Branch("el_pt",      &_el_pt,      "el_pt[3]/F");
  _ttreeGood->Branch("el_px",      &_el_px,      "el_px[3]/F");
  _ttreeGood->Branch("el_py",      &_el_py,      "el_py[3]/F");
  _ttreeGood->Branch("el_pz",      &_el_pz,      "el_pz[3]/F");
  _ttreeGood->Branch("el_E",       &_el_E,       "el_E[3]/F");
  _ttreeGood->Branch("el_opx",     &_el_opx,     "el_opx[3]/F");
  _ttreeGood->Branch("el_opy",     &_el_opy,     "el_opy[3]/F");
  _ttreeGood->Branch("el_opz",     &_el_opz,     "el_opz[3]/F");
  _ttreeGood->Branch("el_opt",     &_el_opt,     "el_opt[3]/F");
  _ttreeGood->Branch("el_oE",     &_el_oE,     "el_oE[3]/F");
  _ttreeGood->Branch("el_emfr",    &_el_emfr,    "el_emfr[3]/F");
  _ttreeGood->Branch("el_iso",     &_el_iso,     "el_iso[3]/F");
  _ttreeGood->Branch("el_hmx8",    &_el_hmx8,    "el_hmx8[3]/F");
  _ttreeGood->Branch("el_eta",     &_el_eta,     "el_eta[3]/F");
  _ttreeGood->Branch("el_deta",    &_el_deta,    "el_deta[3]/F");
  _ttreeGood->Branch("el_phi",     &_el_phi,     "el_phi[3]/F");
  _ttreeGood->Branch("el_dphi",     &_el_dphi,     "el_dphi[3]/F");
  //_ttreeGood->Branch("el_dphmet",  &_el_dphmet,  "el_dphmet[3]/F");
  _ttreeGood->Branch("el_theta",     &_el_theta,     "el_theta[3]/F");
  _ttreeGood->Branch("el_fiducial",     &_el_fiducial,     "el_fiducial[3]/F");
  _ttreeGood->Branch("el_fiducial_eta",     &_el_fiducial_eta,     "el_fiducial_eta[3]/F");
  _ttreeGood->Branch("el_fiducial_phi",     &_el_fiducial_phi,     "el_fiducial_phi[3]/F");
  _ttreeGood->Branch("el_calEM1_E",     &_el_calEM1_E,     "el_calEM1_E[3]/F");
  _ttreeGood->Branch("el_calEM2_E",     &_el_calEM2_E,     "el_calEM2_E[3]/F");
  _ttreeGood->Branch("el_calEM3_E",     &_el_calEM3_E,     "el_calEM3_E[3]/F");
  _ttreeGood->Branch("el_calEM4_E",     &_el_calEM4_E,     "el_calEM4_E[3]/F");
  _ttreeGood->Branch("el_calPS_E",     &_el_calPS_E,     "el_calPS_E[3]/F");
  _ttreeGood->Branch("el_calFH1_E",     &_el_calFH1_E,     "el_calFH1_E[3]/F");
  _ttreeGood->Branch("el_cellE",     &_el_cellE,     "el_cellE[3]/F");
  _ttreeGood->Branch("el_cellEt",     &_el_cellEt,     "el_cellEt[3]/F");
  _ttreeGood->Branch("el_hasTM",     &_el_hasTM,     "el_hasTM[3]/I");
  _ttreeGood->Branch("el_TMchi2",     &_el_TMchi2,     "el_TMchi2[3]/F");
  _ttreeGood->Branch("el_TMspatialchi2",     &_el_TMspatialchi2,     "el_TMspatialchi2[3]/F");
  _ttreeGood->Branch("el_recow", &_el_recow, "el_recow[3]/F");
  _ttreeGood->Branch("el_tmw", &_el_tmw, "el_tmw[3]/F");
  _ttreeGood->Branch("el_trkqpt", &_el_trkqpt, "el_trkqpt[3]/F");
  _ttreeGood->Branch("el_trkphi", &_el_trkphi, "el_trkphi[3]/F");
  _ttreeGood->Branch("el_trktheta", &_el_trktheta, "el_trktheta[3]/F");
  _ttreeGood->Branch("el_trkz0", &_el_trkz0, "el_trkz0[3]/F");
  _ttreeGood->Branch("el_trkr0", &_el_trkr0, "el_trkr0[3]/F");
  _ttreeGood->Branch("inv_mass", &_inv_mass, "inv_mass/F");
  _ttreeGood->Branch("el_l3em", &_el_l3em, "el_l3em[3]/I");
  _ttreeGood->Branch("el_emhi", &_el_emhi, "el_emhi[3]/I");
  _ttreeGood->Branch("el_emhish", &_el_emhish, "el_emhish[3]/I");
  _ttreeGood->Branch("el_emhishtr", &_el_emhishtr, "el_emhishtr[3]/I");
  //----------------------------------------------------------------------
  _ttreeGood->Branch("ev_vtx",     &_ev_vtx,     "ev_vtx[3]/F");
  _ttreeGood->Branch("ev_run"  ,   &_ev_run,     "ev_run/I");
  _ttreeGood->Branch("ev_evt",     &_ev_evt,     "ev_evt/I");
  _ttreeGood->Branch("ev_goodrun",     &_ev_goodrun,     "ev_goodrun/I");
  _ttreeGood->Branch("ev_doubleev",     &_ev_doubleev,     "ev_doubleev/I");
  _ttreeGood->Branch("ev_goodLBN",     &_ev_goodLBN,     "ev_goodLBN/I");
  _ttreeGood->Branch("ev_dataMC",     &_ev_dataMC,     "ev_dataMC/I");
  //----------------------------------------------------------------------
  _ttreeGood->Branch("trg_emhi",    &_trg_emhi,    "trg_emhi/I");
  _ttreeGood->Branch("trg_emhish",    &_trg_emhish,    "trg_emhish/I");
  _ttreeGood->Branch("trg_emhitr",    &_trg_emhitr,    "trg_emhitr/I");
  _ttreeGood->Branch("trg_emhishtr",    &_trg_emhishtr,    "trg_emhishtr/I"); 
  _ttreeGood->Branch("trg_emmx",    &_trg_emmx,    "trg_emmx/I");
  _ttreeGood->Branch("trg_emmxsh",    &_trg_emmxsh,    "trg_emmxsh/I");
  _ttreeGood->Branch("trg_emmxtr",    &_trg_emmxtr,    "trg_emmxtr/I");
  _ttreeGood->Branch("trg_emmxshtr",    &_trg_emmxshtr,    "trg_emmxshtr/I");
  //----------------------------------------------------------------------
  _ttreeGood->Branch("l3e_n",    &_l3e_n,    "l3e_n/I");
  _ttreeGood->Branch("l3e_et",    &_l3e_et,    "l3e_et[40]/F");
  _ttreeGood->Branch("l3e_deteta",    &_l3e_deteta,    "l3e_deteta[40]/F");
  _ttreeGood->Branch("l3e_phi",    &_l3e_phi,    "l3e_phi[40]/F");
  _ttreeGood->Branch("l3e_emfrac",    &_l3e_emfrac,    "l3e_emfrac[40]/F");
  _ttreeGood->Branch("l3e_chi2",    &_l3e_chi2,    "l3e_chi2[40]/F");
  _ttreeGood->Branch("l3e_isolation",    &_l3e_isolation,    "l3e_isolation[40]/F");
  _ttreeGood->Branch("l3e_eoverp",    &_l3e_eoverp,    "l3e_eoverp[40]/F");
  _ttreeGood->Branch("l3e_EM1width",    &_l3e_EM1width,    "l3e_EM1width[40]/F");
  _ttreeGood->Branch("l3e_EM2width",    &_l3e_EM2width,    "l3e_EM2width[40]/F");
  _ttreeGood->Branch("l3e_EM3width",    &_l3e_EM3width,    "l3e_EM3width[40]/F");
  _ttreeGood->Branch("l3e_EM1Rwidth",    &_l3e_EM1Rwidth,    "l3e_EM1Rwidth[40]/F");
  _ttreeGood->Branch("l3e_EM2Rwidth",    &_l3e_EM2Rwidth,    "l3e_EM2Rwidth[40]/F");
  _ttreeGood->Branch("l3e_EM3Rwidth",    &_l3e_EM3Rwidth,    "l3e_EM3Rwidth[40]/F");
  _ttreeGood->Branch("l3e_L3Tool",    &_l3e_L3Tool,    "l3e_L3Tool[40]/F");
  _ttreeGood->Branch("l3e_PSmatch",    &_l3e_PSmatch,    "l3e_PSmatch[40]/F");
  _ttreeGood->Branch("l3e_TRmatch",    &_l3e_TRmatch,    "l3e_TRmatch[40]/F");
  _ttreeGood->Branch("l3e_CALTRmatch",    &_l3e_CALTRmatch,    "l3e_CALTRmatch[40]/F");
  _ttreeGood->Branch("l3e_CPSTRmatch",    &_l3e_CPSTRmatch,    "l3e_CPSTRmatch[40]/F");
  //----------------------------------------------------------------------
  _ttreeGood->Branch("l3trk_n",    &_l3trk_n,    "l3trk_n/I");
  _ttreeGood->Branch("l3trk_pt",    &_l3trk_pt,    "l3trk_pt[40]/F");
  _ttreeGood->Branch("l3trk_eta",    &_l3trk_eta,    "l3trk_eta[40]/F");
  _ttreeGood->Branch("l3trk_tanl",    &_l3trk_tanl,    "l3trk_tanl[40]/F");
  _ttreeGood->Branch("l3trk_phi",    &_l3trk_phi,    "l3trk_phi[40]/F");
  //----------------------------------------------------------------------

  Int_t nbytes = 0, nb = 0;
  //for (Int_t jentry=0; jentry<nentries;jentry++) {
    for (Int_t jentry=0; jentry<100;jentry++) {
    Int_t ientry = LoadTree(jentry); //in case of a TChain, ientry is the entry number in the current file
    if (ientry < 0) break;
    nb = fChain->GetEntry(jentry);   nbytes += nb;
    
    if (fmod(jentry,10000)==0)
      cout << "event counter: " << jentry << endl;

    bool debug = false;

    if(evt_flag==1 || evt_flag==3) {
      ++event_cnt;
      
      /*
      cout << "+++++++++++++++++++++++++++++++++++++++++++++++++"<<endl;
      for (int k=0; k< LETM_Nelectron;++k) {
	cout << LETM_pT[k]<<":"<<LETM_caldeteta[k]<<endl;
      }
      for (int k=0; k< AE_Nelectron;++k) {
	if(AE_isolation[k]<0.15 
	   && AE_emfrac[k]>0.9
	   && AE_Hmx8[k]<20.
	   && (abs(AE_id[k])==10 || abs(AE_id[k])==11)) {
	  cout << AE_pT[k]<<":"<<AE_caldeteta[k]<<":"<<AE_hasTM[k]<<endl;
	}
      }
      */

      //////////////////////////
      // REMOVE DOUBLE EVENTS //
      //////////////////////////
      bool is_doubleevent = false;
      TString eventstr(Form("%i_%i",RunNumber,EventNumber));
      if(StoreEvents.find(eventstr) != StoreEvents.end()) {
	cout << " Found double events !!!!! " << RunNumber << " " << EventNumber << endl;
	is_doubleevent = true;
	doubleev++;
      }
      else {
	++StoreEvents[eventstr];
      }
    
    /////////////////////
    // CHECK GOOD RUNS //
    /////////////////////
    bool goodrun = true;
    if(Data_MC==0)
      goodrun = is_good_cal(RunNumber);


    ////////////////////////
    // CHECK FOR TRIGGERS //
    ////////////////////////
    //bool trigger_fired = false;
    //if (TR_EM_HI || TR_EM_HI_SH_TR || TR_EM_HI_SH ) {
    //trigger_fired = true;
    bool good_LBN=true;
    if(Data_MC==0 && pbadlbnlist_v11->isBadLBN(LumiBlock)) {
      if(debug) cout << "WARNING: Bad Lumi Block " << LumiBlock << endl;
      good_LBN=false;
      //trigger_fired = false;
      badlumi11++;
    }
    //}

    vector<int> goodEM;
    vector<int> goodEM_notm;
    vector<int> iEMs;
     
    goodEM.clear();  
    goodEM_notm.clear();  
    iEMs.clear();  
    
    for (int l=0;l<AE_Nelectron; ++l) {
      if(AE_isolation[l]<0.15 
	 && AE_emfrac[l]>0.9
	 && AE_pT[l]>8.
	 && AE_Hmx8[l]<20.
	 && (AE_id[l]==10 || abs(AE_id[l])==11)) {
	if (AE_hasTM[l]) {
	  goodEM.push_back(l);
	} else {
	  goodEM_notm.push_back(l);
	}
      }
    }  //for loop

    //fill up goodEM vector -> sorting in tm, pt
    if (goodEM_notm.size()>0 && goodEM.size()<3) {
      for (int k=0;k<goodEM_notm.size();++k) {
	if (goodEM.size()<3) {
	  goodEM.push_back(goodEM_notm[k]);
	}
      }
    }


    //if (goodEM.size()>=2 && goodrun && !is_doubleevent && good_LBN) {
    if (goodEM.size()>=2 && goodrun && !is_doubleevent) {
      
      _el_n=goodEM.size();
      _ev_run=RunNumber;
      _ev_evt=EventNumber;
      _ev_vtx[0]=PVtx_x;
      _ev_vtx[1]=PVtx_y;
      _ev_vtx[2]=PVtx_z;
      _ev_goodrun=goodrun;
      _ev_doubleev=is_doubleevent;
      _ev_goodLBN=good_LBN;
      _ev_dataMC=Data_MC;

      _trg_emhi=TR_EM_HI;
      _trg_emhish=TR_EM_HI_SH;
      _trg_emhitr=TR_EM_HI_TR;
      _trg_emhishtr=TR_EM_HI_SH_TR;
      _trg_emmx=TR_EM_MX;
      _trg_emmxsh=TR_EM_MX_SH;
      _trg_emmxtr=TR_EM_MX_TR;
      _trg_emmxshtr=TR_EM_MX_SH_TR;

      //////////////
      // Smearing //
      //////////////
      float smear_em1=0.;
      float em1_spx=0.;
      float em1_spy=0.;
      float em1_spz=0.;
      float em1_sE=0.;
      float em1_spt=0.;

      float smear_em2=0.;
      float em2_spx=0.;
      float em2_spy=0.;
      float em2_spz=0.;
      float em2_sE=0.;
      float em2_spt=0.;

      /*
      if(Data_MC>0) {
	if (AE_fiducial[goodEM[0]]) {
	  //19-Mar-04: smearing in fiducial: 0.035
	  //           smearing in non-fiducial: 0.07
	  smear_em1=rand1->Gaus(0.,0.03);
	} else {
	  //smear_em1=rand1->Gaus(0.,0.075);
	  smear_em1=rand1->Gaus(0.,0.06);
	}
	if (AE_fiducial[goodEM[1]]) {
	  smear_em2=rand1->Gaus(0.,0.03);
	} else {
	  smear_em2=rand1->Gaus(0.,0.06);
	}
      } //if(Data_MC>0)
      */

      em1_spx = AE_px[goodEM[0]]*(1+smear_em1);
      em1_spy = AE_py[goodEM[0]]*(1+smear_em1);
      em1_spz = AE_pz[goodEM[0]]*(1+smear_em1);
      em1_sE  = AE_E[goodEM[0]]*(1+smear_em1);
      em1_spt = sqrt(pow(em1_spx,2)+pow(em1_spy,2));

      em2_spx = AE_px[goodEM[1]]*(1+smear_em2);
      em2_spy = AE_py[goodEM[1]]*(1+smear_em2);
      em2_spz = AE_pz[goodEM[1]]*(1+smear_em2);
      em2_sE  = AE_E[goodEM[1]]*(1+smear_em2);
      em2_spt = sqrt(pow(em2_spx,2)+pow(em2_spy,2));

      if(em1_spt>em2_spt) {
	iEMs.push_back(goodEM[0]);
	iEMs.push_back(goodEM[1]);
	_el_px[0]=em1_spx;
	_el_py[0]=em1_spy;
	_el_pz[0]=em1_spz;
	_el_E[0]=em1_sE;
	_el_pt[0]=em1_spt;
	_el_px[1]=em2_spx;
	_el_py[1]=em2_spy;
	_el_pz[1]=em2_spz;
	_el_E[1]=em2_sE;
	_el_pt[1]=em2_spt;
      } else {
	iEMs.push_back(goodEM[1]);
	iEMs.push_back(goodEM[0]);
	_el_px[0]=em2_spx;
	_el_py[0]=em2_spy;
	_el_pz[0]=em2_spz;
	_el_E[0]=em2_sE;
	_el_pt[0]=em2_spt;
	_el_px[1]=em1_spx;
	_el_py[1]=em1_spy;
	_el_pz[1]=em1_spz;
	_el_E[1]=em1_sE;
	_el_pt[1]=em1_spt;
      }
    
      for (int k=0; k<iEMs.size();++k) {
	//cout<<"filled root tuple with:" <<endl;
	//cout << AE_pT[iEMs[k]] << ":"<<AE_caldeteta[iEMs[k]]<<":"<< AE_hasTM[iEMs[k]]<<endl;

	_el_oE[k]=AE_E[iEMs[k]];
	_el_opt[k]=AE_pT[iEMs[k]];
	_el_opx[k]=AE_px[iEMs[k]];
	_el_opy[k]=AE_py[iEMs[k]];
	_el_opz[k]=AE_pz[iEMs[k]];
	_el_emfr[k]=AE_emfrac[iEMs[k]];
	_el_iso[k]=AE_isolation[iEMs[k]];
	_el_hmx8[k]=AE_Hmx8[iEMs[k]];
	_el_eta[k]=AE_eta[iEMs[k]];
	_el_deta[k]=AE_caldeteta[iEMs[k]];
	_el_phi[k]=AE_phi[iEMs[k]];
	_el_dphi[k]=AE_caldetphi[iEMs[k]];
	_el_theta[k]=AE_theta[iEMs[k]];
	_el_fiducial[k]=AE_fiducial[iEMs[k]];
	_el_fiducial_eta[k]=AE_fiducial_eta[iEMs[k]];
	_el_fiducial_phi[k]=AE_fiducial_phi[iEMs[k]];
	_el_calEM1_E[k]=AE_calEM1_E[iEMs[k]];
	_el_calEM2_E[k]=AE_calEM2_E[iEMs[k]];
	_el_calEM3_E[k]=AE_calEM3_E[iEMs[k]];
	_el_calEM4_E[k]=AE_calEM4_E[iEMs[k]];
	_el_calPS_E[k]=AE_calPS_E[iEMs[k]];
	_el_calFH1_E[k]=AE_calFH1_E[iEMs[k]];
	_el_cellE[k]=AE_CellE[iEMs[k]];
	_el_cellEt[k]=AE_CellEt[iEMs[k]];
	_el_hasTM[k]=AE_hasTM[iEMs[k]];
	_el_TMchi2[k]=AE_TMchi2[iEMs[k]];
	_el_TMspatialchi2[k]=AE_TMspatialchi2[iEMs[k]];


	_el_l3em[k]=-1;
	_el_emhi[k]=0;
	_el_emhish[k]=0;
	_el_emhishtr[k]=0;

	//identify L3 track:
	bool track_found=false;
	if (Data_MC==0) {
	  for(int t=0;t<L3Trk_Ntrack;++t) {
	    if(L3Trk_pT[t]>=12.)
	      track_found=true;
	  }
	  
	  //identify L3 electron
	  for (int l=0;l<L3E_Nelectron; ++l) {
	    float dR=sqrt(pow(L3E_deteta[l]-AE_caldeteta[iEMs[k]],2)+pow(kinem::delta_phi(L3E_phi[l],AE_phi[iEMs[k]]),2));
	    if (dR<0.5) {
	      _el_l3em[k]=l;
	      cout << "found match between probe em and L3 em" << endl;
	      cout << "l3 et,eta,phi,L3: " << L3E_Et[l] <<","<< L3E_deteta[l] <<","<< L3E_phi[l] <<","<< L3E_L3Tool[l] << endl;
	      if(L3E_L3Tool[l]==1 && L3E_Et[l]>30.) _el_emhi[k]=1;
	      if(L3E_L3Tool[l]==2 && L3E_Et[l]>20.) _el_emhish[k]=1;
	      if(L3E_L3Tool[l]==2 && L3E_Et[l]>12. && track_found) _el_emhishtr[k]=1;
	      cout << "em_hi/em_hi_sh/em_hi_sh_tr: " << _el_emhi[k] <<"/"<< _el_emhish[k] <<"/"<< _el_emhishtr[k] << endl;
	    }
	  }
	} //if(Data_MC==0)

	// Set event weight 
	_el_recow[k]=1.;
	_el_tmw[k]=1.;
	/*
	  if (Data_MC>0) {
	  //p14 reconstructions correction (from Marc):
	  if (fabs(AE_caldeteta[iEMs[k]])<1.3)  _el_recow[k]=0.87;//weight=0.96;
	  if (fabs(AE_caldeteta[iEMs[k]])>=1.3 && fabs(AE_caldeteta[iEMs[k]])<2.5)  _el_recow[k]=0.94;//weight=0.94;
	  
	  //p14 track match efficiency correction:
	  if (fabs(AE_caldeteta[iEMs[k]])<2.5 && AE_hasTM[iEMs[k]]) 
	  _el_tmw[k] = MC_TrkMatchEffCorr(AE_caldeteta[iEMs[k]]);
	    } //set event weight
	*/
      } //for loop
	
      //diem mass
      _inv_mass = sqrt(pow(em1_sE+em2_sE,2)
		      - pow(em1_spx+em2_spx,2)
		      - pow(em1_spy+em2_spy,2)
		      - pow(em1_spz+em2_spz,2));
      
      //find track variables if electron has track match
      int index_em1=-1;
      int index_em2=-1;
      for (int j=0;j<LETM_Nelectron;++j) {
	if (AE_hasTM[iEMs[0]] && fabs(AE_pT[iEMs[0]]-LETM_pT[j])<0.01) {
	  index_em1=j;
	}
	if (AE_hasTM[iEMs[1]] && fabs(AE_pT[iEMs[1]]-LETM_pT[j])<0.01) {
	  index_em2=j;
	}
      }
      
      if (index_em1>-1) {
	_el_trkqpt[0]=LETM_trkqpt[index_em1];
	_el_trkphi[0]=LETM_trkphi[index_em1];
	_el_trktheta[0]=LETM_trktheta[index_em1];
	_el_trkz0[0]=LETM_trkz0[index_em1];
	_el_trkr0[0]=LETM_trkr0[index_em1];
      } else {
	_el_trkqpt[0]=-100;
	_el_trkphi[0]=-100;
	_el_trktheta[0]=-100;
	_el_trkz0[0]=-100;
	_el_trkr0[0]=-100;
      }
      if (index_em2>-1) {
	_el_trkqpt[1]=LETM_trkqpt[index_em2];
	_el_trkphi[1]=LETM_trkphi[index_em2];
	_el_trktheta[1]=LETM_trktheta[index_em2];
	_el_trkz0[1]=LETM_trkz0[index_em2];
	_el_trkr0[1]=LETM_trkr0[index_em2];
      } else {
	_el_trkqpt[1]=-100;
	_el_trkphi[1]=-100;
	_el_trktheta[1]=-100;
	_el_trkz0[1]=-100;
	_el_trkr0[1]=-100;
      }


      //fill l3 variables:
      _l3e_n=L3E_Nelectron;
      _l3trk_n=L3Trk_Ntrack;
      for (int k=0;k<L3E_Nelectron;++k) {
	_l3e_et[k]=L3E_Et[k];
	_l3e_et[k]=L3E_deteta[k];
	_l3e_et[k]=L3E_phi[k];
	_l3e_et[k]=L3E_emfrac[k];
	_l3e_et[k]=L3E_chi2[k];
	_l3e_et[k]=L3E_isolation[k];
	_l3e_et[k]=L3E_eoverp[k];
	_l3e_et[k]=L3E_EM1width[k];
	_l3e_et[k]=L3E_EM2width[k];
	_l3e_et[k]=L3E_EM3width[k];
	_l3e_et[k]=L3E_EM1Rwidth[k];
	_l3e_et[k]=L3E_EM2Rwidth[k];
	_l3e_et[k]=L3E_EM3Rwidth[k];
	_l3e_et[k]=L3E_L3Tool[k];
	_l3e_et[k]=L3E_PSmatch[k];
	_l3e_et[k]=L3E_TRmatch[k];
	_l3e_et[k]=L3E_CALTRmatch[k];
	_l3e_et[k]=L3E_CPSTRmatch[k];
      }
      for (int k=0;k<L3Trk_Ntrack;++k){
	_l3trk_eta[k]=L3Trk_eta[k];
	_l3trk_tanl[k]=L3Trk_tanl[k];
	_l3trk_phi[k]=L3Trk_phi[k];
	_l3trk_pt[k]=L3Trk_pT[k];
      }
      //////////////////////
    
      _ttreeGood->Fill();
    } //if (goodEM.size()>=2)


    /*
    //determine trigger eff
    double halfpoint = 20.79;
    double slope = 1.99;
    double plateau = 1.0;
    double offset = 0;
    double arg1 = (em1_spt - halfpoint)/(TMath::Sqrt(em1_spt)*slope);
    double fitval_1 = offset+0.5*plateau*(1+TMath::Erf(arg1));
    double arg2 = (em2_spt - halfpoint)/(TMath::Sqrt(em2_spt)*slope);
    double fitval_2 = offset+0.5*plateau*(1+TMath::Erf(arg2));
    trig_prob=fitval_1+fitval_2-fitval_1*fitval_2;
    if(trig_prob>1.) {
    cout << "SET trig_prob TO 1" << endl;
    trig_prob=1.;
    }
    }
    */

    } //++event_cnt;
  } // loop over all events
  
  //TFile f("diem_tree.root","recreate");
  TFile f(event_name,"recreate");
  _ttreeGood->Write();

  cout << "Event counter             : " << event_cnt << endl;
  cout << "Found " << doubleev << " double events!!! " << endl;
  cout << "Bad Lumi Blocks 11 " << badlumi11 << endl;
}


bool data::is_good_cal(int run)
{
  bool goodcal =  true;
  int badrun = 0;
  if(run ==180949) badrun = 1;
  else if(run ==180678) badrun = 1;
  else if(run ==180677) badrun = 1;
  else if(run ==180676) badrun = 1;
  else if(run ==180675) badrun = 1;
  else if(run ==180667) badrun = 1;
  else if(run ==180625) badrun = 1;
  else if(run ==180407) badrun = 1;
  else if(run ==180268) badrun = 1;
  else if(run ==179908) badrun = 1;
  else if(run ==179907) badrun = 1;
  else if(run ==179906) badrun = 1;
  else if(run ==179905) badrun = 1;
  else if(run ==179904) badrun = 1;
  else if(run ==179903) badrun = 1;
  else if(run ==179902) badrun = 1;
  else if(run ==179901) badrun = 1;
  else if(run ==179900) badrun = 1;
  else if(run ==179899) badrun = 1;
  else if(run ==179898) badrun = 1;
  else if(run ==179897) badrun = 1;
  else if(run ==179738) badrun = 1;
  else if(run ==179470) badrun = 1;
  else if(run ==179468) badrun = 1;
  else if(run ==179192) badrun = 1;
  else if(run ==179054) badrun = 1;
  else if(run ==179053) badrun = 1;
  else if(run ==179051) badrun = 1;
  else if(run ==179048) badrun = 1;
  else if(run ==178734) badrun = 1;
  else if(run ==178645) badrun = 1;
  else if(run ==178619) badrun = 1;
  else if(run ==178618) badrun = 1;
  else if(run ==178225) badrun = 1;
  else if(run ==178224) badrun = 1;
  else if(run ==178221) badrun = 1;
  else if(run ==178220) badrun = 1;
  else if(run ==178150) badrun = 1;
  else if(run ==178093) badrun = 1;
  else if(run ==178091) badrun = 1;
  else if(run ==178070) badrun = 1;
  else if(run ==178069) badrun = 1;
  else if(run ==177740) badrun = 1;
  else if(run ==177640) badrun = 1;
  else if(run ==177638) badrun = 1;
  else if(run ==176896) badrun = 1;
  else if(run ==176821) badrun = 1;
  else if(run ==176818) badrun = 1;
  else if(run ==176711) badrun = 1;
  else if(run ==176591) badrun = 1;
  else if(run ==176590) badrun = 1;
  else if(run ==176589) badrun = 1;
  else if(run ==176588) badrun = 1;
  else if(run ==176587) badrun = 1;
  else if(run ==176585) badrun = 1;
  else if(run ==176571) badrun = 1;
  else if(run ==176570) badrun = 1;
  else if(run ==176569) badrun = 1;
  else if(run ==176214) badrun = 1;
  else if(run ==176213) badrun = 1;
  else if(run ==176202) badrun = 1;
  else if(run ==176200) badrun = 1;
  else if(run ==176199) badrun = 1;
  else if(run ==176198) badrun = 1;
  else if(run ==175665) badrun = 1;
  else if(run ==175651) badrun = 1;
  else if(run ==175650) badrun = 1;
  else if(run ==175649) badrun = 1;
  else if(run ==175647) badrun = 1;
  else if(run ==174846) badrun = 1;
  else if(run ==174845) badrun = 1;
  else if(run ==174752) badrun = 1;
  else if(run ==174589) badrun = 1;
  else if(run ==174588) badrun = 1;
  else if(run ==174587) badrun = 1;
  else if(run ==174586) badrun = 1;
  else if(run ==174511) badrun = 1;
  else if(run ==174510) badrun = 1;
  else if(run ==174509) badrun = 1;
  else if(run ==174508) badrun = 1;
  else if(run ==174507) badrun = 1;
  else if(run ==174506) badrun = 1;
  else if(run ==174505) badrun = 1;
  else if(run ==174504) badrun = 1;
  else if(run ==174503) badrun = 1;
  else if(run ==174502) badrun = 1;
  else if(run ==174338) badrun = 1;
  else if(run ==173974) badrun = 1;
  else if(run ==173968) badrun = 1;
  else if(run ==173966) badrun = 1;
  else if(run ==173965) badrun = 1;
  else if(run ==173923) badrun = 1;
  else if(run ==173922) badrun = 1;
  else if(run ==173921) badrun = 1;
  else if(run ==173920) badrun = 1;
  else if(run ==173919) badrun = 1;
  else if(run ==173483) badrun = 1;
  else if(run ==173352) badrun = 1;
  else if(run ==173007) badrun = 1;
  else if(run ==172850) badrun = 1;
  else if(run ==172811) badrun = 1;
  else if(run ==172810) badrun = 1;
  else if(run ==172708) badrun = 1;
  else if(run ==172526) badrun = 1;
  else if(run ==172359) badrun = 1;
  else if(run ==172217) badrun = 1;
  else if(run ==172179) badrun = 1;
  else if(run ==172178) badrun = 1;
  else if(run ==172177) badrun = 1;
  else if(run ==172176) badrun = 1;
  else if(run ==172175) badrun = 1;
  else if(run ==172174) badrun = 1;
  else if(run ==171899) badrun = 1;
  else if(run ==171885) badrun = 1;
  else if(run ==171829) badrun = 1;
  else if(run ==171818) badrun = 1;
  else if(run ==171755) badrun = 1;
  else if(run ==171659) badrun = 1;
  else if(run ==171634) badrun = 1;
  else if(run ==171611) badrun = 1;
  else if(run ==171605) badrun = 1;
  else if(run ==171604) badrun = 1;
  else if(run ==170175) badrun = 1;
  else if(run ==169929) badrun = 1;
  else if(run ==169867) badrun = 1;
  else if(run ==169861) badrun = 1;
  else if(run ==169860) badrun = 1;
  else if(run ==169512) badrun = 1;
  else if(run ==169511) badrun = 1;
  else if(run ==169510) badrun = 1;
  else if(run ==169055) badrun = 1;
  else if(run ==168870) badrun = 1;
  else if(run ==168514) badrun = 1;
  else if(run ==168414) badrun = 1;
  else if(run ==168388) badrun = 1;
  else if(run ==168033) badrun = 1;
  else if(run ==168032) badrun = 1;
  else if(run ==168031) badrun = 1;
  else if(run ==168030) badrun = 1;
  else if(run ==167393) badrun = 1;
  else if(run ==167390) badrun = 1;
  else if(run ==167387) badrun = 1;
  else if(run ==167262) badrun = 1;
  else if(run ==167257) badrun = 1;
  else if(run ==167100) badrun = 1;
  else if(run ==167099) badrun = 1;
  else if(run ==167098) badrun = 1;
  else if(run ==167094) badrun = 1;
  else if(run ==167092) badrun = 1;
  else if(run ==166788) badrun = 1;
  else if(run ==166785) badrun = 1;
  else if(run ==166784) badrun = 1;
  else if(run ==166315) badrun = 1;
  else if(run ==165060) badrun = 1;
  else if(run ==164504) badrun = 1;
  else if(run ==164503) badrun = 1;
  else if(run ==164502) badrun = 1;
  else if(run ==164500) badrun = 1;
  else if(run ==164497) badrun = 1;
  else if(run ==164496) badrun = 1;
  else if(run ==163972) badrun = 1;
  else if(run ==163307) badrun = 1;
  else if(run ==161972) badrun = 1;
  else if(run ==161928) badrun = 1;
  else if(run ==161920) badrun = 1;
  else if(run ==161916) badrun = 1;
  else if(run ==161914) badrun = 1;
  else if(run ==161912) badrun = 1;
  else if(run ==161878) badrun = 1;
  else if(run ==161877) badrun = 1;
  else if(run ==161875) badrun = 1;
  else if(run ==161874) badrun = 1;
  else if(run ==161873) badrun = 1;
  else if(run ==161872) badrun = 1;
  else if(run ==161871) badrun = 1;
  else if(run ==161865) badrun = 1;
  else if(run ==161864) badrun = 1;
  else if(run ==161863) badrun = 1;
  else if(run ==161664) badrun = 1;
  else if(run ==161663) badrun = 1;
  else if(run ==161662) badrun = 1;
  else if(run ==161661) badrun = 1;
  else if(run ==161660) badrun = 1;
  else if(run ==161659) badrun = 1;
  else if(run ==161658) badrun = 1;
  else if(run ==161657) badrun = 1;
  else if(run ==161656) badrun = 1;
  else if(run ==161654) badrun = 1;
  else if(run ==161648) badrun = 1;
  else if(run ==161647) badrun = 1;
  else if(run ==161646) badrun = 1;
  else if(run ==161645) badrun = 1;
  else if(run ==161644) badrun = 1;
  else if(run ==161643) badrun = 1;
  else if(run ==161630) badrun = 1;
  else if(run ==161629) badrun = 1;
  else if(run ==161628) badrun = 1;
  else if(run ==161627) badrun = 1;
  else if(run ==161624) badrun = 1;
  else if(run ==161623) badrun = 1;
  else if(run ==161622) badrun = 1;
  else if(run ==161621) badrun = 1;
  else if(run ==161620) badrun = 1;
  else if(run ==161619) badrun = 1;
  else if(run ==161618) badrun = 1;
  else if(run ==161617) badrun = 1;
  else if(run ==161616) badrun = 1;
  else if(run ==161311) badrun = 1;
  else if(run ==161310) badrun = 1;
  else if(run ==161309) badrun = 1;
  else if(run ==161308) badrun = 1;
  else if(run ==161307) badrun = 1;
  else if(run ==161303) badrun = 1;
  else if(run ==161301) badrun = 1;
  else if(run ==161300) badrun = 1;
  else if(run ==161299) badrun = 1;
  else if(run ==161297) badrun = 1;
  else if(run ==161296) badrun = 1;
  else if(run ==161294) badrun = 1;
  else if(run ==161293) badrun = 1;
  else if(run ==161202) badrun = 1;
  else if(run ==161201) badrun = 1;
  else if(run ==161200) badrun = 1;
  else if(run ==161199) badrun = 1;
  else if(run ==161198) badrun = 1;
  else if(run ==161197) badrun = 1;
  else if(run ==161196) badrun = 1;
  else if(run ==161195) badrun = 1;
  else if(run ==161194) badrun = 1;
  else if(run ==161193) badrun = 1;
  else if(run ==161192) badrun = 1;
  else if(run ==161191) badrun = 1;
  else if(run ==161189) badrun = 1;
  else if(run ==161188) badrun = 1;
  else if(run ==161187) badrun = 1;
  else if(run ==161186) badrun = 1;
  else if(run ==161110) badrun = 1;
  else if(run ==161106) badrun = 1;
  else if(run ==161101) badrun = 1;
  else if(run ==161099) badrun = 1;
  else if(run ==161098) badrun = 1;
  else if(run ==161083) badrun = 1;
  else if(run ==161053) badrun = 1;
  else if(run ==161049) badrun = 1;
  else if(run ==161048) badrun = 1;
  else if(run ==161041) badrun = 1;
  else if(run ==161031) badrun = 1;
  else if(run ==161030) badrun = 1;
  else if(run ==161029) badrun = 1;
  else if(run ==161028) badrun = 1;
  else if(run ==161026) badrun = 1;
  else if(run ==161025) badrun = 1;
  else if(run ==161023) badrun = 1;
  else if(run ==161020) badrun = 1;
  else if(run ==160704) badrun = 1;
  else if(run ==160700) badrun = 1;
  else if(run ==160699) badrun = 1;
  else if(run ==160694) badrun = 1;
  else if(run ==160693) badrun = 1;
  else if(run ==160692) badrun = 1;
  else if(run ==160690) badrun = 1;
  else if(run ==160689) badrun = 1;
  else if(run ==160688) badrun = 1;
  else if(run ==160686) badrun = 1;
  else if(run ==160588) badrun = 1;
  else if(run ==160586) badrun = 1;
  else if(run ==160585) badrun = 1;
  else if(run ==160582) badrun = 1;
  else if(run ==160554) badrun = 1;
  else if(run ==160553) badrun = 1;
  else if(run ==160552) badrun = 1;
  else if(run ==160551) badrun = 1;
  else if(run ==160550) badrun = 1;
  else if(run ==160549) badrun = 1;
  else if(run ==160548) badrun = 1;
  else if(run ==160546) badrun = 1;
  else if(run ==160545) badrun = 1;
  else if(run ==160540) badrun = 1;
  else if(run ==160539) badrun = 1;
  else if(run ==160538) badrun = 1;
  else if(run ==160472) badrun = 1;
  else if(run ==160471) badrun = 1;
  else if(run ==160470) badrun = 1;
  else if(run ==160459) badrun = 1;
  else if(run ==160449) badrun = 1;
  else if(run ==160433) badrun = 1;
  else if(run ==160432) badrun = 1;
  else if(run ==160431) badrun = 1;
  else if(run ==160430) badrun = 1;
  else if(run ==160429) badrun = 1;
  else if(run ==160428) badrun = 1;
  else if(run ==160426) badrun = 1;
  else if(run ==160425) badrun = 1;
  else if(run ==160424) badrun = 1;
  else if(run ==160423) badrun = 1;
  else if(run ==160422) badrun = 1;
  else if(run ==160421) badrun = 1;
  else if(run ==160420) badrun = 1;
  else if(run ==160419) badrun = 1;
  else if(run ==160362) badrun = 1;
  else if(run ==160361) badrun = 1;
  else if(run ==160360) badrun = 1;
  else if(run ==160229) badrun = 1;
  else if(run ==160228) badrun = 1;
  else if(run ==160227) badrun = 1;
  else if(run ==160226) badrun = 1;
  else if(run ==160225) badrun = 1;
  else if(run ==160224) badrun = 1;
  else if(run ==160223) badrun = 1;
  else if(run ==160221) badrun = 1;
  else if(run ==160220) badrun = 1;
  else if(run ==160219) badrun = 1;
  else if(run ==160218) badrun = 1;
  else if(run ==160217) badrun = 1;
  else if(run ==160216) badrun = 1;
  else if(run ==160215) badrun = 1;
  else if(run ==160214) badrun = 1;
  else if(run ==160210) badrun = 1;
  else if(run ==160209) badrun = 1;
  else if(run ==160208) badrun = 1;
  else if(run ==160207) badrun = 1;
  else if(run ==160206) badrun = 1;
  else if(run ==160205) badrun = 1;
  else if(run ==160198) badrun = 1;
  else if(run ==160196) badrun = 1;
  else if(run ==160145) badrun = 1;
  else if(run ==160144) badrun = 1;
  else if(run ==160143) badrun = 1;
  else if(run ==160142) badrun = 1;
  else if(run ==160141) badrun = 1;
  else if(run ==159776) badrun = 1;
  else if(run ==159775) badrun = 1;
  else if(run ==159774) badrun = 1;
  else if(run ==159773) badrun = 1;
  else if(run ==159771) badrun = 1;
  else if(run ==159770) badrun = 1;
  else if(run ==159769) badrun = 1;
  else if(run ==159768) badrun = 1;
  else if(run ==159767) badrun = 1;
  else if(run ==159766) badrun = 1;
  else if(run ==159765) badrun = 1;
  else if(run ==159764) badrun = 1;
  else if(run ==159763) badrun = 1;
  else if(run ==159707) badrun = 1;
  else if(run ==159706) badrun = 1;
  else if(run ==159705) badrun = 1;
  else if(run ==159704) badrun = 1;
  else if(run ==159703) badrun = 1;
  else if(run ==159702) badrun = 1;
  else if(run ==159701) badrun = 1;
  else if(run ==159700) badrun = 1;
  else if(run ==159699) badrun = 1;
  else if(run ==159698) badrun = 1;
  else if(run ==159453) badrun = 1;
  else if(run ==159452) badrun = 1;
  else if(run ==159451) badrun = 1;
  else if(run ==159449) badrun = 1;
  else if(run ==159448) badrun = 1;
  else if(run ==159447) badrun = 1;
  else if(run ==159446) badrun = 1;
  else if(run ==159445) badrun = 1;
  else if(run ==159443) badrun = 1;
  else if(run ==159442) badrun = 1;
  else if(run ==159441) badrun = 1;
  else if(run ==159044) badrun = 1;
  else if(run ==159040) badrun = 1;
  else if(run ==159033) badrun = 1;
  else if(run ==159030) badrun = 1;
  else if(run ==159029) badrun = 1;
  else if(run ==159028) badrun = 1;
  else if(run ==159027) badrun = 1;
  else if(run ==159002) badrun = 1;
  else if(run ==159001) badrun = 1;
  else if(run ==158999) badrun = 1;
  else if(run ==158998) badrun = 1;
  else if(run ==158997) badrun = 1;
  else if(run ==158996) badrun = 1;
  else if(run ==158995) badrun = 1;
  else if(run ==158994) badrun = 1;
  else if(run ==158993) badrun = 1;
  else if(run ==158983) badrun = 1;
  else if(run ==158982) badrun = 1;
  else if(run ==158981) badrun = 1;
  else if(run ==158979) badrun = 1;
  else if(run ==158978) badrun = 1;
  else if(run ==158977) badrun = 1;
  else if(run ==158976) badrun = 1;
  else if(run ==158975) badrun = 1;
  else if(run ==158974) badrun = 1;
  else if(run ==158973) badrun = 1;
  else if(run ==158972) badrun = 1;
  else if(run ==158971) badrun = 1;
  else if(run ==158827) badrun = 1;
  else if(run ==158826) badrun = 1;
  else if(run ==158825) badrun = 1;
  else if(run ==158824) badrun = 1;
  else if(run ==158823) badrun = 1;
  else if(run ==158821) badrun = 1;
  else if(run ==158803) badrun = 1;
  else if(run ==158802) badrun = 1;
  else if(run ==158801) badrun = 1;
  else if(run ==158800) badrun = 1;
  else if(run ==158799) badrun = 1;
  else if(run ==158798) badrun = 1;
  else if(run ==158797) badrun = 1;
  else if(run ==158796) badrun = 1;
  else if(run ==158795) badrun = 1;
  else if(run ==158794) badrun = 1;
  else if(run ==158793) badrun = 1;
  else if(run ==158792) badrun = 1;
  else if(run ==158579) badrun = 1;
  else if(run ==158578) badrun = 1;
  else if(run ==158576) badrun = 1;
  else if(run ==158575) badrun = 1;
  else if(run ==158574) badrun = 1;
  else if(run ==158573) badrun = 1;
  else if(run ==158572) badrun = 1;
  else if(run ==158571) badrun = 1;
  else if(run ==158569) badrun = 1;
  else if(run ==158546) badrun = 1;
  else if(run ==158545) badrun = 1;
  else if(run ==158544) badrun = 1;
  else if(run ==158543) badrun = 1;
  else if(run ==158542) badrun = 1;
  else if(run ==158541) badrun = 1;
  else if(run ==158539) badrun = 1;
  else if(run ==158538) badrun = 1;
  else if(run ==158537) badrun = 1;
  else if(run ==158536) badrun = 1;
  else if(run ==158535) badrun = 1;
  else if(run ==158532) badrun = 1;
  else if(run ==158531) badrun = 1;
  else if(run ==158527) badrun = 1;
  else if(run ==158465) badrun = 1;
  else if(run ==158069) badrun = 1;
  else if(run ==158068) badrun = 1;
  else if(run ==158064) badrun = 1;
  else if(run ==158063) badrun = 1;
  else if(run ==158062) badrun = 1;
  else if(run ==158061) badrun = 1;
  else if(run ==158041) badrun = 1;
  else if(run ==158040) badrun = 1;
  else if(run ==158036) badrun = 1;
  else if(run ==158035) badrun = 1;
  else if(run ==158034) badrun = 1;
  else if(run ==158033) badrun = 1;
  else if(run ==158031) badrun = 1;
  else if(run ==157713) badrun = 1;
  else if(run ==157709) badrun = 1;
  else if(run ==157708) badrun = 1;
  else if(run ==157705) badrun = 1;
  else if(run ==157702) badrun = 1;
  else if(run ==157597) badrun = 1;
  else if(run ==157596) badrun = 1;
  else if(run ==157595) badrun = 1;
  else if(run ==157591) badrun = 1;
  else if(run ==157581) badrun = 1;
  else if(run ==157577) badrun = 1;
  else if(run ==157576) badrun = 1;
  else if(run ==157575) badrun = 1;
  else if(run ==157574) badrun = 1;
  else if(run ==157573) badrun = 1;
  else if(run ==157572) badrun = 1;
  else if(run ==157571) badrun = 1;
  else if(run ==157566) badrun = 1;
  else if(run ==157564) badrun = 1;
  else if(run ==157563) badrun = 1;
  else if(run ==157562) badrun = 1;
  else if(run ==157561) badrun = 1;
  else if(run ==157560) badrun = 1;
  else if(run ==157559) badrun = 1;
  else if(run ==157556) badrun = 1;
  else if(run ==157546) badrun = 1;
  else if(run ==157536) badrun = 1;
  else if(run ==157531) badrun = 1;
  else if(run ==157484) badrun = 1;
  else if(run ==157483) badrun = 1;
  else if(run ==157482) badrun = 1;
  else if(run ==157478) badrun = 1;
  else if(run ==157477) badrun = 1;
  else if(run ==157476) badrun = 1;

  if(badrun == 1) goodcal = false;
  return goodcal;
}