#ifndef __JEVENT__
#define __JEVENT__

#include <limits>

#include "antcc/Event.hh" 

#include "JDetector/JDetector.hh"
#include "JDetector/JDetectorToolkit.hh"
#include "JDetector/JK40Simulator.hh"
#include "JDetector/JTimeRange.hh"

#include "JSirene/JAntcc.hh"


namespace JSIRENE {


  /**
   * Auxiliary class to allow modifications of the hit list of a Monte Carlo event.
   */
  class JEvent : 
    public Event
  {
  public:
    /**
     * Default constructor.
     */
    JEvent() :
      Event()
    {}


    /**
     * Copy constructor.
     *
     * \param  event            Monte Carlo event
     */
    JEvent(const Event& event) :
      Event(event)
    {}

    
    /**
     * Remove Hits.
     */
    void removeHits()
    {
      HitList_.clear();
    }


    /**
     * Add Hit.
     *
     * \param  hit              hit
     */
    void Add(const Hit& hit)
    {
      HitList_.push_back(hit);
    }

    
    /**
     * Add Secondary, mainly used for EM showers
     *
     * \param secondary         secondary track
     */
    void Add(const Secondary& secondary)
    {
      secondary_.push_back(secondary);
    }

    
    /**
     * Add random noise.
     * This method removes existing noise hits before adding new noise hits.
     * The PMT simulation is not applied. 
     * N.B. The time window extends the time range of hits in the event, if any.
     *
     * \param  detector         detector
     * \param  k40Simulator     K40 simulator
     * \param  period           time window [ns]
     */
    void addNoise(const JDetector&     detector, 
		  const JK40Simulator& k40Simulator,
		  const JTimeRange&    period)
    {
      using namespace std;

      JTimeRange time_range(JTimeRange::JDEFAULT_RANGE);

      vector<Hit>::iterator __end = HitList_.end();

      for (vector<Hit>::iterator i = HitList_.begin(); i != __end; ) {

        if (is_noise(*i)) {

	  iter_swap(i, --__end);

	} else {

	  time_range.include(getTime(*i));

	  ++i;
        }
      }

      HitList_.erase(__end, HitList_.end());

      if (time_range.is_valid())
	time_range.add(period);
      else
	time_range = period;


      const double npe = 1.0;

      JModuleData buffer;

      for (JDetector::const_iterator module = detector.begin(); module != detector.end(); ++module) {

	buffer.reset(module->size());

	k40Simulator.generateHits(*module, time_range, buffer);

	for (unsigned int pmt = 0; pmt != buffer.size(); ++pmt) {

	  const JModuleData::value_type& frame = buffer[pmt];
	  
	  for (JModuleData::value_type::const_iterator hit = frame.begin(); hit != frame.end(); ++hit)
	    HitList_.push_back(Hit(HitList_.size() + 1, module->getPMT(pmt).getID(), hit->t_ns, npe, HIT_TYPE_NOISE, 0, npe, hit->t_ns));
	}
      }
    }


    /**
     * Less than operator for Monte Carlo hits.
     *
     * First hit is defined as:
     *    1/ smallest PMT identifier 
     *    2/ earliest time if same PMT identifier
     *
     * \param  first            first  hit
     * \param  second           second hit
     * \return                  true if first hit earlier than second hit; else false 
     */
    static bool compare(const Hit& first, const Hit& second)
    {
      if (first.pm_id() == second.pm_id())
	return first.pure_dt() < second.pure_dt();
      else
	return first.pm_id()   < second.pm_id();
    }

    
    /**
     * Merge Hits on same PMT that are close in time.
     *
     * The first hit has sum of amplitudes all hits within Tmax.
     * The obsolete hits are subsequently removed.
     * The hit identifiers are reset in ascending order, starting at one.
     *
     * \param  Tmax             maximal time difference [ns]
     */
    void mergeHits(const double Tmax)
    {
      if (Tmax > 0.0 && !HitList().empty()) {

	typedef std::vector<Hit>::iterator  iterator;

	iterator __begin = HitList_.begin();
	iterator __end   = HitList_.end();

	std::sort(__begin, __end, compare);

	iterator out = __begin;

	out->set_ident(1);                // set first identifier

	while (++__begin != __end) {

	  if (out->pm_id() == __begin->pm_id() && getTime(*__begin) - getTime(*out) <= Tmax) {

	    out->set_pure_npe(out->pure_npe() + __begin->pure_npe());
	    out->set_npe     (out->npe()      + __begin->npe());

	  } else {
	    
	    unsigned int ID = out->ident();
	    
	    ++out;                        // increment desitination address
	    
	    *out = *__begin;              // copy first hit
	  
	    out->set_ident(++ID);         // set identifier
	  }
	}

	HitList_.erase(++out, __end);     // remove obsolete hits
      }
    }
  };
}

#endif