Browse code

feat: pwiz returns peak matrices with column names

- `peaks` method for the proteowizard backend also peaks matrices with colnames
set.

jorainer authored on 18/03/2022 10:38:45
Showing4 changed files

... ...
@@ -2,7 +2,7 @@ Package: mzR
2 2
 Type: Package
3 3
 Title: parser for netCDF, mzXML, mzData and mzML and mzIdentML files
4 4
        (mass spectrometry data)
5
-Version: 2.27.1
5
+Version: 2.27.2
6 6
 Author: Bernd Fischer, Steffen Neumann, Laurent Gatto, Qiang Kou, Johannes Rainer
7 7
 Maintainer: Steffen Neumann <sneumann@ipb-halle.de>,
8 8
             Laurent Gatto <laurent.gatto@uclouvain.be>,
... ...
@@ -1,3 +1,7 @@
1
+CHANGES IN VERSION 2.27.2
2
+-------------------------
3
+ o `peaks` method for pwiz backend sets colnames on returned matrices.
4
+
1 5
 CHANGES IN VERSION 2.27.1
2 6
 -------------------------
3 7
  o Add missing atomic_count_sync.hpp BH file (see PR #248 by vjcitn)
... ...
@@ -11,6 +11,7 @@ test_mzXML <- function() {
11 11
     checkTrue(is.matrix(peaks(mzxml,1)))
12 12
     pks <- peaks(mzxml, 2:3)
13 13
     checkTrue(length(pks) == 2)
14
+    checkEquals(colnames(pks[[1L]]), c("mz", "intensity"))
14 15
     peaksCount(mzxml)
15 16
     hdr <- header(mzxml)
16 17
     checkTrue(any(colnames(hdr) == "spectrumId"))
... ...
@@ -39,6 +40,7 @@ test_mzML <- function() {
39 40
     pks <- peaks(mzml)
40 41
     pks <- peaks(mzml,1)
41 42
     pks <- peaks(mzml,2:3)
43
+    checkEquals(colnames(pks[[1L]]), c("mz", "intensity"))
42 44
     peaksCount(mzml)
43 45
     hdr <- header(mzml)
44 46
     checkTrue(any(colnames(hdr) == "spectrumId"))
... ...
@@ -374,6 +374,7 @@ Rcpp::List RcppPwiz::getPeakList(Rcpp::IntegerVector whichScan) {
374 374
     std::vector<double> data;
375 375
     Rcpp::NumericVector data_matrix;
376 376
     Rcpp::List res(n_want);
377
+    Rcpp::CharacterVector cn = Rcpp::CharacterVector::create("mz", "intensity");
377 378
     for (int i = 0; i < n_want; i++) {
378 379
       current_scan = whichScan[i];
379 380
       if (current_scan < 1 || current_scan > n_scans) {
... ...
@@ -385,6 +386,7 @@ Rcpp::List RcppPwiz::getPeakList(Rcpp::IntegerVector whichScan) {
385 386
       ints = sp->getIntensityArray();
386 387
       if (!mzs.get() || !ints.get()) {
387 388
 	Rcpp::NumericMatrix pks(0, 2);
389
+	Rcpp::colnames(pks) = cn;
388 390
 	res[i] = pks;
389 391
 	continue;
390 392
       }
... ...
@@ -394,6 +396,7 @@ Rcpp::List RcppPwiz::getPeakList(Rcpp::IntegerVector whichScan) {
394 396
       data.insert(data.end(), ints->data.begin(), ints->data.end());
395 397
       data_matrix = Rcpp::wrap(data);
396 398
       data_matrix.attr("dim") = Rcpp::Dimension(ints->data.size(), 2);
399
+      Rcpp::colnames(data_matrix) = cn;
397 400
       res[i] = data_matrix;
398 401
     }
399 402
     return res;