
#include "../GlobalFunctions.hh"

//////////////////////////////////////////////////////////////////////
/// OptimizePreSelection()
/// loops over all .root-files created by Preselection() defined in BDTSelection.cc.
/// More kinematic, PID and some more cuts are applied and saved in a new file.
///


THStack *stack(vector<TH1F*> histos){
    THStack *tmp = new THStack("IncMC","IncMC");
    for(unsigned int h = 0; h < histos.size(); h++){
        histos.at(h)->SetFillColor(2+h);
        tmp->Add(histos.at(h));
    }
    return tmp;

}

int plotIncSample(std::string year = "2012", std::string magnet = "down") {

    bool ReferenceChannel = true;

    TFile* input = 0;
    TFile* histos = 0;
    TTree *tree;
        // Kst2Kpluspi0Resolved case

    if(Kst2Kpluspi0Resolved){
        input = TFile::Open(Form("%s/data/MC/BtoXJpsi/%s%s/%s%s_pi0Resolved%s.root",path_to_output_KplusPizero.c_str(), year.c_str(),magnet.c_str(),year.c_str(),magnet.c_str(),(smallSample ? "_small":"") ));
        histos = TFile::Open(Form("%s/data/MC/BtoXJpsi/%s%s/%s%s_pi0Resolved_histos%s.root",path_to_output_KplusPizero.c_str(), year.c_str(),magnet.c_str(),year.c_str(),magnet.c_str(),(smallSample ? "_small":"") ));
    }
        // Kst2Kpluspi0Merged case
    if(Kst2Kpluspi0Merged){
        //TODO
    }
        // Kst2Kspiplus case
    if(Kst2Kspiplus){
        //TODO
    }


    input->cd();
    cout << "In input" << endl;
    tree = (TTree*)input->Get("DecayTree");

    cout << "Tree has " <<  tree->GetEntries() << " entries" << endl;
    tree->SetBranchStatus("*",1);

    Double_t B_plus_M_DTF;
    //assign variables to TRUEID
    Int_t B_plus_TRUEID;
    Int_t K_star_plus_TRUEID;
    Int_t K_plus_TRUEID;
    Int_t K_short_TRUEID;
    Int_t Ks_pi_plus_TRUEID;
    Int_t Ks_pi_minus_TRUEID;
    Int_t pi_plus_TRUEID;
    Int_t pi_zero_TRUEID;
    Int_t mu_plus_TRUEID;
    Int_t mu_minus_TRUEID;
    Int_t J_psi_TRUEID;
    Int_t gamma1_TRUEID;
    Int_t gamma2_TRUEID;

    //assign variables to TRUE_MOTHER_ID
    Int_t B_plus_MOTHER_ID;
    Int_t K_star_plus_MOTHER_ID;
    Int_t K_plus_MOTHER_ID;
    Int_t K_short_MOTHER_ID;
    Int_t Ks_pi_plus_MOTHER_ID;
    Int_t Ks_pi_minus_MOTHER_ID;
    Int_t pi_plus_MOTHER_ID;
    Int_t pi_zero_MOTHER_ID;
    Int_t mu_plus_MOTHER_ID;
    Int_t mu_minus_MOTHER_ID;
    Int_t J_psi_MOTHER_ID;
    Int_t gamma1_MOTHER_ID;
    Int_t gamma2_MOTHER_ID;

    //assign variables to TRUE_GD_MOTHER_ID
    Int_t K_plus_GD_MOTHER_ID;
    Int_t K_plus_GD_GD_MOTHER_ID;
    Int_t K_short_GD_MOTHER_ID;
    Int_t Ks_pi_plus_GD_MOTHER_ID;
    Int_t Ks_pi_minus_GD_MOTHER_ID;
    Int_t pi_plus_GD_MOTHER_ID;
    Int_t pi_zero_GD_MOTHER_ID;
    Int_t mu_plus_GD_MOTHER_ID;
    Int_t mu_minus_GD_MOTHER_ID;
    Int_t gamma1_GD_MOTHER_ID;
    Int_t gamma2_GD_MOTHER_ID;

    //assign variables to TRUE_GD_MOTHER_ID
    Int_t Ks_pi_plus_GD_GD_MOTHER_ID;
    Int_t Ks_pi_minus_GD_GD_MOTHER_ID;
    Int_t gamma1_GD_GD_MOTHER_ID;
    Int_t gamma2_GD_GD_MOTHER_ID;
    Int_t pi_zero_merged_MC_GD_GD_MOTHER_ID;


    //link variables to branches
    tree -> SetBranchAddress( "B_plus_M_DTF" , &B_plus_M_DTF );

    tree -> SetBranchAddress( "B_plus_TRUEID" , &B_plus_TRUEID );
    tree -> SetBranchAddress( "K_star_plus_TRUEID" , &K_star_plus_TRUEID );
    tree -> SetBranchAddress( "mu_plus_TRUEID" , &mu_plus_TRUEID );
    tree -> SetBranchAddress( "mu_minus_TRUEID" , &mu_minus_TRUEID );

    tree -> SetBranchAddress( "B_plus_MC_MOTHER_ID" , &B_plus_MOTHER_ID );
    tree -> SetBranchAddress( "K_star_plus_MC_MOTHER_ID" , &K_star_plus_MOTHER_ID );
    tree -> SetBranchAddress( "mu_plus_MC_MOTHER_ID" , &mu_plus_MOTHER_ID );
    tree -> SetBranchAddress( "mu_minus_MC_MOTHER_ID" , &mu_minus_MOTHER_ID );
    tree -> SetBranchAddress( "mu_plus_MC_GD_MOTHER_ID" , &mu_plus_GD_MOTHER_ID );
    tree -> SetBranchAddress( "mu_minus_MC_GD_MOTHER_ID" , &mu_minus_GD_MOTHER_ID );


    if(Kst2Kpluspi0Resolved){
        tree -> SetBranchAddress( "K_plus_TRUEID" , &K_plus_TRUEID );
        tree -> SetBranchAddress( "pi_zero_resolved_TRUEID" , &pi_zero_TRUEID );
        tree -> SetBranchAddress( "gamma1_TRUEID" , &gamma1_TRUEID );
        tree -> SetBranchAddress( "gamma2_TRUEID" , &gamma2_TRUEID );

        tree -> SetBranchAddress( "K_plus_MC_MOTHER_ID" , &K_plus_MOTHER_ID );
        tree -> SetBranchAddress( "pi_zero_resolved_MC_MOTHER_ID" , &pi_zero_MOTHER_ID );
        tree -> SetBranchAddress( "gamma1_MC_MOTHER_ID" , &gamma1_MOTHER_ID );
        tree -> SetBranchAddress( "gamma2_MC_MOTHER_ID" , &gamma2_MOTHER_ID );

        tree -> SetBranchAddress( "K_plus_MC_GD_MOTHER_ID" , &K_plus_GD_MOTHER_ID );
        tree -> SetBranchAddress( "K_plus_MC_GD_GD_MOTHER_ID" , &K_plus_GD_GD_MOTHER_ID );
        tree -> SetBranchAddress( "pi_zero_resolved_MC_GD_MOTHER_ID" , &pi_zero_GD_MOTHER_ID );
        tree -> SetBranchAddress( "gamma1_MC_GD_MOTHER_ID" , &gamma1_GD_MOTHER_ID );
        tree -> SetBranchAddress( "gamma2_MC_GD_MOTHER_ID" , &gamma2_GD_MOTHER_ID );

        tree -> SetBranchAddress( "gamma1_MC_GD_GD_MOTHER_ID" , &gamma1_GD_GD_MOTHER_ID );
        tree -> SetBranchAddress( "gamma2_MC_GD_GD_MOTHER_ID" , &gamma2_GD_GD_MOTHER_ID );
    }

    if(Kst2Kpluspi0Merged){
        tree -> SetBranchAddress( "K_plus_TRUEID" , &K_plus_TRUEID );
        tree -> SetBranchAddress( "pi_zero_merged_TRUEID" , &pi_zero_TRUEID );

        tree -> SetBranchAddress( "K_plus_MC_MOTHER_ID" , &K_plus_MOTHER_ID );
        tree -> SetBranchAddress( "pi_zero_merged_MC_MOTHER_ID" , &pi_zero_MOTHER_ID );

        tree -> SetBranchAddress( "K_plus_MC_GD_MOTHER_ID" , &K_plus_GD_MOTHER_ID );
        tree -> SetBranchAddress( "pi_zero_merged_MC_GD_MOTHER_ID" , &pi_zero_GD_MOTHER_ID );

        tree -> SetBranchAddress( "pi_zero_merged_MC_GD_GD_MOTHER_ID" , &pi_zero_merged_MC_GD_GD_MOTHER_ID );
    }

    if(Kst2Kspiplus){
        tree -> SetBranchAddress( "K_short_TRUEID" , &K_short_TRUEID );
        tree -> SetBranchAddress( "pi_plus_TRUEID" , &pi_plus_TRUEID );
        tree -> SetBranchAddress( "Ks_pi_plus_TRUEID" , &Ks_pi_plus_TRUEID );
        tree -> SetBranchAddress( "Ks_pi_minus_TRUEID" , &Ks_pi_minus_TRUEID );
        tree -> SetBranchAddress( "K_short_MC_MOTHER_ID" , &K_short_MOTHER_ID );
        tree -> SetBranchAddress( "pi_plus_MC_MOTHER_ID" , &pi_plus_MOTHER_ID );
        tree -> SetBranchAddress( "Ks_pi_plus_MC_MOTHER_ID" , &Ks_pi_plus_MOTHER_ID );
        tree -> SetBranchAddress( "Ks_pi_minus_MC_MOTHER_ID" , &Ks_pi_minus_MOTHER_ID );

        tree -> SetBranchAddress( "K_short_MC_GD_MOTHER_ID" , &K_short_GD_MOTHER_ID );
        tree -> SetBranchAddress( "pi_plus_MC_GD_MOTHER_ID" , &pi_plus_GD_MOTHER_ID );
        tree -> SetBranchAddress( "Ks_pi_plus_MC_GD_MOTHER_ID" , &Ks_pi_plus_GD_MOTHER_ID );
        tree -> SetBranchAddress( "Ks_pi_minus_MC_GD_MOTHER_ID" , &Ks_pi_minus_GD_MOTHER_ID );

        tree -> SetBranchAddress( "Ks_pi_plus_MC_GD_GD_MOTHER_ID" , &Ks_pi_plus_GD_GD_MOTHER_ID );
        tree -> SetBranchAddress( "Ks_pi_minus_MC_GD_GD_MOTHER_ID" , &Ks_pi_minus_GD_GD_MOTHER_ID );
    }

    tree -> SetBranchAddress( "J_psi_TRUEID" , &J_psi_TRUEID );
    tree -> SetBranchAddress( "J_psi_MC_MOTHER_ID" , &J_psi_MOTHER_ID );

    //counter
    Int_t B_plus_VETO = 0;
    Int_t K_star_plus_VETO = 0;
    Int_t K_plus_VETO = 0;
    Int_t K_short_VETO = 0;
    Int_t Ks_pi_plus_VETO = 0;
    Int_t Ks_pi_minus_VETO = 0;
    Int_t pi_plus_VETO = 0;
    Int_t pi_zero_VETO = 0;
    Int_t mu_plus_VETO = 0;
    Int_t mu_minus_VETO = 0;
    Int_t J_psi_VETO = 0;
    Int_t gamma1_VETO = 0;
    Int_t gamma2_VETO = 0;

    Int_t MOTHER_ID_VETO = 0;
    Int_t GD_MOTHER_ID_VETO = 0;
    Int_t GD_GD_MOTHER_ID_VETO = 0;

    Int_t TRUE_events = 0;

    //count the numbers of Kshort daughters, which are pi_zeros
    Int_t Kshort_pi_zero_daughters = 0;

    //count muons, mis-IDed as pions
    Int_t ReconstructedPionsWhichAreMuons = 0;


    //======================================
    //Get histos of all possible backgrounds
    //======================================


    //control plots for pi0
    Double_t pi_zero_resolved_M = 0.0;
    if (Kst2Kpluspi0Resolved)
        tree -> SetBranchAddress( "pi_zero_resolved_M" , &pi_zero_resolved_M );

    TH1 * h_pi0_resolved_mass= new TH1D("h_pi0_resolved_mass" , "m_{#pi_{0}} resolved", 100, 75, 175);
    h_pi0_resolved_mass->GetXaxis()->SetTitle("m_{pi_{0}} resolved [MeV]");
    h_pi0_resolved_mass->GetYaxis()->SetTitle("Counts / 1MeV");

    TH1 * h_pi0_mismatched_mass= new TH1D("h_pi0_mismatched_mass" , "m_{#pi_{0}} IDed as gamma", 100, 75, 175);
    h_pi0_mismatched_mass->GetXaxis()->SetTitle("m_{pi_{0}}[MeV] MisIDed");
    h_pi0_mismatched_mass->GetYaxis()->SetTitle("Counts / 1MeV");


    TH1 * h_B_plus_TRUEID       = new TH1D("h_B_plus_TRUEID",       "B^{+} TRUEID",   100000, -1, 99999);
    TH1 * h_K_star_plus_TRUEID  = new TH1D("h_K_star_plus_TRUEID",  "K^{*+} TRUEID",  100000, -1, 99999);
    TH1 * h_K_plus_TRUEID       = new TH1D("h_K_plus_TRUEID",       "K^{+} TRUEID",   100000, -1, 99999);
    TH1 * h_pi_zero_TRUEID      = new TH1D("h_pi_zero_TRUEID",      "#pi^{0} TRUEID", 100000, -1, 99999);

    TH1 * h_B_plus_TRUEID_rm    = new TH1D("h_B_plus_TRUEID_rm",       "B^{+} TRUEID",   100000, -1, 99999);
    TH1 * h_K_star_plus_TRUEID_rm= new TH1D("h_K_star_plus_TRUEID_rm", "K^{*+} TRUEID",  100000, -1, 99999);
    TH1 * h_K_plus_TRUEID_rm    = new TH1D("h_K_plus_TRUEID_rm",       "K^{+} TRUEID",   100000, -1, 99999);
    TH1 * h_pi_zero_TRUEID_rm   = new TH1D("h_pi_zero_TRUEID_rm",      "#pi^{0} TRUEID", 100000, -1, 99999);


    TH1 * h_B_plus_K1   = new TH1D("h_B_plus_K1",   "h_B_plus_K1",  100,    cut_B_plus_M_low,   cut_B_plus_M_high);
    TH1 * h_B_plus_Kst  = new TH1D("h_B_plus_Kst",  "h_B_plus_Kst", 100,    cut_B_plus_M_low,   cut_B_plus_M_high);
    TH1 * h_B_plus_Kpl  = new TH1D("h_B_plus_Kpl",  "h_B_plus_Kpl", 100,    cut_B_plus_M_low,   cut_B_plus_M_high);


    vector<std::string> Cuts;
    string tmp_string = "";
    int nBins = 75;
    TH1 * h_BplusMass_tot= new TH1D("h_BplusMass_tot" , "h_BplusMass_tot", nBins ,cut_B_plus_M_low, cut_B_plus_M_high);
    h_BplusMass_tot->GetXaxis()->SetTitle("B^{+} mass [MeV]");

    //======================================
    // Define trees for all channels
    //======================================
    TTree *tree_correct = tree->CloneTree(0);
    tree_correct->SetName("tree_correct");
    TH1 * h_BplusMass_correct = new TH1D("h_BplusMass_correct" , "h_BplusMass_correct", nBins ,cut_B_plus_M_low, cut_B_plus_M_high);
    h_BplusMass_correct->GetXaxis()->SetTitle("B^{+} mass [MeV]");

    TTree *tree_K1_1270 = tree->CloneTree(0);
    tree_K1_1270->SetName("tree_K1_1270");
    TH1 * h_BplusMass_K1_1270 = new TH1D("h_BplusMass_K1_1270" , "h_BplusMass_K1_1270", nBins ,cut_B_plus_M_low, cut_B_plus_M_high);
    h_BplusMass_K1_1270->GetXaxis()->SetTitle("B^{+} mass [MeV]");

    TTree *tree_K1_1410 = tree->CloneTree(0);
    tree_K1_1410->SetName("tree_K1_1410");
    TH1 * h_BplusMass_K1_1410 = new TH1D("h_BplusMass_K1_1410" , "h_BplusMass_K1_1410", nBins ,cut_B_plus_M_low, cut_B_plus_M_high);
    h_BplusMass_K1_1410->GetXaxis()->SetTitle("B^{+} mass [MeV]");

    TTree *tree_K1_1400 = tree->CloneTree(0);
    tree_K1_1400->SetName("tree_K1_1400");
    TH1 * h_BplusMass_K1_1400 = new TH1D("h_BplusMass_K1_1400" , "h_BplusMass_K1_1400", nBins ,cut_B_plus_M_low, cut_B_plus_M_high);
    h_BplusMass_K1_1400->GetXaxis()->SetTitle("B^{+} mass [MeV]");

    TTree *tree_K_plus= tree->CloneTree(0);
    tree_K_plus->SetName("tree_K_plus");
    TH1 * h_BplusMass_K_plus = new TH1D("h_BplusMass_K_plus" , "h_BplusMass_K_plus", nBins ,cut_B_plus_M_low, cut_B_plus_M_high);
    h_BplusMass_K_plus->GetXaxis()->SetTitle("B^{+} mass [MeV]");

    TTree *tree_K_plus_misID= tree->CloneTree(0);
    tree_K_plus_misID->SetName("tree_K_plus_misID");
    TH1 * h_BplusMass_K_plus_misID = new TH1D("h_BplusMass_K_plus_misID" , "h_BplusMass_K_plus_misID", nBins ,cut_B_plus_M_low, cut_B_plus_M_high);
    h_BplusMass_K_plus_misID->GetXaxis()->SetTitle("B^{+} mass [MeV]");

    TTree *tree_B_star= tree->CloneTree(0);
    tree_B_star->SetName("tree_B_star");
    TH1 * h_BplusMass_B_star = new TH1D("h_BplusMass_B_star" , "h_BplusMass_B_star", nBins ,cut_B_plus_M_low, cut_B_plus_M_high);
    h_BplusMass_B_star->GetXaxis()->SetTitle("B^{+} mass [MeV]");

    TTree *tree_Kstar_misID= tree->CloneTree(0);
    tree_Kstar_misID->SetName("tree_Kstar_misID");
    TH1 * h_BplusMass_Kstar_misID = new TH1D("h_BplusMass_Kstar_misID" , "h_BplusMass_Kstar_misID", nBins ,cut_B_plus_M_low, cut_B_plus_M_high);
    h_BplusMass_Kstar_misID->GetXaxis()->SetTitle("B^{+} mass [MeV]");

    TTree *tree_bquark= tree->CloneTree(0);
    tree_bquark->SetName("tree_bquark");
    TH1 * h_BplusMass_bquark = new TH1D("h_BplusMass_bquark" , "h_BplusMass_bquark", nBins ,cut_B_plus_M_low, cut_B_plus_M_high);
    h_BplusMass_bquark->GetXaxis()->SetTitle("B^{+} mass [MeV]");

    TTree *tree_combinatorics= tree->CloneTree(0);
    tree_combinatorics->SetName("tree_combinatorics");
    TH1 * h_BplusMass_combinatorics = new TH1D("h_BplusMass_combinatorics" , "h_BplusMass_combinatorics", nBins ,cut_B_plus_M_low, cut_B_plus_M_high);
    h_BplusMass_combinatorics->GetXaxis()->SetTitle("B^{+} mass [MeV]");


    TTree *tree_others = tree->CloneTree(0);
    tree_others->SetName("tree_others");
    TH1 * h_BplusMass_others = new TH1D("h_BplusMass_others" , "h_BplusMass_others", nBins ,cut_B_plus_M_low, cut_B_plus_M_high);
    h_BplusMass_others->GetXaxis()->SetTitle("B^{+} mass [MeV]");


    TH1 * h_BplusMass_JpsiVETO = new TH1D("h_BplusMass_JpsiVETO" , "h_BplusMass_JpsiVETO", nBins ,cut_B_plus_M_low, cut_B_plus_M_high);
    h_BplusMass_JpsiVETO->GetXaxis()->SetTitle("B^{+} mass [MeV]");

    TH1 * h_BplusMass_GDVETO = new TH1D("h_BplusMass_GDVETO" , "h_BplusMass_GDVETO", nBins ,cut_B_plus_M_low, cut_B_plus_M_high);
    h_BplusMass_GDVETO->GetXaxis()->SetTitle("B^{+} mass [MeV]");

    //======================================
    //check all particles for correct TRUEID
    //======================================

    //MC loop
    cout << "Tree has " <<  tree->GetEntries() << " entries" << endl;
    int nEvents = tree->GetEntries();
    std::cout << "Truth matching of " << nEvents << " INC MC events for " << year << magnet << "!" << std::endl;
    for(int i = 0; i < nEvents; i++){
        if(i%2000==0 && i != 0)std::cout << "Loading MC event " << i << " / " << nEvents << std::endl;
        bool flag_other = true;

        //others
        tree->GetEntry(i);
        h_BplusMass_tot->Fill(B_plus_M_DTF);
        //remove nonsense from B and Jpsi
        if(TMath::Abs(B_plus_TRUEID) != TRUEID.B_PLUS && TMath::Abs(B_plus_TRUEID) != 0 && TMath::Abs(B_plus_TRUEID) != 5){
            B_plus_VETO++;
            h_B_plus_TRUEID_rm->Fill(B_plus_TRUEID);
            continue;
        }        
        if(ReferenceChannel){
            if(TMath::Abs(J_psi_TRUEID) != TRUEID.J_PSI){
                J_psi_VETO++;
                h_BplusMass_JpsiVETO->Fill(B_plus_M_DTF);
                continue;
            }
            if(TMath::Abs(mu_plus_MOTHER_ID) != TRUEID.J_PSI){
                MOTHER_ID_VETO++;
                h_BplusMass_others->Fill(B_plus_M_DTF);
                continue;
            }
            if(TMath::Abs(mu_minus_MOTHER_ID) != TRUEID.J_PSI){
                MOTHER_ID_VETO++;
                h_BplusMass_others->Fill(B_plus_M_DTF);
                continue;
            }
            if(TMath::Abs(mu_plus_GD_MOTHER_ID) != TRUEID.B_PLUS){
                GD_MOTHER_ID_VETO++;
                h_BplusMass_GDVETO->Fill(B_plus_M_DTF);
                continue;
            }
            if(TMath::Abs(mu_minus_GD_MOTHER_ID) != TRUEID.B_PLUS){
                GD_MOTHER_ID_VETO++;
                h_BplusMass_GDVETO->Fill(B_plus_M_DTF);
                continue;
            }
        }
        else{
            if(TMath::Abs(J_psi_TRUEID) == TRUEID.J_PSI){
                J_psi_VETO++;
                h_BplusMass_JpsiVETO->Fill(B_plus_M_DTF);
                continue;
            }
            if(TMath::Abs(mu_plus_MOTHER_ID) != TRUEID.B_PLUS){
                MOTHER_ID_VETO++;
                h_BplusMass_others->Fill(B_plus_M_DTF);
                continue;
            }
            if(TMath::Abs(mu_minus_MOTHER_ID) != TRUEID.B_PLUS){
                MOTHER_ID_VETO++;
                h_BplusMass_others->Fill(B_plus_M_DTF);
                continue;
            }
        }



        //correct tree
        if (TMath::Abs(K_star_plus_TRUEID) == TRUEID.K_STAR_PLUS)
            if(TMath::Abs(K_plus_TRUEID) == TRUEID.K_PLUS){ //Paticle IDs
                if(TMath::Abs(K_star_plus_MOTHER_ID) == TRUEID.B_PLUS)
                    if(TMath::Abs(K_plus_MOTHER_ID) == TRUEID.K_STAR_PLUS){ //Mother IDs
                        if(TMath::Abs(K_plus_GD_MOTHER_ID) == TRUEID.B_PLUS){//Grandmother IDs
                            tree_correct->Fill();
                            h_BplusMass_correct->Fill(B_plus_M_DTF);
                            continue;
                        }
                    }
            }




        //K1 tree
        if (TMath::Abs(K_star_plus_TRUEID) == TRUEID.K_ONE_PLUS || TMath::Abs(K_star_plus_TRUEID) == 0 || TMath::Abs(K_star_plus_TRUEID) == TRUEID.B_PLUS )
            if(TMath::Abs(K_plus_TRUEID) == TRUEID.K_STAR_PLUS || TMath::Abs(K_plus_TRUEID) == TRUEID.K_PLUS){ //Paticle IDs
                if(TMath::Abs(K_star_plus_MOTHER_ID) == TRUEID.B_PLUS || TMath::Abs(K_star_plus_MOTHER_ID) == TRUEID.K_ONE_PLUS || TMath::Abs(K_star_plus_MOTHER_ID) == 0)
                    if(TMath::Abs(K_plus_MOTHER_ID) == TRUEID.K_STAR_PLUS  || TMath::Abs(K_plus_MOTHER_ID) == TRUEID.K_ONE_PLUS){ //Mother IDs
                        if(TMath::Abs(K_plus_GD_MOTHER_ID) == TRUEID.B_PLUS || TMath::Abs(K_plus_GD_MOTHER_ID) == TRUEID.K_ONE_PLUS){//Grandmother IDs
                            tree_K1_1270->Fill();
                            h_BplusMass_K1_1270->Fill(B_plus_M_DTF);
                            continue;
                        }
                    }
            }

        //K1 1400 tree
        if (TMath::Abs(K_star_plus_TRUEID) == TRUEID.K_ONE_PLUS_1400 || TMath::Abs(K_star_plus_TRUEID) == 0|| TMath::Abs(K_star_plus_TRUEID) == TRUEID.B_PLUS  )
            if(TMath::Abs(K_plus_TRUEID) == TRUEID.K_STAR_PLUS || TMath::Abs(K_plus_TRUEID) == TRUEID.K_PLUS){ //Paticle ID
                if(TMath::Abs(K_star_plus_MOTHER_ID) == TRUEID.B_PLUS || TMath::Abs(K_star_plus_MOTHER_ID) == TRUEID.K_ONE_PLUS_1400)
                    if(TMath::Abs(K_plus_MOTHER_ID) == TRUEID.K_STAR_PLUS  || TMath::Abs(K_plus_MOTHER_ID) == TRUEID.K_ONE_PLUS_1400){ //Mother IDs
                        if(TMath::Abs(K_plus_GD_MOTHER_ID) == TRUEID.B_PLUS || TMath::Abs(K_plus_GD_MOTHER_ID) == TRUEID.K_ONE_PLUS_1400){//Grandmother IDs
                            tree_K1_1400->Fill();
                            h_BplusMass_K1_1410->Fill(B_plus_M_DTF);
                            tree_others->Fill();
                            h_BplusMass_others->Fill(B_plus_M_DTF);
                            continue;
                        }
                    }
            }
        //K1 1410 tree
        if (TMath::Abs(K_star_plus_TRUEID) == TRUEID.K_ONE_PLUS_1410 || TMath::Abs(K_star_plus_TRUEID) == 0 )
            if(TMath::Abs(K_plus_TRUEID) == TRUEID.K_PLUS){ //Paticle IDs
                if(TMath::Abs(K_star_plus_MOTHER_ID) == TRUEID.B_PLUS || TMath::Abs(K_star_plus_MOTHER_ID) == 0 || TMath::Abs(K_star_plus_MOTHER_ID) == TRUEID.K_ONE_PLUS_1410)
                    if(TMath::Abs(K_plus_MOTHER_ID) == TRUEID.K_STAR_PLUS  || TMath::Abs(K_plus_MOTHER_ID) == TRUEID.K_ONE_PLUS_1410){ //Mother IDs
                        if(TMath::Abs(K_plus_GD_MOTHER_ID) == 0 || TMath::Abs(K_plus_GD_MOTHER_ID) == TRUEID.K_ONE_PLUS_1410){//Grandmother IDs
                            tree_K1_1410->Fill();
                            h_BplusMass_K1_1410->Fill(B_plus_M_DTF);
                            tree_others->Fill();
                            h_BplusMass_others->Fill(B_plus_M_DTF);
                            continue;
                        }
                    }

            }
        //K1 1270 tree 2.0
        if (TMath::Abs(K_star_plus_TRUEID) == TRUEID.K_STAR_PLUS)
            if(TMath::Abs(K_plus_TRUEID) == TRUEID.K_PLUS){ //Paticle IDs
                if(TMath::Abs(K_star_plus_MOTHER_ID) == TRUEID.K_ONE_PLUS)
                    if(TMath::Abs(K_plus_MOTHER_ID) == TRUEID.K_STAR_PLUS){ //Mother IDs
                        if(TMath::Abs(K_plus_GD_MOTHER_ID) == TRUEID.K_ONE_PLUS){//Grandmother IDs
                            tree_K1_1270->Fill();
                            h_BplusMass_K1_1270->Fill(B_plus_M_DTF);
                            continue;
                        }
                    }
            }


        //tree_K_plus_misID tree
        if (TMath::Abs(K_star_plus_TRUEID) == TRUEID.K_ONE_PLUS || TMath::Abs(K_star_plus_TRUEID) == TRUEID.K_STAR_PLUS|| TMath::Abs(K_star_plus_TRUEID) == 0 )
            if(TMath::Abs(K_plus_TRUEID) == TRUEID.PI_PLUS){ //Paticle IDs
                if(TMath::Abs(K_star_plus_MOTHER_ID) == TRUEID.B_PLUS || TMath::Abs(K_star_plus_MOTHER_ID) == TRUEID.K_ONE_PLUS)
                    if(TMath::Abs(K_plus_MOTHER_ID) == TRUEID.K_STAR_PLUS  || TMath::Abs(K_plus_MOTHER_ID) == TRUEID.K_ONE_PLUS){ //Mother IDs
                        if(TMath::Abs(K_plus_GD_MOTHER_ID) == TRUEID.B_PLUS || TMath::Abs(K_plus_GD_MOTHER_ID) == TRUEID.K_ONE_PLUS){//Grandmother IDs
                            tree_K_plus_misID->Fill();
                            h_BplusMass_K_plus_misID->Fill(B_plus_M_DTF);
                            tree_others->Fill();
                            h_BplusMass_others->Fill(B_plus_M_DTF);
                            continue;
                        }
                    }
            }
        //tree_K_plus_misID tree 2.0 (K+->pi+)
        if (TMath::Abs(K_star_plus_TRUEID) == TRUEID.K_STAR_PLUS)
            if(TMath::Abs(K_plus_TRUEID) == TRUEID.PI_PLUS){ //Paticle IDs
                if(TMath::Abs(K_star_plus_MOTHER_ID) == TRUEID.B_PLUS)
                    if(TMath::Abs(K_plus_MOTHER_ID) == TRUEID.K_PLUS){ //Mother IDs
                        if(TMath::Abs(K_plus_GD_MOTHER_ID) == TRUEID.K_STAR_PLUS)
                            if (TMath::Abs(K_plus_GD_GD_MOTHER_ID) == TRUEID.B_PLUS){//Grandmother IDs
                            tree_K_plus_misID->Fill();
                            h_BplusMass_K_plus_misID->Fill(B_plus_M_DTF);
                            tree_others->Fill();
                            h_BplusMass_others->Fill(B_plus_M_DTF);
                            continue;
                        }
                    }
            }

        //Kplus tree
        if (TMath::Abs(K_star_plus_TRUEID) == TRUEID.K_PLUS || TMath::Abs(K_star_plus_TRUEID) == 0)
            if(TMath::Abs(K_plus_TRUEID) == TRUEID.K_PLUS){ //Paticle IDs
                if(TMath::Abs(K_plus_MOTHER_ID) == TRUEID.B_PLUS){ //Mother IDs
                        if(TMath::Abs(K_plus_GD_MOTHER_ID) == 0){//Grandmother IDs
                            tree_K_plus->Fill();
                            h_BplusMass_K_plus->Fill(B_plus_M_DTF);
                            tree_others->Fill();
                            h_BplusMass_others->Fill(B_plus_M_DTF);
                            continue;
                        }
                    }
            }


        //K star misID tree
        if (TMath::Abs(K_star_plus_TRUEID) == TRUEID.K_PLUS)
            if(TMath::Abs(K_plus_TRUEID) == TRUEID.K_PLUS){ //Paticle IDs
                if(TMath::Abs(K_star_plus_MOTHER_ID) == TRUEID.K_STAR_PLUS)
                    if(TMath::Abs(K_plus_MOTHER_ID) == TRUEID.K_STAR_PLUS){ //Mother IDs
                        if(TMath::Abs(K_plus_GD_MOTHER_ID) == TRUEID.B_PLUS){//Grandmother IDs
                            tree_Kstar_misID->Fill();
                            h_BplusMass_Kstar_misID->Fill(B_plus_M_DTF);
                            tree_others->Fill();
                            h_BplusMass_others->Fill(B_plus_M_DTF);
                            continue;
                        }
                    }
            }


        //B+ <-> b quark
        if (TMath::Abs(K_star_plus_TRUEID) == TRUEID.B_PLUS)
            if(TMath::Abs(K_plus_TRUEID) == TRUEID.K_PLUS){ //Paticle IDs
                if(TMath::Abs(K_star_plus_MOTHER_ID) == 5)
                    if(TMath::Abs(K_plus_MOTHER_ID) == TRUEID.B_PLUS){ //Mother IDs
                        tree_bquark->Fill();
                        h_BplusMass_bquark->Fill(B_plus_M_DTF);
                        continue;
                    }
            }



        //B+ <-> b quark 2.0
        if (TMath::Abs(K_star_plus_TRUEID) ==5)
            if(TMath::Abs(K_plus_TRUEID) == TRUEID.K_PLUS){ //Paticle IDs
                if(TMath::Abs(K_star_plus_MOTHER_ID) == 0)
                    if(TMath::Abs(K_plus_MOTHER_ID) == TRUEID.B_PLUS || TMath::Abs(K_plus_MOTHER_ID) == TRUEID.K_STAR_PLUS ){ //Mother IDs
                        tree_bquark->Fill();
                        h_BplusMass_bquark->Fill(B_plus_M_DTF);
                        continue;
                    }
            }


        //Bstar tree
        if (TMath::Abs(K_star_plus_TRUEID) == TRUEID.B_PLUS)
            if(TMath::Abs(K_plus_TRUEID) == TRUEID.K_PLUS){ //Paticle IDs
                if(TMath::Abs(K_star_plus_MOTHER_ID) == 523)
                    if(TMath::Abs(K_plus_MOTHER_ID) == TRUEID.B_PLUS){ //Mother IDs
                        if(TMath::Abs(K_plus_GD_MOTHER_ID) == 523){//Grandmother IDs
                            tree_B_star->Fill();
                            h_BplusMass_B_star->Fill(B_plus_M_DTF);
                            continue;
                        }
                    }
            }

        //Bstar tree 2.0
        if (TMath::Abs(K_star_plus_TRUEID) == TRUEID.K_PLUS)
            if(TMath::Abs(K_plus_TRUEID) == TRUEID.K_PLUS){ //Paticle IDs
                if(TMath::Abs(K_star_plus_MOTHER_ID) == 523)
                    if(TMath::Abs(K_plus_MOTHER_ID) == TRUEID.B_PLUS){ //Mother IDs
                        if(TMath::Abs(K_plus_GD_MOTHER_ID) == 523){//Grandmother IDs
                            tree_B_star->Fill();
                            h_BplusMass_B_star->Fill(B_plus_M_DTF);
                            continue;
                        }
                    }
            }

        //Bstar tree 3.0 (K* misID B+)
        if (TMath::Abs(K_star_plus_TRUEID) == TRUEID.K_PLUS)
            if(TMath::Abs(K_plus_TRUEID) == TRUEID.K_PLUS){ //Paticle IDs
                if(TMath::Abs(K_star_plus_MOTHER_ID) == TRUEID.B_PLUS)
                    if(TMath::Abs(K_plus_MOTHER_ID) == TRUEID.B_PLUS){ //Mother IDs
                        if(TMath::Abs(K_plus_GD_MOTHER_ID) == 523){//Grandmother IDs
                            tree_B_star->Fill();
                            h_BplusMass_B_star->Fill(B_plus_M_DTF);
                            continue;
                        }
                    }
            }
        //Bstar tree 4.0 (K* misID B+)
        if (TMath::Abs(K_star_plus_TRUEID) == 0)
            if(TMath::Abs(K_plus_TRUEID) == TRUEID.K_PLUS){ //Paticle IDs
                if(TMath::Abs(K_star_plus_MOTHER_ID) == 0)
                    if(TMath::Abs(K_plus_MOTHER_ID) == TRUEID.B_PLUS){ //Mother IDs
                        if(TMath::Abs(K_plus_GD_MOTHER_ID) == 523 ||TMath::Abs(K_plus_GD_MOTHER_ID) == 10521){//Grandmother IDs
                            tree_B_star->Fill();
                            h_BplusMass_B_star->Fill(B_plus_M_DTF);
                            continue;
                        }
                    }
            }
        //combinatorics tree
        if (TMath::Abs(K_star_plus_TRUEID) == 0)
            if(TMath::Abs(K_plus_TRUEID) == TRUEID.K_PLUS ||TMath::Abs(K_plus_TRUEID) == 0 ){ //Paticle IDs
                if(TMath::Abs(K_star_plus_MOTHER_ID) == 0)
                    if(TMath::Abs(K_plus_MOTHER_ID) == TRUEID.K_STAR_PLUS ||TMath::Abs(K_plus_TRUEID) == 0 ){ //Mother IDs
                        if(TMath::Abs(K_plus_GD_MOTHER_ID) == TRUEID.B_PLUS){//Grandmother IDs
                            tree_combinatorics->Fill();
                            h_BplusMass_combinatorics->Fill(B_plus_M_DTF);
                            continue;
                       }
                    }
            }

     //the rest of stuff
        tree_others->Fill();
        h_BplusMass_others->Fill(B_plus_M_DTF);

    }
    //output truthmatching results
    std::cout << "Resulting TruthMatching data: " << TRUE_events << "/" << nEvents << " are found to be correctly reconstructed!" << std::endl;
    std::cout << "The incorrect reconstructed particles are the following:" << std::endl;
    std::cout << "(only the first found-incorrect particle is counted)" << std::endl;
    std::cout << std::endl;
    //std::cout << (B0 ? "B_zero:\t\t" : "B_plus:\t\t") << B_plus_VETO << std::endl;
    //std::cout << (B0 ? "K_star_zero:\t" : ( K1 ? "K_one_plus:\t" : "K_star_plus:\t")) << K_star_plus_VETO << std::endl;
    std::cout << "mu_plus:\t" << mu_plus_VETO << std::endl;
    std::cout << "mu_minus:\t" << mu_minus_VETO << std::endl;
    if(Kst2Kpluspi0Resolved){
        std::cout << "K_plus:\t\t" << K_plus_VETO << std::endl;
        std::cout << "pi_zero:\t" << pi_zero_VETO << std::endl;
        std::cout << "gamma1:\t\t" << gamma1_VETO << std::endl;
        std::cout << "gamma2:\t\t" << gamma2_VETO << std::endl;
    }
    if(Kst2Kpluspi0Merged){
        std::cout << "K_plus:\t\t" << K_plus_VETO << std::endl;
        std::cout << "pi_zero:\t" << pi_zero_VETO << std::endl;
    }
    if(Kst2Kspiplus){
        std::cout << "K_short:\t" << K_short_VETO << std::endl;
        std::cout << "Ks_pi_plus:\t" << Ks_pi_plus_VETO << std::endl;
        std::cout << "Ks_pi_minus:\t" << Ks_pi_minus_VETO << std::endl;
        std::cout << "pi_plus:\t" << pi_plus_VETO << std::endl;
    }
    std::cout << "J_psi:\t\t" <<J_psi_VETO << std::endl;

    if(Kst2Kspiplus){
        std::cout << "Pi0 daughers of the Kshort:\t" << Kshort_pi_zero_daughters << std::endl;
        std::cout << "Muons mis-IDed as pions:\t" << ReconstructedPionsWhichAreMuons << "/" << Ks_pi_plus_VETO + Ks_pi_minus_VETO + pi_plus_VETO << std::endl;
    }

    std::cout << "Total incorrect family members of all particles: " << MOTHER_ID_VETO + GD_MOTHER_ID_VETO + GD_GD_MOTHER_ID_VETO << std::endl;
    std::cout << "MotherID:\t" << MOTHER_ID_VETO << std::endl;
    std::cout << "GDmotherID:\t" << GD_MOTHER_ID_VETO << std::endl;
    std::cout << "GDGDmotherID:\t" << GD_GD_MOTHER_ID_VETO << std::endl;

    //save to file
    TFile * outputFile = nullptr;

    // Kst2Kpluspi0Resolved case
    if(Kst2Kpluspi0Resolved){
        outputFile = new TFile(Form("%s/data/MC/BtoXJpsi/%s%s/%s%s_pi0Resolved_plots.root",path_to_output_KplusPizero.c_str(),year.c_str(),magnet.c_str(),year.c_str(),magnet.c_str()), "RECREATE");
        outputFile->cd();
        tree_correct->Write();
        cout << "Correctly id'd entries " << tree_correct->GetEntries() << endl;
        tree_K1_1270->Write();
        cout << "B+->K1(1270) entries " << tree_K1_1270->GetEntries() << endl;
        tree_K1_1400->Write();
        cout << "B+->K1(1400) entries " << tree_K1_1400->GetEntries() << endl;
        tree_K1_1410->Write();
        cout << "B+->K1(1410) entries " << tree_K1_1410->GetEntries() << endl;
        tree_K_plus->Write();
        cout << "B+->K+pi0 entries " << tree_K_plus->GetEntries() << endl;
        tree_K_plus_misID->Write();
        cout << "B+->K*->pi+pi0 " << tree_K_plus_misID->GetEntries() << endl;
        tree_B_star->Write();
        cout << "B*-> B+->K+ entries " << tree_B_star->GetEntries() << endl;
        tree_Kstar_misID->Write();
        cout << "K* misID entries " << tree_Kstar_misID->GetEntries() << endl;
        tree_combinatorics->Write();
        cout << "b-quark entries " << tree_bquark->GetEntries() << endl;
        tree_bquark->Write();
        cout << "Combinatorics " << tree_combinatorics->GetEntries() << endl;
        tree_others->Write();
        cout << "Other entries " << tree_others->GetEntries() << endl;

        THStack *tmp = new THStack("IncMC","IncMC;B^{+} mass [MeV];Entries");
        //tmp->GetXaxis()->SetTitle("B^{+} mass [MeV]");
        //tmp->GetYaxis()->SetTitle("Entries");

        h_BplusMass_others->SetFillColor(kGreen);
        h_BplusMass_others->SetLineColor(kGreen);
        tmp->Add(h_BplusMass_others);
        h_BplusMass_bquark->SetFillColor(kGreen);
        h_BplusMass_bquark->SetLineColor(kGreen);
        tmp->Add(h_BplusMass_bquark);


        h_BplusMass_K1_1270->SetFillColor(kBlue);
        h_BplusMass_K1_1270->SetLineColor(kBlue);
        tmp->Add(h_BplusMass_K1_1270);

//        h_BplusMass_K1_1400->SetFillColor(0);
//        h_BplusMass_K1_1400->SetLineColor(0);
//        tmp->Add(h_BplusMass_K1_1400);

//        h_BplusMass_K1_1410->SetFillColor(0);
//        h_BplusMass_K1_1410->SetLineColor(0);
//        tmp->Add(h_BplusMass_K1_1410);

//        h_BplusMass_K_plus->SetFillColor(0);
//        h_BplusMass_K_plus->SetLineColor(0);
//        tmp->Add(h_BplusMass_K_plus);

//        h_BplusMass_K_plus_misID->SetFillColor(0);
//        h_BplusMass_K_plus_misID->SetLineColor(0);
//        tmp->Add(h_BplusMass_K_plus_misID);

        h_BplusMass_B_star->SetFillColor(kYellow);
        h_BplusMass_B_star->SetLineColor(kYellow);
        tmp->Add(h_BplusMass_B_star);

//        h_BplusMass_Kstar_misID->SetFillColor(0);
//        h_BplusMass_Kstar_misID->SetLineColor(0);
//        tmp->Add(h_BplusMass_Kstar_misID);
//
      //  h_BplusMass_combinatorics->SetFillColor(kBlack);
      //  h_BplusMass_combinatorics->SetLineColor(kBlack);
      //  tmp->Add(h_BplusMass_combinatorics);


        h_BplusMass_JpsiVETO->SetFillColor(kMagenta-10);
        h_BplusMass_JpsiVETO->SetLineColor(kMagenta-10);
        tmp->Add(h_BplusMass_JpsiVETO);

        h_BplusMass_GDVETO->SetFillColor(kOrange-3);
        h_BplusMass_GDVETO->SetLineColor(kOrange-3);
        tmp->Add(h_BplusMass_GDVETO);

        h_BplusMass_correct->SetFillColor(kRed);
        h_BplusMass_correct->SetLineColor(kRed);
        tmp->Add(h_BplusMass_correct);

        tmp->Write();


        TCanvas *c1 = new TCanvas("c1","all",900,700);
        c1->cd();
        gStyle->SetOptStat(0);
       // h_BplusMass_correct->SetFillColor(0);
        h_BplusMass_correct->Draw();
       //  h_BplusMass_others->SetFillColor(0);
        h_BplusMass_others->Draw("SAME");
       //  h_BplusMass_K1_1270->SetFillColor(0);
        h_BplusMass_K1_1270->Draw("SAME");
        //h_BplusMass_K1_1400->Draw("SAME");
        //h_BplusMass_K1_1410->Draw("SAME");
        //h_BplusMass_K_plus->Draw("SAME");
        //h_BplusMass_K_plus_misID->Draw("SAME");
      //   h_BplusMass_B_star->SetFillColor(0);
        h_BplusMass_B_star->Draw("SAME");
        //h_BplusMass_Kstar_misID->Draw("SAME");
        //h_BplusMass_combinatorics->Draw("SAME");
     //    h_BplusMass_bquark->SetFillColor(0);
        h_BplusMass_bquark->Draw("SAME");
      //   h_BplusMass_JpsiVETO->SetFillColor(0);
        h_BplusMass_JpsiVETO->Draw("SAME");
     //    h_BplusMass_GDVETO->SetFillColor(0);
        h_BplusMass_GDVETO->Draw("SAME");

        c1->Write();

        TCanvas *c2 = new TCanvas("c2","al2",900,700);
        c2->cd();
        tmp->Draw();
        h_BplusMass_tot->SetLineWidth(3);
        h_BplusMass_tot->SetLineColor(kBlack);
        h_BplusMass_tot->Draw("SAMEE");
        c2->Write();
        c2->SaveAs(Form("%s/data/MC/BtoXJpsi/%s%s/%s%s_pi0Resolved%s_BKG.eps",path_to_output_KplusPizero.c_str(), year.c_str(),magnet.c_str(),year.c_str(),magnet.c_str(),(smallSample ? "_small":"") ),"eps");

    }
    // Kst2Kpluspi0Merged case
    else if(Kst2Kpluspi0Merged){
        //TODO
    }

    // Kst2Kspiplus case
    else if(Kst2Kspiplus){
        //TODO
    }

    else return 0;

    if(!outputFile->IsOpen()){
       std::cout << "[ERROR]\t\tFile was not opened succesfully!" << std::endl;
       return 0;
    }

    outputFile->cd();
    //treeTruthMatched->Write("",TObject::kWriteDelete);
    if (Kst2Kpluspi0Resolved){
        h_pi0_mismatched_mass->SetLineColor(kRed);
        h_pi0_mismatched_mass->SetLineWidth(3);
        h_pi0_resolved_mass->SetLineColor(kBlue);
        h_pi0_resolved_mass->SetLineWidth(3);
        h_pi0_mismatched_mass->Write("",TObject::kWriteDelete);
        h_pi0_resolved_mass->Write("",TObject::kWriteDelete);

        h_B_plus_TRUEID         ->Write("",TObject::kWriteDelete);
        h_K_star_plus_TRUEID    ->Write("",TObject::kWriteDelete);
        h_K_plus_TRUEID         ->Write("",TObject::kWriteDelete);
        h_pi_zero_TRUEID        ->Write("",TObject::kWriteDelete);

        h_B_plus_TRUEID_rm      ->Write("",TObject::kWriteDelete);
        h_K_star_plus_TRUEID_rm ->Write("",TObject::kWriteDelete);
        h_K_plus_TRUEID_rm      ->Write("",TObject::kWriteDelete);
        h_pi_zero_TRUEID_rm     ->Write("",TObject::kWriteDelete);

        h_B_plus_K1     ->Write("",TObject::kWriteDelete);
        h_B_plus_Kst    ->Write("",TObject::kWriteDelete);
        h_B_plus_Kpl    ->Write("",TObject::kWriteDelete);
    }

    delete outputFile;//->Close();

    delete h_pi0_mismatched_mass;
    delete h_pi0_resolved_mass;

    delete tree;
    std::cout << "[DONE]\t\tFinished Truth-Matching for INC sample: " << year << " " << magnet  << std::endl;

    return 1;
}



int plotIncSampleAll() {


    if ( plotIncSample("2011","down") == 0) return 0;
    if ( plotIncSample("2011","up") == 0) return 0;
    if ( plotIncSample("2012","down") == 0) return 0;
    if ( plotIncSample("2012","up") == 0) return 0;
    if ( plotIncSample("2015","down") == 0) return 0;
    if ( plotIncSample("2015","up") == 0) return 0;
    if ( plotIncSample("2016","down") == 0) return 0;
    if ( plotIncSample("2016","up") == 0) return 0;
    return 1;

}

