Line data    Source code

       1              : #include "crpropa/module/PhotonOutput1D.h"
       2              : #include "crpropa/Units.h"
       3              : 
       4              : #include <iostream>
       5              : #include <sstream>
       6              : #include <cstdio>
       7              : #include <stdexcept>
       8              : 
       9              : #include "kiss/string.h"
      10              : #include "kiss/logger.h"
      11              : 
      12              : #ifdef CRPROPA_HAVE_ZLIB
      13              : #include <ozstream.hpp>
      14              : #endif
      15              : 
      16              : using namespace std;
      17              : 
      18              : namespace crpropa {
      19              : 
      20            0 : PhotonOutput1D::PhotonOutput1D() : out(&std::cout) {
      21            0 :         KISS_LOG_WARNING << "PhotonOutput1D is deprecated and will be removed in the future. Replace with TextOutput or HDF5Output with features ObserverNucleusVeto + ObserverDetectAll";
      22            0 : }
      23              : 
      24            0 : PhotonOutput1D::PhotonOutput1D(std::ostream &out) : out(&out) {
      25            0 :         KISS_LOG_WARNING << "PhotonOutput1D is deprecated and will be removed in the future. Replace with TextOutput or HDF5Output with features ObserverNucleusVeto + ObserverDetectAll";
      26            0 : }
      27              : 
      28            0 : PhotonOutput1D::PhotonOutput1D(const std::string &filename) : outfile(
      29            0 :         filename.c_str(), std::ios::binary), out(&outfile), filename(filename) {
      30            0 :         KISS_LOG_WARNING << "PhotonOutput1D is deprecated and will be removed in the future. Replace with TextOutput or HDF5Output with features ObserverNucleusVeto + ObserverDetectAll";
      31            0 :         if (kiss::ends_with(filename, ".gz"))
      32            0 :                 gzip();
      33              : 
      34            0 :         *out << "#ID\tE\tD\tpID\tpE\tiID\tiE\tiD\n";
      35            0 :         *out << "#\n";
      36            0 :         *out << "# ID          Id of particle (photon, electron, positron)\n";
      37            0 :         *out << "# E           Energy [EeV]\n";
      38            0 :         *out << "# D           Comoving distance to origin [Mpc]\n";
      39            0 :         *out << "# pID         Id of parent particle\n";
      40            0 :         *out << "# pE          Energy [EeV] of parent particle\n";
      41            0 :         *out << "# iID         Id of source particle\n";
      42            0 :         *out << "# iE          Energy [EeV] of source particle\n";
      43            0 :         *out << "# iD          Comoving distance [Mpc] to source\n";
      44            0 :         *out << "#\n";
      45            0 : }
      46              : 
      47            0 : void PhotonOutput1D::process(Candidate *candidate) const {
      48            0 :         int pid = candidate->current.getId();
      49            0 :         if ((pid != 22) and (abs(pid) != 11))
      50            0 :                 return;
      51              : 
      52              :         char buffer[1024];
      53              :         size_t p = 0;
      54              : 
      55            0 :         p += std::snprintf(buffer + p, 1024 - p, "%4i\t", pid);
      56            0 :         p += std::snprintf(buffer + p, 1024 - p, "%g\t", candidate->current.getEnergy() / EeV);
      57            0 :         p += std::snprintf(buffer + p, 1024 - p, "%8.4f\t", candidate->current.getPosition().getR() / Mpc);
      58              : 
      59            0 :         p += std::snprintf(buffer + p, 1024 - p, "%10i\t", candidate->created.getId());
      60            0 :         p += std::snprintf(buffer + p, 1024 - p, "%8.4f\t", candidate->created.getEnergy() / EeV);
      61              : 
      62            0 :         p += std::snprintf(buffer + p, 1024 - p, "%10i\t", candidate->source.getId());
      63            0 :         p += std::snprintf(buffer + p, 1024 - p, "%8.4f\t", candidate->source.getEnergy() / EeV);
      64            0 :         p += std::snprintf(buffer + p, 1024 - p, "%8.4f\n", candidate->source.getPosition().getR() / Mpc);
      65              : 
      66            0 : #pragma omp critical(FileOutput)
      67              :         {
      68            0 :                 out->write(buffer, p);
      69              :         }
      70              : 
      71            0 :         candidate->setActive(false);
      72              : }
      73              : 
      74            0 : void PhotonOutput1D::close() {
      75              :         #ifdef CRPROPA_HAVE_ZLIB
      76            0 :                 zstream::ogzstream *zs = dynamic_cast<zstream::ogzstream *>(out);
      77            0 :                 if (zs) {
      78            0 :                         zs->close();
      79            0 :                         delete out;
      80            0 :                         out = 0;
      81              :                 }
      82              :         #endif
      83            0 :         outfile.flush();
      84            0 : }
      85              : 
      86            0 : string PhotonOutput1D::getDescription() const {
      87            0 :         std::stringstream s;
      88              :         s << "PhotonOutput1D: Output file = " << filename;
      89            0 :         return s.str();
      90            0 : }
      91              : 
      92            0 : PhotonOutput1D::~PhotonOutput1D() {
      93            0 :         close();
      94            0 : }
      95              : 
      96            0 : void PhotonOutput1D::gzip() {
      97              :         #ifdef CRPROPA_HAVE_ZLIB
      98            0 :                 out = new zstream::ogzstream(*out);
      99              :         #else
     100              :                 throw std::runtime_error("CRPropa was build without Zlib compression!");
     101              :         #endif
     102            0 : }
     103              : 
     104              : } // namespace crpropa