Method findApexRT

src/topp/MetaProSIP.cpp:2475–2550

Extract the mono-isotopic trace and reports the rt of the maximum intensity Used to compensate for slight RT shifts (e.g. important if features of a different map are used) n_scans corresponds to the number of neighboring scan rts that should be extracted n_scan = 2 -> vector size = 1 + 2 + 2

  // n_scans corresponds to the number of neighboring scan rts that should be extracted
  // n_scan = 2 -> vector size = 1 + 2 + 2
  vector<double> findApexRT(const FeatureMap::iterator feature_it, double hit_rt, const PeakMap& peak_map, Size n_scans)
  {
    vector<double> seeds_rt;
    vector<Peak2D> mono_trace;

    if (!feature_it->getConvexHulls().empty())
    {
      // extract elution profile of 12C containing mass trace using a bounding box
      // first convex hull contains the monoisotopic 12C trace
      const DBoundingBox<2>& mono_bb = feature_it->getConvexHulls()[0].getBoundingBox();

      //(min_rt, max_rt, min_mz, max_mz)
      PeakMap::ConstAreaIterator ait = peak_map.areaBeginConst(mono_bb.minPosition()[0], mono_bb.maxPosition()[0], mono_bb.minPosition()[1], mono_bb.maxPosition()[1]);
      for (; ait != peak_map.areaEndConst(); ++ait)
      {
        Peak2D p2d;
        p2d.setRT(ait.getRT()); // get rt of scan
        p2d.setMZ(ait->getMZ()); // get peak 1D mz
        p2d.setIntensity(ait->getIntensity());
        mono_trace.push_back(p2d);
      }
    }

    // if there is no 12C mono trace generate a valid starting point
    if (mono_trace.empty())
    {
      Peak2D p2d;
      double next_valid_scan_rt = peak_map.RTBegin(hit_rt - 0.001)->getRT();
      p2d.setRT(next_valid_scan_rt);
      p2d.setMZ(0); // actually not needed
      p2d.setIntensity(0);
      mono_trace.push_back(p2d);
    }

    // determine trace peak with highest intensity
    double max_trace_int = -1e16;
    Size max_trace_int_idx = 0;

    for (Size j = 0; j != mono_trace.size(); ++j)
    {
      if (mono_trace[j].getIntensity() > max_trace_int)
      {
        max_trace_int = mono_trace[j].getIntensity();
        max_trace_int_idx = j;
      }
    }
    double max_trace_int_rt = mono_trace[max_trace_int_idx].getRT();
    seeds_rt.push_back(max_trace_int_rt);

    for (Size i = 1; i <= n_scans; ++i)
    {
      double rt_after = max_trace_int_rt;
      if (max_trace_int_idx < mono_trace.size() - (Int)i)
      {
        rt_after = mono_trace[max_trace_int_idx + i].getRT();
      }

      double rt_before = max_trace_int_rt;