{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "//************************************************************************//\n",
    "//                                                                        //\n",
    "//  Copyright 2017 Bertram Kopf (bertram@ep1.rub.de)                      //\n",
    "//                 - Ruhr-Universität Bochum                              //\n",
    "//                                                                        //\n",
    "//  This file is part of Pawian.                                          //\n",
    "//                                                                        //\n",
    "//  Pawian is free software: you can redistribute it and/or modify        //\n",
    "//  it under the terms of the GNU General Public License as published by  //\n",
    "//  the Free Software Foundation, either version 3 of the License, or     //\n",
    "//  (at your option) any later version.                                   //\n",
    "//                                                                        //\n",
    "//  Pawian is distributed in the hope that it will be useful,             //\n",
    "//  but WITHOUT ANY WARRANTY; without even the implied warranty of        //\n",
    "//  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the         //\n",
    "//  GNU General Public License for more details.                          //\n",
    "//                                                                        //\n",
    "//  You should have received a copy of the GNU General Public License     //\n",
    "//  along with Pawian.  If not, see <http://www.gnu.org/licenses/>.       //\n",
    "//                                                                        //\n",
    "//************************************************************************//\n",
    "\n",
    "#include <iostream>\n",
    "#include <sstream>\n",
    "#include \"TFile.h\"\n",
    "#include \"TH1F.h\"\n",
    "#include \"TCanvas.h\"\n",
    "#include \"TMath.h\"\n",
    "#include \"TKey.h\"\n",
    "#include <algorithm>    // std::max\n",
    "\n",
    "int rebinFactor = 1;\n",
    "\n",
    "TCanvas* canvasMass = new TCanvas(\"Masses\",\"Masses\",20, 20, 900,900);\n",
    "TCanvas* canvasHeliAngles = new TCanvas(\"Angles helicity stystem\",\"Angles helicity stystem\",50, 70, 900,900);\n",
    "TCanvas* canvasGJAngles = new TCanvas(\"Angles Gottfried Jackson stystem\",\"Angles Gottfried Jackson stystem\",\n",
    "                                        80, 120,900,900);"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "void DataFitRebin(int factor) {\n",
    "    rebinFactor=factor;\n",
    "    return;\n",
    "}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "void DataFitComp2(TH1F* histoData, TH1F* histoFit) {\n",
    "    histoData->SetLineColor(kRed);\n",
    "    histoFit->SetLineColor(kBlack);\n",
    "\n",
    "    histoData->SetLineWidth(2);\n",
    "    histoFit->SetLineWidth(2);\n",
    "\n",
    "    histoFit->SetFillColor(kYellow-10);\n",
    "    \n",
    "    histoFit->SetMinimum(0);\n",
    "    histoFit->SetMaximum(histoFit->GetMaximum()*1.3);\n",
    "    \n",
    "    histoFit->Draw(\"hist\");\n",
    "    histoFit->Draw(\"same\");\n",
    "    histoData->Draw(\"same\");\n",
    "    \n",
    "    return;\n",
    "}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "double CalcHistoChi2(TH1F* histoData, TH1F* histoFit) {\n",
    "\n",
    "   int nbins1 = histoData->GetNbinsX();\n",
    "   int nbins2 = histoFit->GetNbinsX();\n",
    "\n",
    "   if(nbins1 != nbins2){\n",
    "      std::cout << \"ERROR: nbins1!=nbins2\" << std::endl;\n",
    "   }\n",
    "\n",
    "   double chi2=0;\n",
    "   int nonzerobins=0;\n",
    "   for(int i=1;i<=nbins1;i++) {\n",
    "      if(histoData->GetBinContent(i) == 0 || histoFit->GetBinContent(i) == 0)\n",
    "         continue;\n",
    "      // if(histoData->GetBinContent(i) == 0 && histoFit->GetBinContent(i) == 0)\n",
    "      //         continue;\n",
    "\n",
    "\n",
    "      //     double error = sqrt(histoData->GetBinContent(i) + histoFit->GetBinContent(i));\n",
    "      double error =  sqrt(histoData->GetBinError(i)*histoData->GetBinError(i) \n",
    "                           + histoFit->GetBinError(i)*histoFit->GetBinError(i));\n",
    "       \n",
    "      chi2+=pow(histoData->GetBinContent(i) - histoFit->GetBinContent(i), 2)/(error*error);\n",
    "      // std::cout << pow(histoData->GetBinContent(i) - histoFit->GetBinContent(i), 2)/(error*error) << std::endl;\n",
    "      nonzerobins++;\n",
    "   }\n",
    "   double chi2Root=histoData->Chi2Test(histoFit,\"WW\");\n",
    "   std::cout << \"Chi2 for histo \" << histoData->GetName() << \"\\t\" << chi2/nonzerobins \n",
    "             << \"\\tp =\" << TMath::Prob(chi2, nonzerobins) \n",
    "             << \"\\tpRoot=\" << chi2Root  << std::endl;\n",
    "\n",
    "   return chi2;\n",
    "}\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "double CalcKolmogorov(TH1F* histoData, TH1F* histoFit) {\n",
    "  double ksTestResult = histoData->KolmogorovTest(histoFit);\n",
    "  std::cout << \"Kolmogorov-Smirnov test for histo \" << histoData->GetName() \n",
    "            << \"\\tp = \" << ksTestResult << std::endl;\n",
    "  return ksTestResult;  \n",
    "}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "void DataFitCompare(std::string rootFileName, std::string histname, bool createcanvas = true) {\n",
    "    TFile* tFile= new TFile(rootFileName.c_str());\n",
    "    std::string dataName = \"Data\" + histname;\n",
    "    std::string fitName = \"Fit\" + histname;\n",
    "\n",
    "    TH1F* histoData = (TH1F*)tFile->Get(dataName.c_str());\n",
    "    TH1F* histoFit = (TH1F*)tFile->Get(fitName.c_str());\n",
    "    histoData->Rebin(rebinFactor);\n",
    "    histoFit->Rebin(rebinFactor);\n",
    "\n",
    "    std::cout << \"\\n\";\n",
    "    CalcHistoChi2(histoData, histoFit);\n",
    "    CalcKolmogorov(histoData, histoFit);\n",
    "    std::cout << std::endl;\n",
    "    DataFitComp2(histoData, histoFit);\n",
    "    return;\n",
    "}\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "void DataFitCompareAll(std::string rootFileName, bool saveImage = true) {\n",
    "  TFile* tFile= new TFile(rootFileName.c_str());   \n",
    "   std::vector<std::string> histnamesMasses;\n",
    "   std::vector<std::string> histnamesAnglesHeli;\n",
    "   std::vector<std::string> histnamesAnglesGJ;\n",
    "\n",
    "   std::string strTheta(\"Theta\");\n",
    "   std::string strPhi(\"Phi\");\n",
    "   std::string strGJ(\"GJ\");\n",
    "\n",
    "   TIter nextkey(tFile->GetListOfKeys());\n",
    "   TKey *key;\n",
    "\n",
    "   while ((key = (TKey*)nextkey())) {\n",
    "      \n",
    "      TH1F *histo = (TH1F*)key->ReadObj();\n",
    "      const char* name = histo->GetName();\n",
    "      //      std::cout << name << std::endl;\n",
    "      std::string sname(name);\n",
    "      std::string pname = sname.substr(0, 4);\n",
    "      std::string aname = sname.substr(4);\n",
    "      \n",
    "      if(pname == \"Data\"){\n",
    "        //sort angle hists\n",
    "        if ( (aname.find(\"Theta\") != string::npos) || (aname.find(\"Phi\") != string::npos) ){\n",
    "          //angle hists in GJ system\n",
    "          if ( aname.find(\"GJ\") != string::npos) histnamesAnglesGJ.push_back(aname);\n",
    "          else histnamesAnglesHeli.push_back(aname);\n",
    "                    }\n",
    "        else histnamesMasses.push_back(aname);\n",
    "      }\n",
    "   }\n",
    "\n",
    "   int numhistsMasses = histnamesMasses.size();\n",
    "   canvasMass->Divide(std::ceil(sqrt(numhistsMasses)), std::ceil(sqrt(numhistsMasses)));\n",
    "   for(int i=0; i<numhistsMasses; i++){\n",
    "      canvasMass->cd(i+1);\n",
    "      DataFitCompare(rootFileName, histnamesMasses.at(i), false);\n",
    "   }\n",
    "\n",
    "   int numhistsHeliAngles = histnamesAnglesHeli.size();\n",
    "   canvasHeliAngles->Divide(std::ceil(sqrt(numhistsHeliAngles)), std::ceil(sqrt(numhistsHeliAngles)));\n",
    "   for(int i=0; i<numhistsHeliAngles; i++){\n",
    "      canvasHeliAngles->cd(i+1);\n",
    "      DataFitCompare(rootFileName, histnamesAnglesHeli.at(i), false);\n",
    "   }\n",
    "   \n",
    "   int numhistsGJAngles = histnamesAnglesGJ.size();\n",
    "   canvasGJAngles->Divide(std::ceil(sqrt(numhistsGJAngles)), std::ceil(sqrt(numhistsGJAngles)));\n",
    "   for(int i=0; i<numhistsGJAngles; i++){\n",
    "      canvasGJAngles->cd(i+1);\n",
    "      DataFitCompare(rootFileName, histnamesAnglesGJ.at(i), false);\n",
    "   }\n",
    "\n",
    "   if(saveImage){\n",
    "      canvasMass->SaveAs(\"DataFitCompFileMasses.png\");\n",
    "      canvasHeliAngles->SaveAs(\"DataFitCompFileHeliAngles.png\");\n",
    "      canvasGJAngles->SaveAs(\"DataFitCompFileGJAngles.png\");\n",
    "   }\n",
    "}   \n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "Chi2 for histo Datapion0Apion0B\t0.965238\tp =0.578647\tpRoot=0.550302\n",
      "Kolmogorov-Smirnov test for histo Datapion0Apion0B\tp = 0.999998\n",
      "\n",
      "\n",
      "Chi2 for histo Datapion0Aeta\t1.02195\tp =0.421113\tpRoot=0.203972\n",
      "Kolmogorov-Smirnov test for histo Datapion0Aeta\tp = 0.727107\n",
      "\n",
      "\n",
      "Chi2 for histo Datapion0Beta\t1.24239\tp =0.0531359\tpRoot=0.0461877\n",
      "Kolmogorov-Smirnov test for histo Datapion0Beta\tp = 0.982953\n",
      "\n",
      "\n",
      "Chi2 for histo DataThetaHeli_pion0Apion0B_Frompion0Apion0Beta\t0.897458\tp =0.759248\tpRoot=0.737875\n",
      "Kolmogorov-Smirnov test for histo DataThetaHeli_pion0Apion0B_Frompion0Apion0Beta\tp = 0.949935\n",
      "\n",
      "\n",
      "Chi2 for histo DataPhiHeli_pion0Apion0B_Frompion0Apion0Beta\t1.23261\tp =0.0572281\tpRoot=0.0503353\n",
      "Kolmogorov-Smirnov test for histo DataPhiHeli_pion0Apion0B_Frompion0Apion0Beta\tp = 0.904724\n",
      "\n",
      "\n",
      "Chi2 for histo DataThetaHeli_pion0Aeta_Frompion0Apion0Beta\t1.17114\tp =0.116277\tpRoot=0.104303\n",
      "Kolmogorov-Smirnov test for histo DataThetaHeli_pion0Aeta_Frompion0Apion0Beta\tp = 0.110322\n",
      "\n",
      "\n",
      "Chi2 for histo DataPhiHeli_pion0Aeta_Frompion0Apion0Beta\t1.12964\tp =0.177089\tpRoot=0.160556\n",
      "Kolmogorov-Smirnov test for histo DataPhiHeli_pion0Aeta_Frompion0Apion0Beta\tp = 0.0798194\n",
      "\n",
      "\n",
      "Chi2 for histo DataThetaHeli_pion0Beta_Frompion0Apion0Beta\t1.19278\tp =0.091637\tpRoot=0.0816076\n",
      "Kolmogorov-Smirnov test for histo DataThetaHeli_pion0Beta_Frompion0Apion0Beta\tp = 0.312808\n",
      "\n",
      "\n",
      "Chi2 for histo DataPhiHeli_pion0Beta_Frompion0Apion0Beta\t1.19216\tp =0.0922855\tpRoot=0.0822041\n",
      "Kolmogorov-Smirnov test for histo DataPhiHeli_pion0Beta_Frompion0Apion0Beta\tp = 0.523201\n",
      "\n",
      "\n",
      "Chi2 for histo DataThetaHeli_pion0A_Frompion0Apion0B\t0.959323\tp =0.596502\tpRoot=0.568644\n",
      "Kolmogorov-Smirnov test for histo DataThetaHeli_pion0A_Frompion0Apion0B\tp = 0.962177\n",
      "\n",
      "\n",
      "Chi2 for histo DataPhiHeli_pion0A_Frompion0Apion0B\t1.0301\tp =0.398284\tpRoot=0.371249\n",
      "Kolmogorov-Smirnov test for histo DataPhiHeli_pion0A_Frompion0Apion0B\tp = 0.404928\n",
      "\n",
      "\n",
      "Chi2 for histo DataThetaHeli_pion0B_Frompion0Apion0B\t0.959323\tp =0.596502\tpRoot=0.568644\n",
      "Kolmogorov-Smirnov test for histo DataThetaHeli_pion0B_Frompion0Apion0B\tp = 0.962177\n",
      "\n",
      "\n",
      "Chi2 for histo DataPhiHeli_pion0B_Frompion0Apion0B\t1.02834\tp =0.403003\tpRoot=0.375866\n",
      "Kolmogorov-Smirnov test for histo DataPhiHeli_pion0B_Frompion0Apion0B\tp = 0.392032\n",
      "\n",
      "\n",
      "Chi2 for histo DataThetaHeli_pion0A_Frompion0Aeta\t0.700774\tp =0.98997\tpRoot=0.987786\n",
      "Kolmogorov-Smirnov test for histo DataThetaHeli_pion0A_Frompion0Aeta\tp = 0.999958\n",
      "\n",
      "\n",
      "Chi2 for histo DataPhiHeli_pion0A_Frompion0Aeta\t1.30251\tp =0.0226983\tpRoot=0.0192775\n",
      "Kolmogorov-Smirnov test for histo DataPhiHeli_pion0A_Frompion0Aeta\tp = 0.393375\n",
      "\n",
      "\n",
      "Chi2 for histo DataThetaHeli_pion0B_Frompion0Beta\t0.889557\tp =0.777595\tpRoot=0.755119\n",
      "Kolmogorov-Smirnov test for histo DataThetaHeli_pion0B_Frompion0Beta\tp = 0.999987\n",
      "\n",
      "\n",
      "Chi2 for histo DataPhiHeli_pion0B_Frompion0Beta\t0.866442\tp =0.826992\tpRoot=0.807786\n",
      "Kolmogorov-Smirnov test for histo DataPhiHeli_pion0B_Frompion0Beta\tp = 0.179073\n",
      "\n",
      "\n",
      "Chi2 for histo DataThetaGJ_pion0Apion0B_Frompion0Apion0Beta\t0.904403\tp =0.742557\tpRoot=0.718486\n",
      "Kolmogorov-Smirnov test for histo DataThetaGJ_pion0Apion0B_Frompion0Apion0Beta\tp = 0.86651\n",
      "\n",
      "\n",
      "Chi2 for histo DataPhiGJ_pion0Apion0B_Frompion0Apion0Beta\t0.95813\tp =0.599855\tpRoot=0.571851\n",
      "Kolmogorov-Smirnov test for histo DataPhiGJ_pion0Apion0B_Frompion0Apion0Beta\tp = 0.365316\n",
      "\n",
      "\n",
      "Chi2 for histo DataThetaGJ_pion0Aeta_Frompion0Apion0Beta\t1.20517\tp =0.0795146\tpRoot=0.0698091\n",
      "Kolmogorov-Smirnov test for histo DataThetaGJ_pion0Aeta_Frompion0Apion0Beta\tp = 0.655725\n",
      "\n",
      "\n",
      "Chi2 for histo DataPhiGJ_pion0Aeta_Frompion0Apion0Beta\t1.12891\tp =0.178332\tpRoot=0.160633\n",
      "Kolmogorov-Smirnov test for histo DataPhiGJ_pion0Aeta_Frompion0Apion0Beta\tp = 0.160925\n",
      "\n",
      "\n",
      "Chi2 for histo DataThetaGJ_pion0Beta_Frompion0Apion0Beta\t0.881926\tp =0.794631\tpRoot=0.773278\n",
      "Kolmogorov-Smirnov test for histo DataThetaGJ_pion0Beta_Frompion0Apion0Beta\tp = 0.929447\n",
      "\n",
      "\n",
      "Chi2 for histo DataPhiGJ_pion0Beta_Frompion0Apion0Beta\t0.830251\tp =0.890275\tpRoot=0.875842\n",
      "Kolmogorov-Smirnov test for histo DataPhiGJ_pion0Beta_Frompion0Apion0Beta\tp = 0.768697\n",
      "\n",
      "\n",
      "Chi2 for histo DataThetaGJ_pion0A_Frompion0Apion0B\t0.997031\tp =0.489565\tpRoot=0.461791\n",
      "Kolmogorov-Smirnov test for histo DataThetaGJ_pion0A_Frompion0Apion0B\tp = 0.79276\n",
      "\n",
      "\n",
      "Chi2 for histo DataPhiGJ_pion0A_Frompion0Apion0B\t0.936044\tp =0.660869\tpRoot=0.634282\n",
      "Kolmogorov-Smirnov test for histo DataPhiGJ_pion0A_Frompion0Apion0B\tp = 0.472767\n",
      "\n",
      "\n",
      "Chi2 for histo DataThetaGJ_pion0B_Frompion0Apion0B\t0.997031\tp =0.489565\tpRoot=0.461791\n",
      "Kolmogorov-Smirnov test for histo DataThetaGJ_pion0B_Frompion0Apion0B\tp = 0.79276\n",
      "\n",
      "\n",
      "Chi2 for histo DataPhiGJ_pion0B_Frompion0Apion0B\t0.935408\tp =0.662587\tpRoot=0.636043\n",
      "Kolmogorov-Smirnov test for histo DataPhiGJ_pion0B_Frompion0Apion0B\tp = 0.366517\n",
      "\n",
      "\n",
      "Chi2 for histo DataThetaGJ_pion0A_Frompion0Aeta\t0.77408\tp =0.954347\tpRoot=0.946903\n",
      "Kolmogorov-Smirnov test for histo DataThetaGJ_pion0A_Frompion0Aeta\tp = 0.999126\n",
      "\n",
      "\n",
      "Chi2 for histo DataPhiGJ_pion0A_Frompion0Aeta\t1.02418\tp =0.414223\tpRoot=0.386848\n",
      "Kolmogorov-Smirnov test for histo DataPhiGJ_pion0A_Frompion0Aeta\tp = 0.208214\n",
      "\n",
      "\n",
      "Chi2 for histo DataThetaGJ_pion0B_Frompion0Beta\t0.898192\tp =0.757509\tpRoot=0.734099\n",
      "Kolmogorov-Smirnov test for histo DataThetaGJ_pion0B_Frompion0Beta\tp = 0.991458\n",
      "\n",
      "\n",
      "Chi2 for histo DataPhiGJ_pion0B_Frompion0Beta\t0.936502\tp =0.65963\tpRoot=0.633013\n",
      "Kolmogorov-Smirnov test for histo DataPhiGJ_pion0B_Frompion0Beta\tp = 0.281276\n",
      "\n"
     ]
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "Info in <TH1F::Chi2TestX>: There is a bin in h1 with less than 10 effective events.\n",
      "\n",
      "Info in <TH1F::Chi2TestX>: There is a bin in h2 with less than 10 effective events.\n",
      "\n",
      "Info in <TH1F::Chi2TestX>: There is a bin in h2 with less than 10 effective events.\n",
      "\n",
      "Info in <TH1F::Chi2TestX>: There is a bin in h1 with less than 10 effective events.\n",
      "\n",
      "Info in <TH1F::Chi2TestX>: There is a bin in h2 with less than 10 effective events.\n",
      "\n",
      "Info in <TCanvas::Print>: png file DataFitCompFileMasses.png has been created\n",
      "Info in <TCanvas::Print>: png file DataFitCompFileHeliAngles.png has been created\n"
     ]
    },
    {
     "data": {
      "image/png": "[encoded data, 36752 characters]\n",
      "text/plain": [
       "<IPython.core.display.Image object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "image/png": "[encoded data, 62516 characters]\n",
      "text/plain": [
       "<IPython.core.display.Image object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "image/png": "[encoded data, 64296 characters]\n",
      "text/plain": [
       "<IPython.core.display.Image object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "Info in <TCanvas::Print>: png file DataFitCompFileGJAngles.png has been created\n"
     ]
    }
   ],
   "source": [
    "\n",
    "DataFitCompareAll(\"/data/jollyj/bertram/fits/pbarp/coupled_PiPiEta_PiEtaEta_KKPi/paperFinalHyper/final/Ochs2013---/4f2_0orderBgV5a/qaPiPiEta/pawianHists.root\", true);\n",
    "\n",
    "canvasMass->Draw();\n",
    "canvasHeliAngles->Draw();\n",
    "canvasGJAngles->Draw();\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "ROOT C++ (Python 3)",
   "language": "c++",
   "name": "python36-jupyroot"
  },
  "language_info": {
   "codemirror_mode": "text/x-c++src",
   "file_extension": ".C",
   "mimetype": " text/x-c++src",
   "name": "c++"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}
