1. mzR
  2. Commit b4a2f2dad2cc3ae880ab1e8676a6350271d32330
Browse code

Merge pull request #207 from jorainer/master

Rewrite the pwiz header and peaks functions

Laurent Gatto authored on 02/10/2019 06:50:54 • GitHub committed on 02/10/2019 06:50:54
Showing 10 changed files
NEWS
History View file @ b4a2f2dad2
... ... 
@@ -1,3 +1,10 @@
1 
+CHANGES IN VERSION 2.19.5
2 
+-------------------------
3 
+ o header for the pwiz backend returns NA instead of 0 for not defined or
4 
+   missing information <2019-09-24 Tue>.
5 
+ o peaks for pwiz backend rewritten (small performance improvement)
6 
+   <2019-09-26 Thu>.
7 
+
1 8 
 CHANGES IN VERSION 2.19.4
2 9 
 -------------------------
3 10 
  o Add header columns scanWindowLowerLimit and scanWindowUpperLimit
R/methods-mzRpwiz.R
History View file @ b4a2f2dad2
... ... 
@@ -18,7 +18,6 @@ setMethod("chromatogramsInfo", "mzRpwiz",
18 18 
               .Defunct("chromatogram")
19 19 
           })
20 20
21 
-
22 21 
 setMethod("manufacturer", "mzRpwiz",
23 22 
           function(object) {
24 23 
             info <- instrumentInfo(object)
... ... 
@@ -82,19 +81,12 @@ setMethod("spectra", "mzRpwiz",
82 81
83 82 
 setMethod("peaksCount", c("mzRpwiz", "numeric"),
84 83 
           function(object, scans) {
85 
-            if (length(scans) == 1) {
86 
-              return(object@backend$getPeakList(scans)$peaksCount)
87 
-            } else {
88 
-                return(sapply(scans,
89 
-                              function(x)
90 
-                                  object@backend$getPeakList(x)$peaksCount))
91 
-            }
84 
+              lengths(object@backend$getPeakList(scans))/2
92 85 
           })
93 86
94 87 
 setMethod("peaksCount", c("mzRpwiz", "missing"),
95 88 
           function(object) {
96 
-            n <- length(object)
97 
-            return(peaksCount(object, 1:n))
89 
+            peaksCount(object, seq_along(object))
98 90 
           })
99 91
100 92 
 setMethod("runInfo", "mzRpwiz",
R/methods-mzRramp.R
History View file @ b4a2f2dad2
... ... 
@@ -18,18 +18,18 @@ setMethod("length",
18 18 
           function(x) return(x@backend$getLastScan()))
19 19
20 20 
 setMethod("peaks", "mzRramp",
21 
-          function(object, scans) .peaks(object, scans))
21 
+          function(object, scans) .peaks_ramp(object, scans))
22 22
23 23 
 setMethod("spectra", "mzRramp",
24 
-          function(object, scans) .peaks(object, scans))
24 
+          function(object, scans) .peaks_ramp(object, scans))
25 25
26 26 
 setMethod("peaksCount",
27 27 
           signature=c("mzRramp","numeric"),
28 28 
           function(object,scans) {
29 29 
               if (length(scans)==1) {
30 
-                  return(object@backend$getPeakList(scans)$peaksCount)
30 
+                  object@backend$getPeakList(scans)$peaksCount
31 31 
               } else {
32 
-                  return(sapply(scans,function(x) object@backend$getPeakList(x)$peaksCount))
32 
+                  sapply(scans,function(x) object@backend$getPeakList(x)$peaksCount)
33 33 
               }
34 34 
           })
35 35
R/utils.R
History View file @ b4a2f2dad2
... ... 
@@ -1,20 +1,26 @@
1 1 
 .peaks <- function(object, scans) {
2 
+    if (missing(scans))
3 
+        scans <- 1:length(object)
4 
+    res <- object@backend$getPeakList(scans)
5 
+    if (length(res) == 1)
6 
+        res[[1]]
7 
+    else res
8 
+}
9 
+
10 
+.peaks_ramp <- function(object, scans) {
2 11 
     if (missing(scans))
3 12 
         scans <- 1:length(object)
4 13 
     if (length(scans) == 1) {
5 
-        return(object@backend$getPeakList(scans)$peaks)
14 
+        object@backend$getPeakList(scans)$peaks
6 15 
     } else {
7 
-        return(sapply(scans,
8 
-                      function(x) object@backend$getPeakList(x)$peaks,
9 
-                      simplify = FALSE))
16 
+        sapply(scans, function(x) object@backend$getPeakList(x)$peaks,
17 
+               simplify = FALSE)
10 18 
     }
11 19 
 }
12 20
13 
-
14 21 
 setMethod("isolationWindow", "character",
15 22 
           function(object, ...) .isolationWindow(object, ...))
16 23
17 
-
18 24 
 .isolationWindow <- function(x, unique. = TRUE, simplify = TRUE) {
19 25 
     stopifnot(all(file.exists(x)))
20 26 
     if (!requireNamespace("XML"))
inst/unitTests/test_pwiz.R
History View file @ b4a2f2dad2
... ... 
@@ -83,38 +83,53 @@ test_getScanHeaderInfo <- function() {
83 83 
     ramp <- openMSfile(file, backend = "Ramp")
84 84 
     ## Read single scan header.
85 85 
     scan_3 <- header(mzml, scans = 3)
86 
+    scan_3_ramp <- header(ramp, scans = 3)
86 87 
     cn <- names(scan_3)
87 88 
     cn <- cn[!(cn %in% c("spectrumId", "scanWindowLowerLimit",
88 
-                         "scanWindowUpperLimit"))]
89 
-    scan_3_ramp <- header(ramp, scans = 3)
89 
+                         "scanWindowUpperLimit", "centroided"))]
90 
+    scan_3 <- scan_3[cn]
91 
+    scan_3_ramp <- scan_3_ramp[cn]
92 
+    ## Special case for columns that have an NA reported: they might have a
93 
+    ## 0 or NA in Ramp.
94 
+    nas <- is.na(scan_3)
90 95 
     ## Ramp does not read polarity
91 96 
     scan_3$polarity <- 0
92 
-    scan_3$centroided <- NA
93 
-    checkEquals(scan_3[cn], scan_3_ramp[cn])
94 
-
97 
+    checkEquals(scan_3[!nas], scan_3_ramp[!nas])
98 
+    ## All others have to be 0 or NA
99 
+    checkEquals(all(scan_3_ramp[nas] == 0 | is.na(scan_3_ramp[nas])), TRUE)
100 
+
95 101 
     ## Read all scan header
96 102 
     all_scans <- header(mzml)
97 103 
     all_scans_ramp <- header(ramp)
104 
+    all_scans <- all_scans[, cn]
105 
+    all_scans_ramp <- all_scans_ramp[, cn]
98 106 
     all_scans$polarity <- 0
99 
-    all_scans$centroided <- NA
100 
-    checkEquals(all_scans[, cn], all_scans_ramp[, cn])
101 
-
107 
+    nas <- vapply(all_scans, function(z) all(is.na(z)), logical(1))
108 
+    checkEquals(all_scans[, !nas], all_scans_ramp[, !nas])
109 
+    ## All others have to be 0 or NA
110 
+    checkEquals(all(all_scans_ramp[nas] == 0 | is.na(all_scans_ramp[nas])), TRUE)
111 
+
102 112 
     ## passing the index of all scan headers should return the same
103 
-    all_scans_2 <- header(mzml, scans = 1:nrow(all_scans))
104 
-    all_scans_ramp_2 <- header(ramp, scans = 1:nrow(all_scans))
105 
-    all_scans_2$polarity <- 0
106 
-    all_scans_2$centroided <- NA
107 
-    checkEquals(all_scans, all_scans_2)
108 
-    checkEquals(as.list(all_scans[3, ]), scan_3)
109 
-    checkEquals(all_scans_2[, cn], all_scans_ramp_2[, cn])
110 
-
113 
+    all_scans <- header(mzml, scans = 1:nrow(all_scans))
114 
+    all_scans_ramp <- header(ramp, scans = 1:nrow(all_scans))
115 
+    all_scans <- all_scans[, cn]
116 
+    all_scans_ramp <- all_scans_ramp[, cn]
117 
+    all_scans$polarity <- 0
118 
+    nas <- vapply(all_scans, function(z) all(is.na(z)), logical(1))
119 
+    checkEquals(all_scans[, !nas], all_scans_ramp[, !nas])
120 
+    ## All others have to be 0 or NA
121 
+    checkEquals(all(all_scans_ramp[, nas] == 0 | is.na(all_scans_ramp[, nas])), TRUE)
122 
+
111 123 
     ## Some selected scans
112 
-    scan_3 <- header(mzml, scans = c(3, 1, 14))
113 
-    scan_3_ramp <- header(ramp, scans = c(3, 1, 14))
114 
-    ## Ramp does not read polarity
115 
-    scan_3$polarity <- 0
116 
-    scan_3$centroided <- NA
117 
-    checkEquals(scan_3[, cn], scan_3_ramp[, cn])
124 
+    all_scans <- header(mzml, scans = c(3, 1, 14))
125 
+    all_scans_ramp <- header(ramp, scans = c(3, 1, 14))
126 
+    all_scans <- all_scans[, cn]
127 
+    all_scans_ramp <- all_scans_ramp[, cn]
128 
+    all_scans$polarity <- 0
129 
+    nas <- vapply(all_scans, function(z) all(is.na(z)), logical(1))
130 
+    checkEquals(all_scans[, !nas], all_scans_ramp[, !nas])
131 
+    ## All others have to be 0 or NA
132 
+    checkEquals(all(all_scans_ramp[, nas] == 0 | is.na(all_scans_ramp[, nas])), TRUE)
118 133
119 134 
     close(mzml)
120 135 
     close(ramp)
... ... 
@@ -126,33 +141,65 @@ test_getScanHeaderInfo <- function() {
126 141 
     ramp <- openMSfile(file, backend = "Ramp")
127 142 
     ## Read single scan header.
128 143 
     scan_3 <- header(mzml, scans = 3)
129 
-    scan_3$centroided <- NA
130 144 
     scan_3_ramp <- header(ramp, scans = 3)
131 
-    checkEquals(scan_3[cn], scan_3_ramp[cn])
145 
+    scan_3 <- scan_3[cn]
146 
+    scan_3_ramp <- scan_3_ramp[cn]
147 
+    nas <- vapply(scan_3, function(z) all(is.na(z)), logical(1))
148 
+    checkEquals(scan_3[!nas], scan_3_ramp[!nas])
149 
+    checkEquals(all(scan_3_ramp[nas] == 0 | is.na(scan_3_ramp[nas])), TRUE)
132 150
133 151 
     ## Read all scan header
134 152 
     all_scans <- header(mzml)
135 
-    all_scans$centroided <- NA
136 153 
     all_scans_ramp <- header(ramp)
154 
+    all_scans <- all_scans[, cn]
155 
+    all_scans_ramp <- all_scans_ramp[, cn]
156 
+
137 157 
     ## Ramp unable to read precursorScanNum from an mzXML file.
138 158 
     all_scans$precursorScanNum <- 0
139 
-    checkEquals(all_scans[, cn], all_scans_ramp[, cn])
159 
+    ## Replace all NA in all_scans with 0
160 
+    all_scans <- data.frame(lapply(all_scans, function(z) {
161 
+        z[is.na(z)] <- 0
162 
+        z
163 
+    }))
164 
+    all_scans_ramp <- data.frame(lapply(all_scans_ramp, function(z) {
165 
+        z[is.na(z)] <- 0
166 
+        z
167 
+    }))
168 
+    checkEquals(all_scans, all_scans_ramp)
140 169
141 170 
     ## passing the index of all scan headers should return the same
142 
-    all_scans_2 <- header(mzml, scans = 1:nrow(all_scans))
143 
-    all_scans_2$centroided <- NA
144 
-    all_scans_ramp_2 <- header(ramp, scans = 1:nrow(all_scans))
145 
-    all_scans_2$precursorScanNum <- 0
146 
-    checkEquals(all_scans, all_scans_2)
147 
-    checkEquals(as.list(all_scans[3, cn]), scan_3[cn])
148 
-    checkEquals(all_scans_2[, cn], all_scans_ramp_2[, cn])
171 
+    all_scans <- header(mzml, scans = 1:nrow(all_scans))
172 
+    all_scans_ramp <- header(ramp, scans = 1:nrow(all_scans))
173 
+    all_scans$precursorScanNum <- 0
174 
+    all_scans <- all_scans[, cn]
175 
+    all_scans_ramp <- all_scans_ramp[, cn]
176 
+    ## Replace all NA in all_scans with 0
177 
+    all_scans <- data.frame(lapply(all_scans, function(z) {
178 
+        z[is.na(z)] <- 0
179 
+        z
180 
+    }))
181 
+    all_scans_ramp <- data.frame(lapply(all_scans_ramp, function(z) {
182 
+        z[is.na(z)] <- 0
183 
+        z
184 
+    }))
185 
+    checkEquals(all_scans, all_scans_ramp)
149 186
150 187 
     ## Some selected scans
151 
-    scan_3 <- header(mzml, scans = c(3, 1, 14))
152 
-    scan_3_ramp <- header(ramp, scans = c(3, 1, 14))
153 
-    scan_3$precursorScanNum <- 0
154 
-    scan_3$centroided <- NA
155 
-    checkEquals(scan_3[, cn], scan_3_ramp[, cn])
188 
+    all_scans <- header(mzml, scans = c(3, 1, 14))
189 
+    all_scans_ramp <- header(ramp, scans = c(3, 1, 14))
190 
+    all_scans$precursorScanNum <- 0
191 
+    all_scans <- all_scans[, cn]
192 
+    all_scans_ramp <- all_scans_ramp[, cn]
193 
+    ## Replace all NA in all_scans with 0
194 
+    all_scans <- data.frame(lapply(all_scans, function(z) {
195 
+        z[is.na(z)] <- 0
196 
+        z
197 
+    }))
198 
+    all_scans_ramp <- data.frame(lapply(all_scans_ramp, function(z) {
199 
+        z[is.na(z)] <- 0
200 
+        z
201 
+    }))
202 
+    checkEquals(all_scans, all_scans_ramp)
156 203
157 204 
     close(mzml)
158 205 
     close(ramp)
... ... 
@@ -166,55 +213,62 @@ test_getScanHeaderInfo <- function() {
166 213 
     scan_3 <- header(mzml, scans = 3)
167 214 
     scan_3_ramp <- header(ramp, scans = 3)
168 215 
     ## Ramp does not read polarity, injectionTime or filterString
169 
-    scan_3$polarity <- 0
170 
-    scan_3$injectionTime <- 0
171 
-    scan_3$filterString <- NA_character_
172 
-    scan_3$centroided <- NA
173 
-    scan_3$isolationWindowTargetMZ <- NA_real_
174 
-    scan_3$isolationWindowLowerOffset <- NA_real_
175 
-    scan_3$isolationWindowUpperOffset <- NA_real_
176 
-    checkEquals(scan_3[cn], scan_3_ramp[cn])
177 
-
216 
+    cn <- names(scan_3)
217 
+    cn <- cn[!(cn %in% c("spectrumId", "scanWindowLowerLimit",
218 
+                         "scanWindowUpperLimit", "centroided",
219 
+                         "polarity", "filterString",
220 
+                         "isolationWindowTargetMZ",
221 
+                         "isolationWindowLowerOffset",
222 
+                         "isolationWindowUpperOffset",
223 
+                         "injectionTime"))]
224 
+    scan_3 <- scan_3[cn]
225 
+    scan_3_ramp <- scan_3_ramp[cn]
226 
+    scan_3[is.na(scan_3)] <- 0
227 
+    scan_3_ramp[is.na(scan_3_ramp)] <- 0
228 
+    checkEquals(scan_3, scan_3_ramp)
229 
+
178 230 
     ## Read all scan header
179 
-    all_scans <- header(mzml)
180 
-    all_scans_ramp <- header(ramp)
181 
-    ## Ramp does not read polarity, injectionTime or filterString
182 
-    all_scans$polarity <- 0
183 
-    all_scans$injectionTime <- 0
184 
-    all_scans$filterString <- NA_character_
185 
-    all_scans$centroided <- NA
186 
-    all_scans$isolationWindowTargetMZ <- NA_real_
187 
-    all_scans$isolationWindowLowerOffset <- NA_real_
188 
-    all_scans$isolationWindowUpperOffset <- NA_real_
189 
-    checkEquals(all_scans[, cn], all_scans_ramp[, cn])
231 
+    all_scans <- header(mzml)[, cn]
232 
+    all_scans_ramp <- header(ramp)[, cn]
233 
+    ## Replace all NA in all_scans with 0
234 
+    all_scans <- data.frame(lapply(all_scans, function(z) {
235 
+        z[is.na(z)] <- 0
236 
+        z
237 
+    }))
238 
+    all_scans_ramp <- data.frame(lapply(all_scans_ramp, function(z) {
239 
+        z[is.na(z)] <- 0
240 
+        z
241 
+    }))
242 
+    checkEquals(all_scans, all_scans_ramp)
190 243
191 244 
     ## passing the index of all scan headers should return the same
192 
-    all_scans_2 <- header(mzml, scans = 1:nrow(all_scans))
193 
-    all_scans_ramp_2 <- header(ramp, scans = 1:nrow(all_scans))
194 
-    ## Ramp does not read polarity, injectionTime or filterString
195 
-    all_scans_2$polarity <- 0
196 
-    all_scans_2$injectionTime <- 0
197 
-    all_scans_2$filterString <- NA_character_
198 
-    all_scans_2$centroided <- NA
199 
-    all_scans_2$isolationWindowTargetMZ <- NA_real_
200 
-    all_scans_2$isolationWindowLowerOffset <- NA_real_
201 
-    all_scans_2$isolationWindowUpperOffset <- NA_real_
202 
-    checkEquals(all_scans[, cn], all_scans_2[, cn])
203 
-    checkEquals(as.list(all_scans[3, cn]), scan_3[cn])
204 
-    checkEquals(all_scans_2[, cn], all_scans_ramp_2[, cn])
205 
-
245 
+    all_scans_2 <- header(mzml, scans = 1:nrow(all_scans))[, cn]
246 
+    all_scans_ramp <- header(ramp, scans = 1:nrow(all_scans))[, cn]
247 
+    ## Replace all NA in all_scans with 0
248 
+    all_scans_2 <- data.frame(lapply(all_scans_2, function(z) {
249 
+        z[is.na(z)] <- 0
250 
+        z
251 
+    }))
252 
+    all_scans_ramp <- data.frame(lapply(all_scans_ramp, function(z) {
253 
+        z[is.na(z)] <- 0
254 
+        z
255 
+    }))
256 
+    checkEquals(all_scans, all_scans_ramp)
257 
+    checkEquals(all_scans, all_scans_2)
258 
+
206 259 
     ## Some selected scans
207 
-    scan_3 <- header(mzml, scans = c(3, 1, 14))
208 
-    scan_3_ramp <- header(ramp, scans = c(3, 1, 14))
209 
-    ## Ramp does not read polarity, injectionTime or filterString
210 
-    scan_3$polarity <- 0
211 
-    scan_3$injectionTime <- 0
212 
-    scan_3$filterString <- NA_character_
213 
-    scan_3$centroided <- NA
214 
-    scan_3$isolationWindowTargetMZ <- NA_real_
215 
-    scan_3$isolationWindowLowerOffset <- NA_real_
216 
-    scan_3$isolationWindowUpperOffset <- NA_real_
217 
-    checkEquals(scan_3[, cn], scan_3_ramp[, cn])
260 
+    all_scans <- header(mzml, scans = c(3, 1, 14))[, cn]
261 
+    all_scans_ramp <- header(ramp, scans = c(3, 1, 14))[, cn]
262 
+    ## Replace all NA in all_scans with 0
263 
+    all_scans <- data.frame(lapply(all_scans, function(z) {
264 
+        z[is.na(z)] <- 0
265 
+        z
266 
+    }))
267 
+    all_scans_ramp <- data.frame(lapply(all_scans_ramp, function(z) {
268 
+        z[is.na(z)] <- 0
269 
+        z
270 
+    }))
271 
+    checkEquals(all_scans, all_scans_ramp)
218 272
219 273 
     close(mzml)
220 274 
     close(ramp)
inst/unitTests/test_write.R
History View file @ b4a2f2dad2
... ... 
@@ -482,7 +482,7 @@ test_writeMSData <- function() {
482 482 
     ## I don't quite understand that, but the acquisitionNum and the
483 483 
     ## precursorScanNum are different while the spectrumId is the same.
484 484 
     ## Still, check that the precursorScanNum is what we expect:
485 
-    checkEquals(hdr_new$precursorScanNum, c(0, 1, 0, 0, 4, 4))
485 
+    checkEquals(hdr_new$precursorScanNum, c(NA, 1, NA, NA, 4, 4))
486 486 
     hdr_new$acquisitionNum <- as.integer(factor(hdr_new$acquisitionNum))
487 487 
     hdr_new$precursorScanNum <- as.integer(factor(hdr_new$precursorScanNum))
488 488 
     hdr_sub$acquisitionNum <- as.integer(factor(hdr_sub$acquisitionNum))
... ... 
@@ -506,7 +506,7 @@ test_writeMSData <- function() {
506 506 
     ii_new <- mzR::instrumentInfo(mzml_new)
507 507 
     mzR::close(mzml_new)
508 508 
     checkEquals(pks_new, pks_sub)
509 
-    checkEquals(hdr_new$precursorScanNum, c(0, 1, 0, 0, 4, 4))
509 
+    checkEquals(hdr_new$precursorScanNum, c(NA, 1, NA, NA, 4, 4))
510 510 
     rownames(hdr_sub) <- NULL
511 511 
     rownames(hdr_new) <- NULL
512 512 
     hdr_sub$acquisitionNum <- as.integer(factor(hdr_sub$acquisitionNum))
... ... 
@@ -562,16 +562,25 @@ test_writeMSData <- function() {
562 562 
     test_file <- system.file("iontrap", "extracted.mzData", package = "msdata")
563 563 
     in_file <- openMSfile(test_file, backend = "Ramp")
564 564 
     hdr <- header(in_file)
565 
+    hdr_orig <- hdr
565 566 
     pks <- peaks(in_file)
566 567 
     mzR::close(in_file)
567 568
568 569 
     ## mzML
569 570 
     out_file <- paste0(test_folder, "/test_write.mzML")
570 
-    writeMSData(file = out_file, header = hdr, object = pks)
571 
+    writeMSData(file = out_file, header = hdr_orig, object = pks)
571 572 
     in_file <- openMSfile(out_file, backend = "pwiz")
572 573 
     hdr_2 <- header(in_file)
573 574 
     pks_2 <- peaks(in_file)
574 575 
     mzR::close(in_file)
576 
+    hdr <- data.frame(lapply(hdr, function(z) {
577 
+        z[is.na(z)] <- 0
578 
+        z
579 
+    }))
580 
+    hdr_2 <- data.frame(lapply(hdr_2, function(z) {
581 
+        z[is.na(z)] <- 0
582 
+        z
583 
+    }))
575 584 
     checkEquals(hdr, hdr_2)
576 585 
     checkEquals(pks, pks_2)
577 586 
     ## validate mzML:
... ... 
@@ -581,13 +590,17 @@ test_writeMSData <- function() {
581 590
582 591 
     ## mzXML output:
583 592 
     out_file <- paste0(test_folder, "test_write.mzXML")
584 
-    writeMSData(file = out_file, header = hdr, object = pks,
593 
+    writeMSData(file = out_file, header = hdr_orig, object = pks,
585 594 
                 outformat = "mzXML")
586 595 
     in_file <- openMSfile(out_file, backend = "pwiz")
587 596 
     hdr_2 <- header(in_file)
597 
+    hdr_2 <- data.frame(lapply(hdr_2, function(z) {
598 
+        z[is.na(z)] <- 0
599 
+        z
600 
+    }))
588 601 
     pks_2 <- peaks(in_file)
589 602 
     mzR::close(in_file)
590 603 
     checkEquals(pks, pks_2)
591 
-    hdr$centroided <- FALSE
592 
-    checkEquals(hdr, hdr_2)
604 
+    checkEquals(hdr[, colnames(hdr_2) != "centroided"],
605 
+                hdr_2[, colnames(hdr_2) != "centroided"])
593 606 
 }
man/metadata.Rd
History View file @ b4a2f2dad2
... ... 
@@ -42,7 +42,7 @@ psms(object, ...)
42 42 
 substitutions(object)
43 43 
 database(object, ...)
44 44 
 enzymes(object)
45 
-tolerance(object)
45 
+tolerance(object, ...)
46 46 
 score(x, ...)
47 47 
 para(object)
48 48 
 specParams(object)
man/peaks.Rd
History View file @ b4a2f2dad2
... ... 
@@ -128,7 +128,7 @@
128 128 
   \code{isolationWindowUpperOffset}, \code{scanWindowLowerLimit} and
129 129 
   \code{scanWindowUpperLimit}. If multiple scans are queried, a
130 130 
   \code{data.frame} is returned with the scans reported along the
131 
-  rows.
131 
+  rows. For missing or not defined spectrum variables \code{NA} is reported.
132 132
133 133 
   The \code{get3Dmap} function performs a simple resampling between
134 134 
   \code{lowMz} and \code{highMz} with \code{reMz} resolution. A matrix
src/RcppPwiz.cpp
History View file @ b4a2f2dad2
... ... 
@@ -38,43 +38,6 @@ void RcppPwiz::close()
38 38 
     }
39 39 
 }
40 40
41 
-
42 
-// void RcppPwiz::writeMSfile(const string& file, const string& format)
43 
-// {
44 
-//     if (msd != NULL)
45 
-//     {
46 
-//         if(format == "mgf")
47 
-//         {
48 
-//             std::ofstream* mgfOutFileP = new std::ofstream(file.c_str());
49 
-//             Serializer_MGF serializerMGF;
50 
-//             serializerMGF.write(*mgfOutFileP, *msd);
51 
-//             mgfOutFileP->flush();
52 
-//             mgfOutFileP->close();
53 
-//         }
54 
-//         else if(format == "mzxml")
55 
-//         {
56 
-//             std::ofstream mzXMLOutFileP(file.c_str());
57 
-//             Serializer_mzXML::Config config;
58 
-//             config.binaryDataEncoderConfig.compression = BinaryDataEncoder::Compression_Zlib;
59 
-//             Serializer_mzXML serializerMzXML(config);
60 
-//             serializerMzXML.write(mzXMLOutFileP, *msd);
61 
-//         }
62 
-//         else if(format == "mzml")
63 
-//         {
64 
-//             std::ofstream mzXMLOutFileP(file.c_str());
65 
-//             Serializer_mzML::Config config;
66 
-//             config.binaryDataEncoderConfig.compression = BinaryDataEncoder::Compression_Zlib;
67 
-//             Serializer_mzML mzmlSerializer(config);
68 
-//             mzmlSerializer.write(mzXMLOutFileP, *msd);
69 
-//         }
70 
-//         else
71 
-//             Rcpp::Rcerr << format << " format not supported! Please try mgf, mzML, mzXML or mz5." << std::endl;
72 
-//     }
73 
-//     else
74 
-//         Rcpp::Rcerr << "No pwiz object available! Please open a file first!" << std::endl;
75 
-// }
76 
-
77 
-
78 41 
 string RcppPwiz::getFilename() {
79 42 
   return filename;
80 43 
 }
... ... 
@@ -164,249 +127,279 @@ Rcpp::List RcppPwiz::getInstrumentInfo ( )
164 127 
   return instrumentInfo;
165 128 
 }
166 129
167 
-
168 
-Rcpp::DataFrame RcppPwiz::getScanHeaderInfo (Rcpp::IntegerVector whichScan)
169 
-{
170 
-  if (msd != NULL)
171 
-    {
172 
-      SpectrumListPtr slp = msd->run.spectrumListPtr;
173 
-      int N = slp->size();
174 
-
175 
-      int N_scans = whichScan.size();
176 
-
177 
-      ScanHeaderStruct scanHeader;
178 
-      RAMPAdapter * adapter = new  RAMPAdapter(filename);
179 
-      Rcpp::IntegerVector seqNum(N_scans); // number in sequence observed file (1-based)
180 
-      Rcpp::IntegerVector acquisitionNum(N_scans); // scan number as declared in File (may be gaps)
181 
-      Rcpp::IntegerVector msLevel(N_scans);
182 
-      Rcpp::IntegerVector polarity(N_scans);
183 
-      Rcpp::IntegerVector peaksCount(N_scans);
184 
-      Rcpp::NumericVector totIonCurrent(N_scans);
185 
-      Rcpp::NumericVector retentionTime(N_scans);        /* in seconds */
186 
-      Rcpp::NumericVector basePeakMZ(N_scans);
187 
-      Rcpp::NumericVector basePeakIntensity(N_scans);
188 
-      Rcpp::NumericVector collisionEnergy(N_scans);
189 
-      Rcpp::NumericVector ionisationEnergy(N_scans);
190 
-      Rcpp::NumericVector lowMZ(N_scans);
191 
-      Rcpp::NumericVector highMZ(N_scans);
192 
-      Rcpp::IntegerVector precursorScanNum(N_scans); /* only if MS level > 1 */
193 
-      Rcpp::NumericVector precursorMZ(N_scans);  /* only if MS level > 1 */
194 
-      Rcpp::IntegerVector precursorCharge(N_scans);  /* only if MS level > 1 */
195 
-      Rcpp::NumericVector precursorIntensity(N_scans);  /* only if MS level > 1 */
196 
-      //char scanType[SCANTYPE_LENGTH];
197 
-      //char activationMethod[SCANTYPE_LENGTH];
198 
-      //char possibleCharges[SCANTYPE_LENGTH];
199 
-      //int numPossibleCharges;
200 
-      //bool possibleChargesArray[CHARGEARRAY_LENGTH]; /* NOTE: does NOT include "precursorCharge" information; only from "possibleCharges" */
201 
-      Rcpp::IntegerVector mergedScan(N_scans);  /* only if MS level > 1 */
202 
-      Rcpp::IntegerVector mergedResultScanNum(N_scans); /* scan number of the resultant merged scan */
203 
-      Rcpp::IntegerVector mergedResultStartScanNum(N_scans); /* smallest scan number of the scanOrigin for merged scan */
204 
-      Rcpp::IntegerVector mergedResultEndScanNum(N_scans); /* largest scan number of the scanOrigin for merged scan */
205 
-      Rcpp::NumericVector ionInjectionTime(N_scans); /* The time spent filling an ion trapping device*/
206 
-      Rcpp::StringVector filterString(N_scans);
207 
-      Rcpp::StringVector spectrumId(N_scans);
208 
-      Rcpp::LogicalVector centroided(N_scans);
209 
-      Rcpp::NumericVector ionMobilityDriftTime(N_scans);
210 
-      Rcpp::NumericVector isolationWindowTargetMZ(N_scans);
211 
-      Rcpp::NumericVector isolationWindowLowerOffset(N_scans);
212 
-      Rcpp::NumericVector isolationWindowUpperOffset(N_scans);
213 
-      Rcpp::NumericVector scanWindowLowerLimit(N_scans);
214 
-      Rcpp::NumericVector scanWindowUpperLimit(N_scans);
215 
-
216 
-      for (int i = 0; i < N_scans; i++)
217 
-	{
218 
-	  int current_scan = whichScan[i];
219 
-	  adapter->getScanHeader(current_scan - 1, scanHeader, false);
220 
-	  seqNum[i] = scanHeader.seqNum;
221 
-	  acquisitionNum[i] = scanHeader.acquisitionNum;
222 
-	  msLevel[i] = scanHeader.msLevel;
223 
-
224 
-	  SpectrumPtr sp = slp->spectrum(current_scan-1, false); // Is TRUE neccessary here ?
225 
-	  Scan dummy;
226 
-	  Scan& scan = sp->scanList.scans.empty() ? dummy : sp->scanList.scans[0];
227 
-	  CVParam param = sp->cvParamChild(MS_scan_polarity);
228 
-	  polarity[i] = (param.cvid==MS_negative_scan ? 0 : (param.cvid==MS_positive_scan ? +1 : -1 ) );
229 
-	  // ionInjectionTime[i] = scan.cvParam(MS_ion_injection_time).valueAs<double>();
230 
-	  ionInjectionTime[i] = (scan.cvParam(MS_ion_injection_time).timeInSeconds() * 1000);
231 
-	  filterString[i] = scan.cvParam(MS_filter_string).value.empty() ? NA_STRING : Rcpp::String(scan.cvParam(MS_filter_string).value);
232 
-	  ionMobilityDriftTime[i] = scan.cvParam(MS_ion_mobility_drift_time).value.empty() ? NA_REAL : (scan.cvParam(MS_ion_mobility_drift_time).timeInSeconds() * 1000);
233 
-
234 
-	  if (!scan.scanWindows.empty()) {
235 
-	    scanWindowLowerLimit[i] = scan.scanWindows[0].cvParam(MS_scan_window_lower_limit).valueAs<double>();
236 
-	    scanWindowUpperLimit[i] = scan.scanWindows[0].cvParam(MS_scan_window_upper_limit).valueAs<double>();
237 
-	  } else {
238 
-	    scanWindowLowerLimit[i] = NA_REAL;
239 
-	    scanWindowUpperLimit[i] = NA_REAL;
240 
-	  }
241 
-
242 
-	  if (!sp->precursors.empty()) {
243 
-	    IsolationWindow iwin = sp->precursors[0].isolationWindow;
244 
-	    if (!iwin.empty()) {
245 
-	      isolationWindowTargetMZ[i] = iwin.cvParam(MS_isolation_window_target_m_z).value.empty() ? NA_REAL : iwin.cvParam(MS_isolation_window_target_m_z).valueAs<double>();
246 
-	      isolationWindowLowerOffset[i] = iwin.cvParam(MS_isolation_window_lower_offset).value.empty() ? NA_REAL : iwin.cvParam(MS_isolation_window_lower_offset).valueAs<double>();
247 
-	      isolationWindowUpperOffset[i] = iwin.cvParam(MS_isolation_window_upper_offset).value.empty() ? NA_REAL : iwin.cvParam(MS_isolation_window_upper_offset).valueAs<double>();
248 
-	    } else {
249 
-	      isolationWindowTargetMZ[i] = NA_REAL;
250 
-	      isolationWindowLowerOffset[i] = NA_REAL;
251 
-	      isolationWindowUpperOffset[i] = NA_REAL;
252 
-	    }
130 
+Rcpp::DataFrame RcppPwiz::getScanHeaderInfo (Rcpp::IntegerVector whichScan) {
131 
+  if (msd != NULL) {
132 
+    SpectrumListPtr slp = msd->run.spectrumListPtr;
133 
+    int N = slp->size();
134 
+    int N_scans = whichScan.size();
135 
+    CVID nativeIdFormat = id::getDefaultNativeIDFormat(*msd);
136 
+    Rcpp::IntegerVector seqNum(N_scans); // number in sequence observed file (1-based)
137 
+    Rcpp::IntegerVector acquisitionNum(N_scans); // scan number as declared in File (may be gaps)
138 
+    Rcpp::IntegerVector msLevel(N_scans);
139 
+    Rcpp::IntegerVector polarity(N_scans);
140 
+    Rcpp::IntegerVector peaksCount(N_scans);
141 
+    Rcpp::NumericVector totIonCurrent(N_scans);
142 
+    Rcpp::NumericVector retentionTime(N_scans);        /* in seconds */
143 
+    Rcpp::NumericVector basePeakMZ(N_scans);
144 
+    Rcpp::NumericVector basePeakIntensity(N_scans);
145 
+    Rcpp::NumericVector collisionEnergy(N_scans);
146 
+    Rcpp::NumericVector ionisationEnergy(N_scans);
147 
+    Rcpp::NumericVector lowMZ(N_scans);
148 
+    Rcpp::NumericVector highMZ(N_scans);
149 
+    Rcpp::IntegerVector precursorScanNum(N_scans); /* only if MS level > 1 */
150 
+    Rcpp::NumericVector precursorMZ(N_scans);  /* only if MS level > 1 */
151 
+    Rcpp::IntegerVector precursorCharge(N_scans);  /* only if MS level > 1 */
152 
+    Rcpp::NumericVector precursorIntensity(N_scans);  /* only if MS level > 1 */
153 
+    Rcpp::IntegerVector mergedScan(N_scans);  /* only if MS level > 1 */
154 
+    Rcpp::IntegerVector mergedResultScanNum(N_scans); /* scan number of the resultant merged scan */
155 
+    Rcpp::IntegerVector mergedResultStartScanNum(N_scans); /* smallest scan number of the scanOrigin for merged scan */
156 
+    Rcpp::IntegerVector mergedResultEndScanNum(N_scans); /* largest scan number of the scanOrigin for merged scan */
157 
+    Rcpp::NumericVector ionInjectionTime(N_scans); /* The time spent filling an ion trapping device*/
158 
+    Rcpp::StringVector filterString(N_scans);
159 
+    Rcpp::StringVector spectrumId(N_scans);
160 
+    Rcpp::LogicalVector centroided(N_scans);
161 
+    Rcpp::NumericVector ionMobilityDriftTime(N_scans);
162 
+    Rcpp::NumericVector isolationWindowTargetMZ(N_scans);
163 
+    Rcpp::NumericVector isolationWindowLowerOffset(N_scans);
164 
+    Rcpp::NumericVector isolationWindowUpperOffset(N_scans);
165 
+    Rcpp::NumericVector scanWindowLowerLimit(N_scans);
166 
+    Rcpp::NumericVector scanWindowUpperLimit(N_scans);
167 
+
168 
+    for (int i = 0; i < N_scans; i++) {
169 
+      int current_scan = whichScan[i];
170 
+      SpectrumPtr sp = slp->spectrum(current_scan - 1, false);
171 
+      Scan dummy;
172 
+      Scan& scan = sp->scanList.scans.empty() ? dummy : sp->scanList.scans[0];
173 
+      // seqNum
174 
+      seqNum[i] = current_scan;
175 
+      // acquisitionNum
176 
+      string id = sp->id;
177 
+      string scanNumber = id::translateNativeIDToScanNumber(nativeIdFormat, id);
178 
+      if (scanNumber.empty()) {
179 
+	acquisitionNum[i] = current_scan;
180 
+      } else {
181 
+	acquisitionNum[i] = lexical_cast<int>(scanNumber);
182 
+      }
183 
+      // spectrumId
184 
+      spectrumId[i] = sp->id;
185 
+      // msLevel
186 
+      msLevel[i] = sp->cvParam(MS_ms_level).valueAs<int>();
187 
+      // peaksCount
188 
+      peaksCount[i] = static_cast<int>(sp->defaultArrayLength);
189 
+      // totIonCurrent
190 
+      totIonCurrent[i] = sp->cvParam(MS_total_ion_current).valueAs<double>();
191 
+      // basePeakMZ
192 
+      basePeakMZ[i] = sp->cvParam(MS_base_peak_m_z).valueAs<double>();
193 
+      // basePeakIntensity
194 
+      basePeakIntensity[i] = sp->cvParam(MS_base_peak_intensity).valueAs<double>();
195 
+      // ionisationEnerty
196 
+      ionisationEnergy[i] = sp->cvParam(MS_ionization_energy_OBSOLETE).valueAs<double>();
197 
+      // lowMZ
198 
+      lowMZ[i] = sp->cvParam(MS_lowest_observed_m_z).valueAs<double>();
199 
+      // highMZ
200 
+      highMZ[i] = sp->cvParam(MS_highest_observed_m_z).valueAs<double>();
201 
+      // polarity
202 
+      CVParam param = sp->cvParamChild(MS_scan_polarity);
203 
+      polarity[i] = (param.cvid==MS_negative_scan ? 0 : (param.cvid==MS_positive_scan ? +1 : -1 ) );
204 
+      // centroided
205 
+      param = sp->cvParamChild(MS_spectrum_representation);
206 
+      centroided[i] = (param.cvid==MS_centroid_spectrum ? TRUE : (param.cvid==MS_profile_spectrum ? FALSE : NA_LOGICAL));
207 
+      // retentionTime
208 
+      retentionTime[i] = scan.cvParam(MS_scan_start_time).timeInSeconds();
209 
+      // ionInjectionTime
210 
+      ionInjectionTime[i] = (scan.cvParam(MS_ion_injection_time).timeInSeconds() * 1000);
211 
+      // filterString
212 
+      filterString[i] = scan.cvParam(MS_filter_string).value.empty() ? NA_STRING : Rcpp::String(scan.cvParam(MS_filter_string).value);
213 
+      // ionMobilityDriftTime
214 
+      ionMobilityDriftTime[i] = scan.cvParam(MS_ion_mobility_drift_time).value.empty() ? NA_REAL : (scan.cvParam(MS_ion_mobility_drift_time).timeInSeconds() * 1000);
215 
+      // scanWindowLowerLimit and scanWindowUpperLimit
216 
+      if (!scan.scanWindows.empty()) {
217 
+	scanWindowLowerLimit[i] = scan.scanWindows[0].cvParam(MS_scan_window_lower_limit).valueAs<double>();
218 
+	scanWindowUpperLimit[i] = scan.scanWindows[0].cvParam(MS_scan_window_upper_limit).valueAs<double>();
219 
+      } else {
220 
+	scanWindowLowerLimit[i] = NA_REAL;
221 
+	scanWindowUpperLimit[i] = NA_REAL;
222 
+      }
223 
+      // mergedScan - also not supported by RAMPAdapter
224 
+      mergedScan[i] = NA_INTEGER;
225 
+      mergedResultScanNum[i] = NA_INTEGER;
226 
+      mergedResultStartScanNum[i] = NA_INTEGER;
227 
+      mergedResultEndScanNum[i] = NA_INTEGER;
228 
+      if (!sp->precursors.empty()) {
229 
+	const Precursor& precursor = sp->precursors[0];
230 
+	// collisionEnergy
231 
+	collisionEnergy[i] = precursor.activation.cvParam(MS_collision_energy).valueAs<double>();
232 
+	// precursorScanNum
233 
+	size_t precursorIndex = slp->find(precursor.spectrumID);
234 
+	if (precursorIndex < slp->size()) {
235 
+	  const SpectrumIdentity& precursorSpectrum = slp->spectrumIdentity(precursorIndex);
236 
+	  string precursorScanNumber = id::translateNativeIDToScanNumber(nativeIdFormat, precursorSpectrum.id);
237 
+	  if (precursorScanNumber.empty()) {
238 
+	    precursorScanNum[i] = precursorIndex + 1;
253 239 
 	  } else {
254 
-	      isolationWindowTargetMZ[i] = NA_REAL;
255 
-	      isolationWindowLowerOffset[i] = NA_REAL;
256 
-	      isolationWindowUpperOffset[i] = NA_REAL;
240 
+	    precursorScanNum[i] = lexical_cast<int>(precursorScanNumber);
257 241 
 	  }
258 
-	  peaksCount[i] = scanHeader.peaksCount;
259 
-	  totIonCurrent[i] = scanHeader.totIonCurrent;
260 
-	  retentionTime[i] = scanHeader.retentionTime;
261 
-	  basePeakMZ[i] = scanHeader.basePeakMZ;
262 
-	  basePeakIntensity[i] = scanHeader.basePeakIntensity;
263 
-	  collisionEnergy[i] = scanHeader.collisionEnergy;
264 
-	  ionisationEnergy[i] = scanHeader.ionisationEnergy;
265 
-	  lowMZ[i] = scanHeader.lowMZ;
266 
-	  highMZ[i] = scanHeader.highMZ;
267 
-	  precursorScanNum[i] = scanHeader.precursorScanNum;
268 
-	  precursorMZ[i] = scanHeader.precursorMZ;
269 
-	  precursorCharge[i] = scanHeader.precursorCharge;
270 
-	  precursorIntensity[i] = scanHeader.precursorIntensity;
271 
-	  mergedScan[i] = scanHeader.mergedScan;
272 
-	  mergedResultScanNum[i] = scanHeader.mergedResultScanNum;
273 
-	  mergedResultStartScanNum[i] = scanHeader.mergedResultStartScanNum;
274 
-	  mergedResultEndScanNum[i] = scanHeader.mergedResultEndScanNum;
275 
-	  spectrumId[i] = sp->id;
276 
-	  CVParam prm = sp->cvParamChild(MS_spectrum_representation);
277 
-	  centroided[i] = (prm.cvid==MS_centroid_spectrum ? TRUE : (prm.cvid==MS_profile_spectrum ? FALSE : NA_LOGICAL));
242 
+	} else {
243 
+	  precursorScanNum[i] = NA_INTEGER;
278 244 
 	}
279 
-      // delete adapter issue #64
280 
-      delete adapter;
281 
-      adapter = NULL;
282 
-
283 
-      Rcpp::List header(31);
284 
-      std::vector<std::string> names;
285 
-      int i = 0;
286 
-      names.push_back("seqNum");
287 
-      header[i++] = Rcpp::wrap(seqNum);
288 
-      names.push_back("acquisitionNum");
289 
-      header[i++] = Rcpp::wrap(acquisitionNum);
290 
-      names.push_back("msLevel");
291 
-      header[i++] = Rcpp::wrap(msLevel);
292 
-      names.push_back("polarity");
293 
-      header[i++] = Rcpp::wrap(polarity);
294 
-      names.push_back("peaksCount");
295 
-      header[i++] = Rcpp::wrap(peaksCount);
296 
-      names.push_back("totIonCurrent");
297 
-      header[i++] = Rcpp::wrap(totIonCurrent);
298 
-      names.push_back("retentionTime");
299 
-      header[i++] = Rcpp::wrap(retentionTime);
300 
-      names.push_back("basePeakMZ");
301 
-      header[i++] = Rcpp::wrap(basePeakMZ);
302 
-      names.push_back("basePeakIntensity");
303 
-      header[i++] = Rcpp::wrap(basePeakIntensity);
304 
-      names.push_back("collisionEnergy");
305 
-      header[i++] = Rcpp::wrap(collisionEnergy);
306 
-      names.push_back("ionisationEnergy");
307 
-      header[i++] = Rcpp::wrap(ionisationEnergy);
308 
-      names.push_back("lowMZ");
309 
-      header[i++] = Rcpp::wrap(lowMZ);
310 
-      names.push_back("highMZ");
311 
-      header[i++] = Rcpp::wrap(highMZ);
312 
-      names.push_back("precursorScanNum");
313 
-      header[i++] = Rcpp::wrap(precursorScanNum);
314 
-      names.push_back("precursorMZ");
315 
-      header[i++] = Rcpp::wrap(precursorMZ);
316 
-      names.push_back("precursorCharge");
317 
-      header[i++] = Rcpp::wrap(precursorCharge);
318 
-      names.push_back("precursorIntensity");
319 
-      header[i++] = Rcpp::wrap(precursorIntensity);
320 
-      names.push_back("mergedScan");
321 
-      header[i++] = Rcpp::wrap(mergedScan);
322 
-      names.push_back("mergedResultScanNum");
323 
-      header[i++] = Rcpp::wrap(mergedResultScanNum);
324 
-      names.push_back("mergedResultStartScanNum");
325 
-      header[i++] = Rcpp::wrap(mergedResultStartScanNum);
326 
-      names.push_back("mergedResultEndScanNum");
327 
-      header[i++] = Rcpp::wrap(mergedResultEndScanNum);
328 
-      names.push_back("injectionTime");
329 
-      header[i++] = Rcpp::wrap(ionInjectionTime);
330 
-      names.push_back("filterString");
331 
-      header[i++] = Rcpp::wrap(filterString);
332 
-      names.push_back("spectrumId");
333 
-      header[i++] = Rcpp::wrap(spectrumId);
334 
-      names.push_back("centroided");
335 
-      header[i++] = Rcpp::wrap(centroided);
336 
-      names.push_back("ionMobilityDriftTime");
337 
-      header[i++] = Rcpp::wrap(ionMobilityDriftTime);
338 
-      names.push_back("isolationWindowTargetMZ");
339 
-      header[i++] = Rcpp::wrap(isolationWindowTargetMZ);
340 
-      names.push_back("isolationWindowLowerOffset");
341 
-      header[i++] = Rcpp::wrap(isolationWindowLowerOffset);
342 
-      names.push_back("isolationWindowUpperOffset");
343 
-      header[i++] = Rcpp::wrap(isolationWindowUpperOffset);
344 
-      names.push_back("scanWindowLowerLimit");
345 
-      header[i++] = Rcpp::wrap(scanWindowLowerLimit);
346 
-      names.push_back("scanWindowUpperLimit");
347 
-      header[i++] = Rcpp::wrap(scanWindowUpperLimit);
348 
-      header.attr("names") = names;
349 
-
350 
-      return header;
245 
+	// precursorMZ, precursorCharge, precursorIntensity
246 
+	if (!precursor.selectedIons.empty()) {
247 
+	  precursorMZ[i] = precursor.selectedIons[0].cvParam(MS_selected_ion_m_z).value.empty() ? precursor.selectedIons[0].cvParam(MS_m_z).valueAs<double>() : precursor.selectedIons[0].cvParam(MS_selected_ion_m_z).valueAs<double>();
248 
+	  precursorCharge[i] = precursor.selectedIons[0].cvParam(MS_charge_state).valueAs<int>();
249 
+	  precursorIntensity[i] = precursor.selectedIons[0].cvParam(MS_peak_intensity).valueAs<double>();
250 
+	}
251 
+	// isolationWindowTargetMZ, ...
252 
+	IsolationWindow iwin = sp->precursors[0].isolationWindow;
253 
+	if (!iwin.empty()) {
254 
+	  isolationWindowTargetMZ[i] = iwin.cvParam(MS_isolation_window_target_m_z).value.empty() ? NA_REAL : iwin.cvParam(MS_isolation_window_target_m_z).valueAs<double>();
255 
+	  isolationWindowLowerOffset[i] = iwin.cvParam(MS_isolation_window_lower_offset).value.empty() ? NA_REAL : iwin.cvParam(MS_isolation_window_lower_offset).valueAs<double>();
256 
+	  isolationWindowUpperOffset[i] = iwin.cvParam(MS_isolation_window_upper_offset).value.empty() ? NA_REAL : iwin.cvParam(MS_isolation_window_upper_offset).valueAs<double>();
257 
+	} else {
258 
+	  isolationWindowTargetMZ[i] = NA_REAL;
259 
+	  isolationWindowLowerOffset[i] = NA_REAL;
260 
+	  isolationWindowUpperOffset[i] = NA_REAL;
261 
+	}
262 
+      } else {
263 
+	collisionEnergy[i] = NA_REAL;
264 
+	precursorScanNum[i] = NA_INTEGER;
265 
+	precursorMZ[i] = NA_REAL;
266 
+	precursorCharge[i] = NA_INTEGER;
267 
+	precursorIntensity[i] = NA_REAL;
268 
+	mergedScan[i] = NA_INTEGER;
269 
+	mergedResultScanNum[i] = NA_INTEGER;
270 
+	mergedResultStartScanNum[i] = NA_INTEGER;
271 
+	mergedResultEndScanNum[i] = NA_INTEGER;
272 
+	isolationWindowTargetMZ[i] = NA_REAL;
273 
+	isolationWindowLowerOffset[i] = NA_REAL;
274 
+	isolationWindowUpperOffset[i] = NA_REAL;
275 
+      }
351 276 
     }
277 
+
278 
+    Rcpp::List header(31);
279 
+    std::vector<std::string> names;
280 
+    int i = 0;
281 
+    names.push_back("seqNum");
282 
+    header[i++] = Rcpp::wrap(seqNum);
283 
+    names.push_back("acquisitionNum");
284 
+    header[i++] = Rcpp::wrap(acquisitionNum);
285 
+    names.push_back("msLevel");
286 
+    header[i++] = Rcpp::wrap(msLevel);
287 
+    names.push_back("polarity");
288 
+    header[i++] = Rcpp::wrap(polarity);
289 
+    names.push_back("peaksCount");
290 
+    header[i++] = Rcpp::wrap(peaksCount);
291 
+    names.push_back("totIonCurrent");
292 
+    header[i++] = Rcpp::wrap(totIonCurrent);
293 
+    names.push_back("retentionTime");
294 
+    header[i++] = Rcpp::wrap(retentionTime);
295 
+    names.push_back("basePeakMZ");
296 
+    header[i++] = Rcpp::wrap(basePeakMZ);
297 
+    names.push_back("basePeakIntensity");
298 
+    header[i++] = Rcpp::wrap(basePeakIntensity);
299 
+    names.push_back("collisionEnergy");
300 
+    header[i++] = Rcpp::wrap(collisionEnergy);
301 
+    names.push_back("ionisationEnergy");
302 
+    header[i++] = Rcpp::wrap(ionisationEnergy);
303 
+    names.push_back("lowMZ");
304 
+    header[i++] = Rcpp::wrap(lowMZ);
305 
+    names.push_back("highMZ");
306 
+    header[i++] = Rcpp::wrap(highMZ);
307 
+    names.push_back("precursorScanNum");
308 
+    header[i++] = Rcpp::wrap(precursorScanNum);
309 
+    names.push_back("precursorMZ");
310 
+    header[i++] = Rcpp::wrap(precursorMZ);
311 
+    names.push_back("precursorCharge");
312 
+    header[i++] = Rcpp::wrap(precursorCharge);
313 
+    names.push_back("precursorIntensity");
314 
+    header[i++] = Rcpp::wrap(precursorIntensity);
315 
+    names.push_back("mergedScan");
316 
+    header[i++] = Rcpp::wrap(mergedScan);
317 
+    names.push_back("mergedResultScanNum");
318 
+    header[i++] = Rcpp::wrap(mergedResultScanNum);
319 
+    names.push_back("mergedResultStartScanNum");
320 
+    header[i++] = Rcpp::wrap(mergedResultStartScanNum);
321 
+    names.push_back("mergedResultEndScanNum");
322 
+    header[i++] = Rcpp::wrap(mergedResultEndScanNum);
323 
+    names.push_back("injectionTime");
324 
+    header[i++] = Rcpp::wrap(ionInjectionTime);
325 
+    names.push_back("filterString");
326 
+    header[i++] = Rcpp::wrap(filterString);
327 
+    names.push_back("spectrumId");
328 
+    header[i++] = Rcpp::wrap(spectrumId);
329 
+    names.push_back("centroided");
330 
+    header[i++] = Rcpp::wrap(centroided);
331 
+    names.push_back("ionMobilityDriftTime");
332 
+    header[i++] = Rcpp::wrap(ionMobilityDriftTime);
333 
+    names.push_back("isolationWindowTargetMZ");
334 
+    header[i++] = Rcpp::wrap(isolationWindowTargetMZ);
335 
+    names.push_back("isolationWindowLowerOffset");
336 
+    header[i++] = Rcpp::wrap(isolationWindowLowerOffset);
337 
+    names.push_back("isolationWindowUpperOffset");
338 
+    header[i++] = Rcpp::wrap(isolationWindowUpperOffset);
339 
+    names.push_back("scanWindowLowerLimit");
340 
+    header[i++] = Rcpp::wrap(scanWindowLowerLimit);
341 
+    names.push_back("scanWindowUpperLimit");
342 
+    header[i++] = Rcpp::wrap(scanWindowUpperLimit);
343 
+    header.attr("names") = names;
344 
+
345 
+    return header;
346 
+  }
352 347 
   Rf_warningcall(R_NilValue, "pwiz not yet initialized.");
353 348 
   return Rcpp::DataFrame::create( );
354 349 
 }
355 350
356 
-Rcpp::DataFrame RcppPwiz::getAllScanHeaderInfo ( )
357 
-{
358 
-  if (msd != NULL)
359 
-    {
360 
-      if (!isInCacheAllScanHeaderInfo)
361 
-        {
362 
-	  SpectrumListPtr slp = msd->run.spectrumListPtr;
363 
-	  int N = slp->size();
364 
-
365 
-	  allScanHeaderInfo = getScanHeaderInfo(Rcpp::seq(1, N));
366 
-	  isInCacheAllScanHeaderInfo = TRUE;
367 
-        }
368 
-      return allScanHeaderInfo ;
351 
+Rcpp::DataFrame RcppPwiz::getAllScanHeaderInfo ( ) {
352 
+  if (msd != NULL) {
353 
+    if (!isInCacheAllScanHeaderInfo) {
354 
+      SpectrumListPtr slp = msd->run.spectrumListPtr;
355 
+      int N = slp->size();
356 
+
357 
+      allScanHeaderInfo = getScanHeaderInfo(Rcpp::seq(1, N));
358 
+      isInCacheAllScanHeaderInfo = TRUE;
369 359 
     }
360 
+    return allScanHeaderInfo ;
361 
+  }
370 362 
   Rf_warningcall(R_NilValue, "pwiz not yet initialized.");
371 363 
   return Rcpp::DataFrame::create( );
372 364 
 }
373 365
374 
-Rcpp::List RcppPwiz::getPeakList ( int whichScan )
375 
-{
376 
-  if (msd != NULL)
377 
-    {
378 
-      SpectrumListPtr slp = msd->run.spectrumListPtr;
379 
-
380 
-      if ((whichScan <= 0) || (whichScan > slp->size()))
381 
-        {
382 
-	  Rprintf("Index whichScan out of bounds [1 ... %d].\n", slp->size());
383 
-	  return Rcpp::List::create( );
384 
-        }
385 
-
386 
-      SpectrumPtr s = slp->spectrum(whichScan - 1, true);
387 
-      vector<MZIntensityPair> pairs;
388 
-      s->getMZIntensityPairs(pairs);
389 
-
390 
-      Rcpp::NumericMatrix peaks(pairs.size(), 2);
391 
-
392 
-      if(pairs.size()!=0)
393 
-        {
394 
-	  for (int i = 0; i < pairs.size(); i++)
395 
-            {
396 
-	      MZIntensityPair p = pairs.at(i);
397 
-	      peaks(i,0) = p.mz;
398 
-	      peaks(i,1) = p.intensity;
399 
-            }
400 
-
401 
-        }
402 
-
403 
-      return Rcpp::List::create(
404 
-				Rcpp::_["peaksCount"]  = pairs.size(),
405 
-				Rcpp::_["peaks"]  = peaks
406 
-				) ;
366 
+Rcpp::List RcppPwiz::getPeakList(Rcpp::IntegerVector whichScan) {
367 
+  if (msd != NULL) {
368 
+    SpectrumListPtr slp = msd->run.spectrumListPtr;
369 
+    int n_scans = slp->size();
370 
+    int n_want = whichScan.size();
371 
+    int current_scan;
372 
+    SpectrumPtr sp;
373 
+    BinaryDataArrayPtr mzs,ints;
374 
+    std::vector<double> data;
375 
+    Rcpp::NumericVector data_matrix;
376 
+    Rcpp::List res(n_want);
377 
+    for (int i = 0; i < n_want; i++) {
378 
+      current_scan = whichScan[i];
379 
+      if (current_scan < 1 || current_scan > n_scans) {
380 
+	Rprintf("Index whichScan out of bounds [1 ... %d].\n", n_scans);
381 
+	return Rcpp::List::create( );
382 
+      }
383 
+      sp = slp->spectrum(current_scan - 1, true);
384 
+      mzs = sp->getMZArray();
385 
+      ints = sp->getIntensityArray();
386 
+      if (!mzs.get() || !ints.get()) {
387 
+	Rcpp::NumericMatrix pks(0, 2);
388 
+	res[i] = pks;
389 
+	continue;
390 
+      }
391 
+      if (mzs->data.size() != ints->data.size())
392 
+	Rcpp::Rcerr << "Sizes of mz and intensity arrays don't match." << std::endl;
393 
+      data = mzs->data;
394 
+      data.insert(data.end(), ints->data.begin(), ints->data.end());
395 
+      data_matrix = Rcpp::wrap(data);
396 
+      data_matrix.attr("dim") = Rcpp::Dimension(ints->data.size(), 2);
397 
+      res[i] = data_matrix;
407 398 
     }
399 
+    return res;
400 
+  }
408 401 
   Rf_warningcall(R_NilValue, "pwiz not yet initialized.");
409 
-  return Rcpp::List::create( );
402 
+  return Rcpp::List::create();
410 403 
 }
411 404
412 405 
 /**
src/RcppPwiz.h
History View file @ b4a2f2dad2
... ... 
@@ -98,7 +98,7 @@ public:
98 98
99 99 
     Rcpp::DataFrame getAllChromatogramHeaderInfo();
100 100
101 
-    Rcpp::List getPeakList(int whichScan);
101 
+    Rcpp::List getPeakList(Rcpp::IntegerVector whichScan);
102 102
103 103 
     Rcpp::NumericMatrix get3DMap(std::vector<int> scanNumbers, double whichMzLow, double whichMzHigh, double resMz);
104 104