LoadEventPreNexus2.cpp

        local_num_error_events++;
        local_num_bad_events ++;
        continue;
      }

      //Covert the pixel ID from DAS pixel to our pixel ID
      // downstream monitor pixel for SNAP
      if(pid ==1073741843) pid = 1179648;
      else if (this->using_mapping_file)
      {
        PixelType unmapped_pid = pid % this->numpixel;
        pid = this->pixelmap[unmapped_pid];
      }

      // Wrong pixel IDs
      if (pid > static_cast<PixelType>(detid_max))
      {
        iswrongdetid = true;

        local_num_error_events++;
        local_num_wrongdetid_events++;
        local_wrongdetids.insert(pid);
      }

      // Now check if this pid we want to load.
      if (loadOnlySomeSpectra && !iswrongdetid)
      {
        std::map<int64_t, bool>::iterator it;
        it = spectraLoadMap.find(pid);
        if (it == spectraLoadMap.end())
        {
          //Pixel ID was not found, so the event is being ignored.
          local_num_ignored_events++;
          continue;
        }
      }

      // Upon this point, only 'good' events are left to work on

      // Pulse: Find the pulse time for this event index
      if (pulse_i < numPulses-1)
      {
        //This is the total offset into the file
        size_t total_i = i + fileOffset;
        // Go through event_index until you find where the index increases to encompass the current index.
        // Your pulse = the one before.
        while (!((total_i >= event_indices[pulse_i]) && (total_i < event_indices[pulse_i+1])) )
        {
          pulse_i++;
          if (pulse_i >= (numPulses-1))
            break;
        }

        //Save the pulse time at this index for creating those events
        pulsetime = pulsetimes[pulse_i];
      } // Find pulse time

      // TOF
      double tof = static_cast<double>(temp.tof) * TOF_CONVERSION;

      if (!iswrongdetid)
      {
        // Regular event that is belonged to a defined detector
        // Find the overall max/min tof
        if (tof < local_shortest_tof)
          local_shortest_tof = tof;
        if (tof > local_longest_tof)
          local_longest_tof = tof;

        // This is equivalent to workspace->getEventList(this->pixel_to_wkspindex[pid]).addEventQuickly(event);
        // But should be faster as a bunch of these calls were cached.
#if defined(__GNUC__) && !(defined(__INTEL_COMPILER)) && !(defined(__clang__))
        // This avoids a copy constructor call but is only available with GCC (requires variadic templates)
        arrayOfVectors[pid]->emplace_back( tof, pulsetime );
#else
        arrayOfVectors[pid]->push_back(TofEvent(tof, pulsetime));
#endif

        ++ local_num_good_events;
      }
      else
      {
        // Special events/Wrong detector id
        // i.  get/add index of the entry in map
        std::map<PixelType, size_t>::iterator it;
        it = local_pidindexmap.find(pid);
        size_t theindex = 0;
        if (it == local_pidindexmap.end())
        {
          // Initialize it!
          size_t newindex = local_pulsetimes.size();
          local_pidindexmap[pid] = newindex;

          std::vector<Kernel::DateAndTime> tempvectime;
          std::vector<double> temptofs;
          local_pulsetimes.push_back(tempvectime);
          local_tofs.push_back(temptofs);

          theindex = newindex;

          ++ numwrongpid;

          g_log.debug() << "Find New Wrong Pixel ID = " << pid << "\n";
        }
        else
        {
          // existing
          theindex = it->second;
        }

        // ii. calculate and add absolute time
        // int64_t abstime = (pulsetime.totalNanoseconds()+int64_t(tof*1000));
        local_pulsetimes[theindex].push_back(pulsetime);
        local_tofs[theindex].push_back(tof);

      } // END-IF-ELSE: On Event's Pixel's Nature

    } // ENDFOR each event

    // Update local statistics to their global counterparts
    PARALLEL_CRITICAL( LoadEventPreNexus2_global_statistics )
    {
      this->num_good_events += local_num_good_events;
      this->num_ignored_events += local_num_ignored_events;
      this->num_error_events += local_num_error_events;

      this->num_bad_events += local_num_bad_events;
      this->num_wrongdetid_events += local_num_wrongdetid_events;

      std::set<PixelType>::iterator it;
      for (it = local_wrongdetids.begin(); it != local_wrongdetids.end(); ++it)
      {
        PixelType tmpid = *it;
        this->wrongdetids.insert(*it);

        // Create class map entry if not there
        size_t mindex = 0;
        std::map<PixelType, size_t>::iterator git = this->wrongdetidmap.find(tmpid);
        if (git == this->wrongdetidmap.end())
        {
          // create entry
          size_t newindex = this->wrongdetid_pulsetimes.size();
          this->wrongdetidmap[tmpid] = newindex;

          std::vector<Kernel::DateAndTime> temppulsetimes;
          std::vector<double> temptofs;
          this->wrongdetid_pulsetimes.push_back(temppulsetimes);
          this->wrongdetid_tofs.push_back(temptofs);

          mindex = newindex;
        } else {
          mindex = git->second;
        }

        // 2. Find local
        std::map<PixelType, size_t>::iterator lit = local_pidindexmap.find(tmpid);
        size_t localindex= lit->second;

        for (size_t iv = 0; iv < local_pulsetimes[localindex].size(); iv ++)
        {
          this->wrongdetid_pulsetimes[mindex].push_back(local_pulsetimes[localindex][iv]);
          this->wrongdetid_tofs[mindex].push_back(local_tofs[localindex][iv]);
        }
        // std::sort(this->wrongdetid_abstimes[mindex].begin(), this->wrongdetid_abstimes[mindex].end());

      }


      if (local_shortest_tof < shortest_tof)
        shortest_tof = local_shortest_tof;
      if (local_longest_tof > longest_tof)
        longest_tof = local_longest_tof;
    } // END_CRITICAL

    return;
  }

  //----------------------------------------------------------------------------------------------
  /** Comparator for sorting dasevent lists
    */
  bool vzintermediatePixelIDComp(IntermediateEvent x, IntermediateEvent y)
  {
    return (x.pid < y.pid);
  }

  //-----------------------------------------------------------------------------
  /**
    * Add a sample environment log for the proton chage (charge of the pulse in picoCoulombs)
    * and set the scalar value (total proton charge, microAmps*hours, on the sample)
    *
    * @param workspace :: Event workspace to set the proton charge on
    */
  void LoadEventPreNexus2::setProtonCharge(DataObjects::EventWorkspace_sptr & workspace)
  {
    if (this->proton_charge.empty()) // nothing to do
      return;

    Run& run = workspace->mutableRun();

    //Add the proton charge entries.
    TimeSeriesProperty<double>* log = new TimeSeriesProperty<double>("proton_charge");
    log->setUnits("picoCoulombs");

    //Add the time and associated charge to the log
    log->addValues(this->pulsetimes, this->proton_charge);


    /// TODO set the units for the log
    run.addLogData(log);
    double integ = run.integrateProtonCharge();

    g_log.information() << "Total proton charge of " << integ << " microAmp*hours found by integrating.\n";

    return;
  }

  //-----------------------------------------------------------------------------
  /** Load a pixel mapping file
    * @param filename :: Path to file.
    */
  void LoadEventPreNexus2::loadPixelMap(const std::string &filename)
  {
    this->using_mapping_file = false;
    this->pixelmap.clear();

    // check that there is a mapping file
    if (filename.empty()) {
      this->g_log.information("NOT using a mapping file");
      return;
    }

    // actually deal with the file
    this->g_log.debug("Using mapping file \"" + filename + "\"");

    //Open the file; will throw if there is any problem
    BinaryFile<PixelType> pixelmapFile(filename);
    PixelType max_pid = static_cast<PixelType>(pixelmapFile.getNumElements());
    //Load all the data
    pixelmapFile.loadAllInto( this->pixelmap );

    //Check for funky file
    if (std::find_if(pixelmap.begin(), pixelmap.end(), std::bind2nd(std::greater<PixelType>(), max_pid))
            != pixelmap.end())
    {
      this->g_log.warning("Pixel id in mapping file was out of bounds. Loading without mapping file");
      this->numpixel = 0;
      this->pixelmap.clear();
      this->using_mapping_file = false;
      return;
    }

    //If we got here, the mapping file was loaded correctly and we'll use it
    this->using_mapping_file = true;
    //Let's assume that the # of pixels in the instrument matches the mapping file length.
    this->numpixel = static_cast<uint32_t>(pixelmapFile.getNumElements());
  }



  //-----------------------------------------------------------------------------
  /** Open an event file
    * @param filename :: file to open.
    */
  void LoadEventPreNexus2::openEventFile(const std::string &filename)
  {
    //Open the file
    eventfile = new BinaryFile<DasEvent>(filename);
    num_events = eventfile->getNumElements();
    g_log.debug() << "File contains " << num_events << " event records.\n";

    // Check if we are only loading part of the event file
    const int chunk = getProperty("ChunkNumber");
    if ( isEmpty(chunk) ) // We are loading the whole file
    {
      first_event = 0;
      max_events = num_events;
    }
    else // We are loading part - work out the event number range
    {
      const int totalChunks = getProperty("TotalChunks");
      max_events = num_events/totalChunks;
      first_event = (chunk - 1) * max_events;
      // Need to add any remainder to the final chunk
      if ( chunk == totalChunks ) max_events += num_events%totalChunks;
    }

    g_log.information()<< "Reading " <<  max_events << " event records\n";

    return;
  }


  //-----------------------------------------------------------------------------
  /** Read a pulse ID file
   * @param filename :: file to load.
   * @param throwError :: Flag to trigger error throwing instead of just logging
   */
  void LoadEventPreNexus2::readPulseidFile(const std::string &filename, const bool throwError)
  {
    this->proton_charge_tot = 0.;
    this->num_pulses = 0;
    this->pulsetimesincreasing = true;

    // jump out early if there isn't a filename
    if (filename.empty()) {
      this->g_log.information("NOT using a pulseid file");
      return;
    }

    std::vector<Pulse> * pulses;

    // set up for reading
    //Open the file; will throw if there is any problem
    try {
      BinaryFile<Pulse> pulseFile(filename);

      //Get the # of pulse
      this->num_pulses = pulseFile.getNumElements();
      this->g_log.information() << "Using pulseid file \"" << filename << "\", with " << num_pulses
          << " pulses.\n";

      //Load all the data
      pulses = pulseFile.loadAll();
    } catch (runtime_error &e) {
      if (throwError)
      {
        throw;
      }
      else
      {
        this->g_log.information() << "Encountered error in pulseidfile (ignoring file): " << e.what() << "\n";
        return;
      }
    }

    if (num_pulses > 0)
    {
      double temp;
      DateAndTime lastPulseDateTime(0, 0);
      this->pulsetimes.reserve(num_pulses);
      for (size_t i=0; i < num_pulses; i++)
      {
        Pulse & it = (*pulses)[i];
        DateAndTime pulseDateTime( (int64_t) it.seconds, (int64_t) it.nanoseconds);
        this->pulsetimes.push_back(pulseDateTime);
        this->event_indices.push_back(it.event_index);

        if (pulseDateTime < lastPulseDateTime)
          this->pulsetimesincreasing = false;
        else
          lastPulseDateTime = pulseDateTime;

        temp = it.pCurrent;
        this->proton_charge.push_back(temp);
        if (temp < 0.)
          this->g_log.warning("Individual proton charge < 0 being ignored");
        else
          this->proton_charge_tot += temp;
      }
    }

    this->proton_charge_tot = this->proton_charge_tot * CURRENT_CONVERSION;

    //Clear the vector
    delete pulses;

  }

} // namespace DataHandling
} // namespace Mantid