JDrawPDF.cc File Reference

Auxiliary program to draw PDF of Cherenkov light from muon. More...

#include <string>
#include <iostream>
#include <iomanip>
#include "TROOT.h"
#include "TFile.h"
#include "TH1D.h"
#include "JPhysics/JPDF.hh"
#include "JPhysics/Antares.hh"
#include "JPhysics/KM3NeT.hh"
#include "JTools/JElement.hh"
#include "JTools/JFunction1D_t.hh"
#include "JTools/JQuantiles.hh"
#include "JTools/JAbstractHistogram.hh"
#include "JGeometry3D/JAngle3D.hh"
#include "Jeep/JProperties.hh"
#include "Jeep/JPrint.hh"
#include "Jeep/JParser.hh"
#include "Jeep/JMessage.hh"

Go to the source code of this file.

Functions
double	getAbsorptionLength (const double lambda)
double	getScatteringLength (const double lambda)
int	main (int argc, char **argv)

Variables
double	absorptionLengthFactor
	Scaling of absorption and scattering length. More...
double	scatteringLengthFactor

Detailed Description

Auxiliary program to draw PDF of Cherenkov light from muon.

Author

mdejong

Definition in file JDrawPDF.cc.

Function Documentation

double getAbsorptionLength ( const double  lambda )

inline

Definition at line 30 of file JDrawPDF.cc.

{
  return absorptionLengthFactor * NAMESPACE::getAbsorptionLength(lambda);
}

double getScatteringLength ( const double  lambda )

inline

Definition at line 36 of file JDrawPDF.cc.

{
  return scatteringLengthFactor * NAMESPACE::getScatteringLength(lambda);
}

int main	(	int	argc,
		char **	argv
	)

Definition at line 48 of file JDrawPDF.cc.

{
  using namespace std;
  using namespace JPP;

  typedef JAbstractHistogram<double> JHistogram_t;

  string        outputFile;
  int           numberOfPoints;
  double        epsilon;
  double        E;
  double        R;
  JAngle3D      dir;
  vector<int>   function;
  JHistogram_t  histogram;
  int           debug;

  try { 

    JProperties properties;

    properties.insert(gmake_property(JPHYSICS::MODULE_RADIUS_M));
  
    JParser<> zap("Auxiliary program to draw PDF of Cherenkov light from muon.");

    zap['@'] = make_field(properties)             = JPARSER::initialised();
    zap['o'] = make_field(outputFile)                                           = "";
    zap['n'] = make_field(numberOfPoints,         "points for integration")     = 25;
    zap['e'] = make_field(epsilon,                "precision for integration")  = 1.0e-10;
    zap['A'] = make_field(absorptionLengthFactor, "scaling factor")             = 1.0;
    zap['S'] = make_field(scatteringLengthFactor, "scaling factor")             = 1.0;
    zap['E'] = make_field(E,                      "muon energy [GeV]")          = 1.0;
    zap['R'] = make_field(R,                      "distance of approach [m]");
    zap['D'] = make_field(dir,                    "(theta, phi) of PMT [rad]");
    zap['F'] = make_field(function,               "PDF type");
    zap['H'] = make_field(histogram,              "histogram binning")          = JHistogram_t();
    zap['d'] = make_field(debug)                  = 0;

    zap(argc, argv);
  }
  catch(const exception &error) {
    FATAL(error.what() << endl);
  }
  

  const JPDF_C
    pdf(NAMESPACE::getPhotocathodeArea(),
        NAMESPACE::getQE,
        NAMESPACE::getAngularAcceptance,
        getAbsorptionLength,
        getScatteringLength,
        NAMESPACE::getScatteringProbability,
        NAMESPACE::getAmbientPressure(),
        getMinimalWavelength(),
        getMaximalWavelength(),
        numberOfPoints,
        epsilon);


  if (outputFile == "") {

    for (double dt; cin >> dt; ) {

      for (vector<int>::const_iterator F = function.begin(); F != function.end(); ++F) {

        cout << setw(2)         << *F             << ' '
             << SCIENTIFIC(7,1) << E              << ' '
             << FIXED(5,1)      << R              << ' '
             << FIXED(5,2)      << dir.getTheta() << ' '
             << FIXED(5,2)      << dir.getPhi()   << ' '
             << FIXED(5,1)      << dt             << ' '
             << SCIENTIFIC(9,3) << pdf.getLightFromMuon(*F, E, R, dir.getTheta(), dir.getPhi(), dt) * E << endl;
      }
    }

    return 0;
  }
  
 
  TFile out(outputFile.c_str(), "recreate");

  //const double z0   = -100.0;                                 // [m]
  //const double t0   = (-z0 +  R * getTanThetaC()) / C;        // time [ns]
  const double t0   =  0.0;                                   // time [ns]

  if (!histogram.is_valid()) {

    if (function.size() == 1 && function[0] == DIRECT_LIGHT_FROM_MUON) {

      histogram = JHistogram_t(t0 - 20.0, t0 + 50.0);

      histogram.setBinWidth(0.1);

    } else {

      histogram = JHistogram_t(t0 - 20.0, t0 + 500.0);

      histogram.setBinWidth(0.5);
    }
  }

  TH1D h0("h0", NULL, histogram.getNumberOfBins(), histogram.getLowerLimit(), histogram.getUpperLimit());

  JSplineFunction1S_t f1;
 
  for (int i = 1; i <= h0.GetNbinsX(); ++i) {

    const double dt = h0.GetBinCenter(i) - t0;

    double value = 0.0;

    for (vector<int>::const_iterator F = function.begin(); F != function.end(); ++F) {
      value += pdf.getLightFromMuon(*F, E, R, dir.getTheta(), dir.getPhi(), dt);
    }

    h0.SetBinContent(i, value);

    f1[dt] = value;
  }

  f1.setExceptionHandler(new JSplineFunction1S_t::JDefaultResult(JMATH::zero));
  f1.compile();

  const double T = 5;  // [ns]

  JQuantiles quantiles(f1);

  const double t1 = quantiles.getX();
  const double y  = f1(t1 + T).v - f1(t1 - T).v;

  DEBUG("E     " << E                       << endl);
  DEBUG("R     " << R                       << endl);
  DEBUG("theta " << dir.getTheta()          << endl);
  DEBUG("phi   " << dir.getPhi()            << endl);
  DEBUG("int   " << quantiles.getIntegral() << endl);
  DEBUG("t1    " << t1                      << endl);
  DEBUG("max   " << quantiles.getY()        << endl);
  DEBUG("FWHM  " << quantiles.getFWHM()     << endl);
  DEBUG("int[] " << y                       << endl);  

  
  out.Write();
  out.Close();
}

Variable Documentation

double absorptionLengthFactor

Scaling of absorption and scattering length.

Definition at line 27 of file JDrawPDF.cc.

double scatteringLengthFactor

Definition at line 28 of file JDrawPDF.cc.