/* This file is part of MAUS: http://micewww.pp.rl.ac.uk:8080/projects/maus
 *
 * MAUS is free software: you can redistribute it and/or modify
 * it under the terms of the GNU General Public License as published by
 * the Free Software Foundation, either version 3 of the License, or
 * (at your option) any later version.
 *
 * MAUS is distributed in the hope that it will be useful,
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 * GNU General Public License for more details.
 *
 * You should have received a copy of the GNU General Public License
 * along with MAUS.  If not, see <http://www.gnu.org/licenses/>.
 *
 */


#include "src/reduce/ReduceCppEVExport/EVHepRepExporter.hh"

#include <unistd.h>
#include <string>
#include <sstream>
#include <fstream>
#include <vector>

namespace EventViewer {

EVHepRepExporter::EVHepRepExporter(EVEvent evEvent) {
  event = evEvent;
}

EVHepRepExporter::~EVHepRepExporter() {
  delete mm; // is this really necessary?
  // Do nothing
}

void EVHepRepExporter::Export(int display) {

  HepRepXMLWriter* writer = new HepRepXMLWriter();
  std::stringstream ssOutFileName;

  // - retrieve destination directory
  char *destDir = getenv("EV_DEST_DIR");
  if (destDir != NULL)
    ssOutFileName << destDir << "/";
  else
    ssOutFileName << "./";

  // - open output heprep file
  ssOutFileName << "MICEEvent_" << event.runNumber << "_" << event.spillNumber
                << "_" << event.eventNumber << ".heprep";
  std::string dataFileName = "data.heprep";
  writer->open(dataFileName.c_str());

  // - define tof points
  writer->addType("TOFPoints", 0);
  writer->addInstance();
  writer->addAttValue("DrawAs", "Point");
  writer->addAttValue("MarkName", "Box");
  writer->addAttValue("MarkColor", 0., 1., 1.);
  for (int point = 0; point < 3; ++point) {
    double x = event.tofPoints[point][0];
    double y = event.tofPoints[point][1];
    double z = event.tofPoints[point][2];
    if (z != 5673390) {
      writer->addPrimitive();
      writer->addPoint(x, y, z);
    }
  }

  // - define track points
  writer->addType("TrackPoints", 0);
  writer->addInstance();
  writer->addAttValue("DrawAs", "Point");

  // - plot TKU track points in HepRApp
  writer->addType("TKUTrackPoints", 1);
  writer->addInstance();
  writer->addAttValue("MarkName", "Cross");
  writer->addAttValue("MarkColor", 1., 1., 0.);
  for (int point = 0; point < 5; ++point) {
    writer->addType("TKUTrackPoint", 2); // add layer of depth to apply attibutes
    writer->addInstance();

    double x = event.scifiUSTrackerPoints[point][0];
    double y = event.scifiUSTrackerPoints[point][1];
    double z = event.scifiUSTrackerPoints[point][2];

    // attributes - momentum components
    writer->addAttDef("Px", "Track momentum x component", "Physics", "GeV");
    writer->addAttDef("Py", "Track momentum y component", "Physics", "GeV");
    writer->addAttDef("Pz", "Track momentum z component", "Physics", "GeV");
    writer->addAttValue("Px", event.scifiUSTrackerPointsMomenta[point][0]);
    writer->addAttValue("Py", event.scifiUSTrackerPointsMomenta[point][1]);
    writer->addAttValue("Pz", event.scifiUSTrackerPointsMomenta[point][2]);

    if (z != 5673390) {
      writer->addPrimitive();
      writer->addPoint(x, y, z);
    }
  }

  // - plot TKD track points in HepRApp
  writer->addType("TKDTrackPoints", 1);
  writer->addInstance();
  writer->addAttValue("MarkName", "Cross");
  writer->addAttValue("MarkColor", 1., 1., 0.);
  for (int point = 0; point < 5; ++point) {
    writer->addType("TKDTrackPoint", 2); // add layer of depth to apply attibutes
    writer->addInstance();

    double x = event.scifiDSTrackerPoints[point][0];
    double y = event.scifiDSTrackerPoints[point][1];
    double z = event.scifiDSTrackerPoints[point][2];

    // attributes - momentum components
    writer->addAttDef("Px", "Track momentum x component", "Physics", "GeV");
    writer->addAttDef("Py", "Track momentum y component", "Physics", "GeV");
    writer->addAttDef("Pz", "Track momentum z component", "Physics", "GeV");
    writer->addAttValue("Px", event.scifiUSTrackerPointsMomenta[point][0]);
    writer->addAttValue("Py", event.scifiUSTrackerPointsMomenta[point][1]);
    writer->addAttValue("Pz", event.scifiUSTrackerPointsMomenta[point][2]);

    if (z != 5673390) {
      writer->addPrimitive();
      writer->addPoint(x, y, z);
    }
  }

  // - define space points
  writer->addType("SpacePoints", 0);
  writer->addInstance();
  writer->addAttValue("DrawAs", "Point");

  // - plot TKU space points in HepRApp
  writer->addType("TKUSpacePoints", 1);
  writer->addInstance();
  writer->addAttValue("MarkName", "Cross");
  writer->addAttValue("MarkColor", 1., 0.3, 0.);
  for (std::vector<MAUS::ThreeVector>::iterator it = event.scifiUSTrackerSpacePoints.begin(),
       end = event.scifiUSTrackerSpacePoints.end(); it != end; ++it) {
    writer->addType("TKUSpacePoint", 2); // add layer of depth to apply attibutes
    writer->addInstance();

    double x = it->x();
    double y = it->y();
    double z = it->z();

    writer->addAttDef("alpha P0", "Cluster info for plane 0 - alpha", "Physics", "channel");
    writer->addAttDef("z P0", "Cluster info for plane 0 - z", "Physics", "mm");
    writer->addAttDef("alpha P1", "Cluster info for plane 1", "Physics", "channel");
    writer->addAttDef("z P1", "Cluster info for plane 1 - z", "Physics", "mm");
    writer->addAttDef("alpha P2", "Cluster info for plane 2", "Physics", "channel");
    writer->addAttDef("z P2", "Cluster info for plane 2 - z", "Physics", "mm");
    unsigned int station = it - event.scifiUSTrackerSpacePoints.begin();
    if (station < event.scifiUSTrackerClusters[0][0].size()) {
      writer->addAttValue("alpha P0", event.scifiUSTrackerClusters[0][0][station]);
      writer->addAttValue("z P0", event.scifiUSTrackerClusters[0][1][station]);
    }
    if (station < event.scifiUSTrackerClusters[1][0].size()) {
      writer->addAttValue("alpha P1", event.scifiUSTrackerClusters[1][0][station]);
      writer->addAttValue("z P1", event.scifiUSTrackerClusters[1][1][station]);
    }
    if (station < event.scifiUSTrackerClusters[2][0].size()) {
      writer->addAttValue("alpha P2", event.scifiUSTrackerClusters[2][0][station]);
      writer->addAttValue("z P2", event.scifiUSTrackerClusters[2][1][station]);
    }

    if (z != 5673390) { // should be 5673390
      writer->addPrimitive();
      writer->addPoint(x, y, z);
    }
  }

  // - plot TKD space points in HepRApp
  writer->addType("TKDSpacePoints", 1);
  writer->addInstance();
  writer->addAttValue("MarkName", "Cross");
  writer->addAttValue("MarkColor", 1., 0.3, 0.);
  for (std::vector<MAUS::ThreeVector>::iterator it = event.scifiDSTrackerSpacePoints.begin(),
       end = event.scifiDSTrackerSpacePoints.end(); it != end; ++it) {
    writer->addType("TKDSpacePoint", 2); // add layer of depth to apply attibutes
    writer->addInstance();

    double x = it->x();
    double y = it->y();
    double z = it->z();

    writer->addAttDef("alpha P0", "Cluster info for plane 0 - alpha", "Physics", "channel");
    writer->addAttDef("z P0", "Cluster info for plane 0 - z", "Physics", "mm");
    writer->addAttDef("alpha P1", "Cluster info for plane 1", "Physics", "channel");
    writer->addAttDef("z P1", "Cluster info for plane 1 - z", "Physics", "mm");
    writer->addAttDef("alpha P2", "Cluster info for plane 2", "Physics", "channel");
    writer->addAttDef("z P2", "Cluster info for plane 2 - z", "Physics", "mm");
    unsigned int station = it - event.scifiDSTrackerSpacePoints.begin();
    if (station < event.scifiDSTrackerClusters[0][0].size()) {
      writer->addAttValue("alpha P0", event.scifiDSTrackerClusters[0][0][station]);
      writer->addAttValue("z P0", event.scifiDSTrackerClusters[0][1][station]);
    }
    if (station < event.scifiDSTrackerClusters[1][0].size()) {
      writer->addAttValue("alpha P1", event.scifiDSTrackerClusters[1][0][station]);
      writer->addAttValue("z P1", event.scifiDSTrackerClusters[1][1][station]);
    }
    if (station < event.scifiDSTrackerClusters[2][0].size()) {
      writer->addAttValue("alpha P2", event.scifiDSTrackerClusters[2][0][station]);
      writer->addAttValue("z P2", event.scifiDSTrackerClusters[2][1][station]);
    }

    if (z != 5673390) { // should be 5673390
      writer->addPrimitive();
      writer->addPoint(x, y, z);
    }
  }

  // - define straight track
  writer->addType("StraightTrack", 0);
  writer->addInstance();
  writer->addAttValue("DrawAs", "Line");
  writer->addAttValue("LineColor", 1., 0., 0.);

  // - plot TKU straight track(s)
  writer->addType("TKUStraightTrack", 1);
  writer->addInstance();
  writer->addPrimitive();
  for (int point = 0; point < 2; ++point) {
    double x = event.scifiUSTrackerStraightTrackPoints[point][0];
    double y = event.scifiUSTrackerStraightTrackPoints[point][1];
    double z = event.scifiUSTrackerStraightTrackPoints[point][2];
    if (z != 5673390) { // ToDo: define dummy value at one place and use through out
      writer->addPoint(x, y, z);
    }
  }


  // - load parent geometry mice module (needed for helical track transformations)
  std::stringstream modulePath;
  if (getenv("EV_GEOM_HOME"))
      modulePath << getenv("EV_GEOM_HOME") << "/ParentGeometryFile.dat";
  else
      modulePath << getenv("MAUS_ROOT_DIR") << "/files/geometry/download/ParentGeometryFile.dat";

  // - Adam's fix (redirecting stdout/stderr)
  std::ostream* stdout = new std::ofstream();
  stdout->rdbuf(std::cout.rdbuf());
  std::ostream* stderr = new std::ofstream();
  stderr->rdbuf(std::cerr.rdbuf());

  // - leaving this commented out code in case OnRec crashes debugging is resumed
  /*try {
      std::cout << "Run: " << event.runNumber << "  Spill: " << event.spillNumber
                << "  Event: " << event.eventNumber << " Try MM" << std::endl;
      mm = new MiceModule(modulePath.str());
      std::cout << "Just MM init" << std::endl;
      std::cout.rdbuf(stdout->rdbuf());
      std::cerr.rdbuf(stderr->rdbuf());
      std::cout << "Just before catch" << std::endl;
  } catch (MAUS::Exceptions::Exception e) {
      std::cerr << "[EVHepRepExporter]Caught exception while building MiceModules! " 
                << "Possibly ParentGeometryFile.dat missing! Using default geometry!" << std::endl;
      e.what();
      e.Print();
      modulePath.str("");
      modulePath << getenv("MAUS_ROOT_DIR") << "/src/utilities/event-viewer/ParentGeometryFile.dat";
      try {
          mm = new MiceModule(modulePath.str());
          std::cout.rdbuf(stdout->rdbuf());
          std::cerr.rdbuf(stderr->rdbuf());
      } catch (MAUS::Exceptions::Exception e) {
          std::cerr << "[EVHepRepExporter]Caught exception while building MiceModules! " 
                << "Possibly ParentGeometryFile.dat missing! Using default geometry!" << std::endl;
          e.what();
          e.Print();
          return;
      }
  }
  
  if (mm) {// substitute this with simple check
    std::cout << "MM path: " << modulePath.str() << "  MM name: " << mm->name() << std::endl;
  } else {
    std::cerr << "Bad mm pointer: " << mm << std::endl;
    return;
  }*/

  try {
      mm = new MiceModule(modulePath.str());
      std::cout.rdbuf(stdout->rdbuf());
      std::cerr.rdbuf(stderr->rdbuf());
  } catch (MAUS::Exceptions::Exception e) {
      std::cerr << "Caught exception while building MiceModules" << std::endl;
      e.what();
      e.Print();
      modulePath.str("");
      // this file should be included in the release so no try/catch neccessary
      // (file in this dir not referenced? remove?)
      modulePath << getenv("MAUS_ROOT_DIR") << "/src/utilities/event-viewer/ParentGeometryFile.dat";
      mm = new MiceModule(modulePath.str());
      std::cout.rdbuf(stdout->rdbuf());
      std::cerr.rdbuf(stderr->rdbuf());
  }/**/


  // --- fetch US tracker global geometry
  mmTrackerUS = mm->findModulesByName("ParentGeometryFile.dat/Tracker0.dat");

  // - get TrackerRef0 position
  double trackerRef0posX = mmTrackerUS[0]->daughter(0)->globalPosition().x();
  double trackerRef0posY = mmTrackerUS[0]->daughter(0)->globalPosition().y();
  double trackerRef0posZ = mmTrackerUS[0]->daughter(0)->globalPosition().z();

  // - get rotation matrix
  double xx = mmTrackerUS[0]->globalRotation().xx();
  double xy = mmTrackerUS[0]->globalRotation().xy();
  double xz = mmTrackerUS[0]->globalRotation().xz();
  double yx = mmTrackerUS[0]->globalRotation().yx();
  double yy = mmTrackerUS[0]->globalRotation().yy();
  double yz = mmTrackerUS[0]->globalRotation().yz();
  double zx = mmTrackerUS[0]->globalRotation().zx();
  double zy = mmTrackerUS[0]->globalRotation().zy();
  double zz = mmTrackerUS[0]->globalRotation().zz();

  // - fill matrix elements (for VectorTransform)
  std::vector<double> matrixElements;
  matrixElements.push_back(xx);
  matrixElements.push_back(xy);
  matrixElements.push_back(xz);
  matrixElements.push_back(yx);
  matrixElements.push_back(yy);
  matrixElements.push_back(yz);
  matrixElements.push_back(zx);
  matrixElements.push_back(zy);
  matrixElements.push_back(zz);/**/

  // - dimensions (y coordinate is length)
  Hep3Vector dimensions = mmTrackerUS[0]->dimensions();
  double trackerLength = dimensions.y();

  // - plot TKU helical track(s)
  writer->addType("TKUHelicalTrackTest", 1);
  writer->addInstance();
  writer->addPrimitive();
  double circle_x0 = event.scifiUSTrackerHelicalTrackParameters[0];
  double circle_y0 = event.scifiUSTrackerHelicalTrackParameters[1];
  double radius = event.scifiUSTrackerHelicalTrackParameters[2];
  double dsdz = event.scifiUSTrackerHelicalTrackParameters[3];
  double line_dz_c = event.scifiUSTrackerHelicalTrackParameters[4];
  double posZ = event.scifiUSTrackerHelicalTrackParameters[7];

  double zLocal = posZ;
  int zStep = trackerLength/300;
  while (zLocal < trackerLength) {
      zLocal += zStep; // z value in tracker local reference system
      // x value in tracker local reference system
      double upstream_helix_x_value
              = circle_x0 + (radius*cos((1/radius)*(dsdz*zLocal+line_dz_c)));
      // y value in tracker local reference system
      double upstream_helix_y_value
              = circle_y0 + (radius*sin((1/radius)*(dsdz*zLocal+line_dz_c)));

      // transformation to global coordinate system
      std::vector<double> vSpacePoint
              = {upstream_helix_x_value, upstream_helix_y_value, zLocal};
      VectorTransform trans(vSpacePoint);
      trans.Rotate(matrixElements, "reg");
      trans.Translate(trackerRef0posX, trackerRef0posY, trackerRef0posZ, "reg");

      double x = trans.GetX();
      double y = trans.GetY();
      double z = trans.GetZ();

      writer->addPoint(x, y, z);
  }


  // - plot TKD straight track(s)
  writer->addType("TKDStraightTrack", 1);
  writer->addInstance();
  writer->addPrimitive();
  for (int point = 0; point < 2; ++point) {
    double x = event.scifiDSTrackerStraightTrackPoints[point][0];
    double y = event.scifiDSTrackerStraightTrackPoints[point][1];
    double z = event.scifiDSTrackerStraightTrackPoints[point][2];
    if (z != 5673390) {
      writer->addPoint(x, y, z);
    }
  }


  // --- fetch DS tracker global geometry
  mmTrackerDS = mm->findModulesByName("ParentGeometryFile.dat/Tracker1.dat");

  // - get TrackerRef1 position
  double trackerRef1posX = mmTrackerDS[0]->daughter(0)->globalPosition().x();
  double trackerRef1posY = mmTrackerDS[0]->daughter(0)->globalPosition().y();
  double trackerRef1posZ = mmTrackerDS[0]->daughter(0)->globalPosition().z();

  // - get rotation matrix
  xx = mmTrackerDS[0]->globalRotation().xx();
  xy = mmTrackerDS[0]->globalRotation().xy();
  xz = mmTrackerDS[0]->globalRotation().xz();
  yx = mmTrackerDS[0]->globalRotation().yx();
  yy = mmTrackerDS[0]->globalRotation().yy();
  yz = mmTrackerDS[0]->globalRotation().yz();
  zx = mmTrackerDS[0]->globalRotation().zx();
  zy = mmTrackerDS[0]->globalRotation().zy();
  zz = mmTrackerDS[0]->globalRotation().zz();

  // - fill matrix elements (for VectorTransform)
  matrixElements.clear();
  matrixElements.push_back(xx);
  matrixElements.push_back(xy);
  matrixElements.push_back(xz);
  matrixElements.push_back(yx);
  matrixElements.push_back(yy);
  matrixElements.push_back(yz);
  matrixElements.push_back(zx);
  matrixElements.push_back(zy);
  matrixElements.push_back(zz);


  // - plot TKD helical track(s)
  writer->addType("TKDHelicalTrackTest", 1);
  writer->addInstance();
  writer->addPrimitive();
  circle_x0 = event.scifiDSTrackerHelicalTrackParameters[0];
  circle_y0 = event.scifiDSTrackerHelicalTrackParameters[1];
  radius = event.scifiDSTrackerHelicalTrackParameters[2];
  dsdz = event.scifiDSTrackerHelicalTrackParameters[3];
  line_dz_c = event.scifiDSTrackerHelicalTrackParameters[4];
  posZ = event.scifiDSTrackerHelicalTrackParameters[7];
  zLocal = posZ;

  while (zLocal < trackerLength) {
      zLocal += zStep; // z value in tracker local reference system
      // x value in tracker local reference system
      double downstream_helix_x_value
              = circle_x0 + (radius*cos((1/radius)*(dsdz*zLocal+line_dz_c)));
      // y value in tracker local reference system
      double downstream_helix_y_value
              = circle_y0 + (radius*sin((1/radius)*(dsdz*zLocal+line_dz_c)));

      // transformation to global coordinate system
      std::vector<double> vSpacePoint
              = {downstream_helix_x_value, downstream_helix_y_value, zLocal};
      VectorTransform trans(vSpacePoint);
      trans.Rotate(matrixElements, "reg");
      trans.Translate(trackerRef1posX, trackerRef1posY, trackerRef1posZ, "reg");

      double x = trans.GetX();
      double y = trans.GetY();
      double z = trans.GetZ();

      writer->addPoint(x, y, z);
  }

  writer->close();

  char *detGeomFileName = getenv("EV_DETGEOM_FILE");
  // disable writing to HepRep options if geometry file is not selected
  if (detGeomFileName == NULL) {
    std::cerr << "WARNING: Detector geometry file not defined! Modify and source configure.sh or "
              << "choose geometry file using selection button (if using GUI). "
              << "Not writting output heprep file!" << std::endl;
  } else {
    std::string sOutFileName = ssOutFileName.str();
    std::string sDetGeomFileName(detGeomFileName);
    Concatenate(sDetGeomFileName, dataFileName, sOutFileName);
  }

  // - remove MAUS data heprep file
  std::stringstream cmdRemove;
  cmdRemove << "rm " << dataFileName;
  system(cmdRemove.str().c_str());

  // - open the file for viewing
  if (display) {
    std::cout << display << std::endl;
    std::stringstream cmdDisplay;
    int displayCheck = 1;

    char *javaDir = getenv("EV_JAVA_DIR");
    if (javaDir != NULL) {
      cmdDisplay << javaDir << "/bin/java -jar ";
    } else {
      std::cerr << "WARNING: No java directory selected! HepRApp display unavailable! "
                << "Select directory in which your JRE can be found!" << std::endl;
      displayCheck = 0;
    }

    char *heprappDir = getenv("EV_HEPRAPP_DIR");
    if (heprappDir != NULL) {
      cmdDisplay << heprappDir << "/HepRApp.jar -file ";
    } else {
      std::cerr << "WARNING: No HepRapp directory selected! HepRApp display unavailable! "
                << "Select directory in which your HepRApp.jar can be found!" << std::endl;
      displayCheck = 0;
    }

    cmdDisplay << ssOutFileName.str();
    std::cout << cmdDisplay.str() << std::endl;
    if (displayCheck)
      system(cmdDisplay.str().c_str());
  }
}

int EVHepRepExporter::Concatenate(std::string& geometryFileName, std::string& dataFileName,
                                  std::string& outFileName) {
  std::ifstream geomFile;
  geomFile.open(geometryFileName.c_str(), std::ifstream::in);
  if (!geomFile.is_open()) {
    std::cerr << "ERROR: Could not open geometry file needed for output to heprep!" << std::endl;
    return 0;
  }

  std::ifstream dataFile;
  dataFile.open(dataFileName.c_str(), std::fstream::in);
  if (!dataFile.is_open()) {
    std::cerr << "ERROR: Could not open data file needed for output to heprep!" << std::endl;
    return 0;
  }

  std::ofstream outFile;
  outFile.open(outFileName.c_str(), std::ifstream::out);
  if (!outFile.is_open()) {
    std::cerr << "ERROR: Could not open temporary file needed for output to heprep!" << std::endl;
    return 0;
  }


  while (!geomFile.eof()) {
    std::string line, target;
    std::getline(geomFile, line);
    target = line.substr(0, 16);
    if (target != "</heprep:heprep>")
      outFile << line << std::endl;
  }

  int counter = 0;
  while (!dataFile.eof()) {
    ++counter;
    std::string line;
    std::getline(dataFile, line);
    if (counter > 3)
      outFile << line << std::endl;
  }

  return 1;
}

}  // ~namespace EventViewer