Browse code

restored old ramp Interface

Steffen Neumann authored on 25/06/2014 04:55:42
Showing 7 changed files

... ...
@@ -1,27 +1,27 @@
1
-Package: mzR
2
-Type: Package
3
-Title: parser for netCDF, mzXML, mzData and mzML files (mass
4
-        spectrometry data)
5
-Version: 1.11.3
6
-Author: Bernd Fischer, Steffen Neumann, Laurent Gatto, Qiang Kou
7
-Maintainer: Bernd Fischer <bernd.fischer@embl.de>, Steffen Neumann
8
- <sneumann@ipb-halle.de>, Laurent Gatto <lg390@cam.ac.uk>,
9
- Qiang Kou <qkou@umail.iu.edu>
10
-Description: mzR provides a unified API to the common file formats and
11
-        parsers available for mass spectrometry data. It comes with a
12
-        wrapper for the ISB random access parser for mass spectrometry
13
-        mzXML, mzData and mzML files. The package contains the original
14
-        code written by the ISB, and a subset of the proteowizard
15
-        library for mzML. The netCDF reading code has previously been
16
-        used in XCMS.
17
-License: Artistic-2.0
18
-LazyLoad: yes
19
-Depends: Rcpp (>= 0.10.1), methods, utils
20
-Imports: Biobase
21
-Suggests: msdata (>= 0.1.9), RUnit, faahKO
22
-LinkingTo: Rcpp
23
-RcppModules: Ramp
24
-SystemRequirements: GNU make, NetCDF, zlib
25
-biocViews: Infrastructure, DataImport, Proteomics, Metabolomics,
26
-        MassSpectrometry
27
-Packaged: 2014-06-21 21:42:58 UTC; Administrator
1
+Package: mzR
2
+Type: Package
3
+Title: parser for netCDF, mzXML, mzData and mzML files (mass
4
+        spectrometry data)
5
+Version: 1.11.4
6
+Author: Bernd Fischer, Steffen Neumann, Laurent Gatto, Qiang Kou
7
+Maintainer: Bernd Fischer <bernd.fischer@embl.de>, Steffen Neumann
8
+ <sneumann@ipb-halle.de>, Laurent Gatto <lg390@cam.ac.uk>,
9
+ Qiang Kou <qkou@umail.iu.edu>
10
+Description: mzR provides a unified API to the common file formats and
11
+        parsers available for mass spectrometry data. It comes with a
12
+        wrapper for the ISB random access parser for mass spectrometry
13
+        mzXML, mzData and mzML files. The package contains the original
14
+        code written by the ISB, and a subset of the proteowizard
15
+        library for mzML. The netCDF reading code has previously been
16
+        used in XCMS.
17
+License: Artistic-2.0
18
+LazyLoad: yes
19
+Depends: Rcpp (>= 0.10.1), methods, utils
20
+Imports: Biobase
21
+Suggests: msdata (>= 0.1.9), RUnit, faahKO
22
+LinkingTo: Rcpp
23
+RcppModules: Ramp
24
+SystemRequirements: GNU make, NetCDF, zlib
25
+biocViews: Infrastructure, DataImport, Proteomics, Metabolomics,
26
+        MassSpectrometry
27
+Packaged: 2014-06-25; sneumann
... ...
@@ -1,6 +1,6 @@
1
-CHANGES IN VERSION 1.11.3
1
+CHANGES IN VERSION 1.11.4
2 2
 -------------------------
3
- o Updated pwiz and fixed building on all platforms [2014-06-22 Sun]
3
+ o restored the old ramp Interface
4 4
 
5 5
 CHANGES IN VERSION 1.11.2
6 6
 -------------------------
7 7
new file mode 100644
... ...
@@ -0,0 +1,144 @@
1
+rampInit <- function() {
2
+
3
+    result <- .C("RampRInit",
4
+                 PACKAGE = "mzR")
5
+}
6
+
7
+rampPrintFiles <- function() {
8
+
9
+    result <- .C("RampRPrintFiles",
10
+                 PACKAGE = "mzR")
11
+}
12
+
13
+rampIsFile <- function(filename) {
14
+
15
+    # The C version doesn't do anything extra
16
+    #.C("RampRIsFile",
17
+    #   as.character(filename),
18
+    #   isfile = logical(1),
19
+    #   PACKAGE = "mzR")$isfile
20
+
21
+    if (!file.exists(filename))
22
+        return(FALSE)
23
+    text <- readChar(filename, 1024)
24
+
25
+    length(text) > 0
26
+}
27
+
28
+rampOpen <- function(filename) {
29
+
30
+    result <- .C("RampROpen",
31
+                 as.character(filename),
32
+                 rampid = integer(1),
33
+                 status = integer(1),
34
+                 PACKAGE = "mzR")
35
+
36
+    if (result$status)
37
+        return(result$status)
38
+
39
+    return(result$rampid)
40
+}
41
+
42
+rampClose <- function(rampid) {
43
+
44
+    result <- .C("RampRClose",
45
+                 as.integer(rampid),
46
+                 PACKAGE = "mzR")
47
+}
48
+
49
+rampCloseAll <- function() {
50
+
51
+    result <- .C("RampRCloseAll",
52
+                 PACKAGE = "mzR")
53
+}
54
+
55
+rampNumScans <- function(rampid) {
56
+
57
+    result <- .C("RampRNumScans",
58
+                 as.integer(rampid),
59
+                 numscans = integer(1),
60
+                 status = integer(1),
61
+                 PACKAGE = "mzR")
62
+
63
+    if (result$status)
64
+        return(NA)
65
+
66
+    return(result$numscans)
67
+}
68
+
69
+rampScanHeaders <- function(rampid) {
70
+
71
+    .Call("RampRScanHeaders",
72
+          as.integer(rampid),
73
+          PACKAGE = "mzR")
74
+}
75
+
76
+rampSIPeaks <- function(rampid, seqNum, peaksCount) {
77
+
78
+    if (!is.integer(seqNum))
79
+        seqNum <- as.integer(seqNum)
80
+    if (!is.integer(peaksCount))
81
+        peaksCount <- as.integer(peaksCount)
82
+    .Call("RampRSIPeaks",
83
+          as.integer(rampid),
84
+          seqNum,
85
+          peaksCount,
86
+          PACKAGE = "mzR")
87
+}
88
+
89
+rampRawData <- function(rampid) {
90
+
91
+    scanHeaders <- rampScanHeaders(rampid)
92
+
93
+    # Some of these checks work around buggy RAMP indexing code
94
+    scans <- scanHeaders$msLevel == 1 & scanHeaders$seqNum > 0 &
95
+             !duplicated(scanHeaders$acquisitionNum) &
96
+             scanHeaders$peaksCount > 0
97
+    if ("Full" %in% levels(scanHeaders$scanType))
98
+        scans <- scans & scanHeaders$scanType == "Full"
99
+
100
+    scans <- which(scans)
101
+
102
+    sipeaks <- rampSIPeaks(rampid, scans, scanHeaders$peaksCount[scans])
103
+
104
+    return(list(rt = scanHeaders$retentionTime[scans],
105
+                acquisitionNum = scanHeaders$acquisitionNum[scans],
106
+                tic = scanHeaders$totIonCurrent[scans],
107
+                scanindex = sipeaks$scanindex,
108
+                mz = sipeaks$mz,
109
+                intensity = sipeaks$intensity,
110
+                polarity = scanHeaders$polarity[scans]))
111
+}
112
+
113
+rampRawDataMSn <- function(rampid) {
114
+
115
+    # Check if we have MSn at all
116
+    scanHeaders <- rampScanHeaders(rampid)
117
+    if (max(scanHeaders[,"msLevel"]) < 2) {
118
+        warning("MSn spectra requested but not found")
119
+        return (NULL);
120
+    }
121
+
122
+    # Some of these checks work around buggy RAMP indexing code
123
+    scans <- ( scanHeaders$msLevel >= 2 & scanHeaders$seqNum > 0
124
+              & !duplicated(scanHeaders$acquisitionNum)
125
+              & scanHeaders$peaksCount > 0)
126
+
127
+    scans <- which(scans)
128
+
129
+    sipeaks <- rampSIPeaks(rampid, scans, scanHeaders$peaksCount[scans])
130
+
131
+    retdata <- list(rt = scanHeaders$retentionTime[scans],
132
+                    acquisitionNum = scanHeaders$acquisitionNum[scans],
133
+                    precursorNum=scanHeaders$precursorScanNum[scans],
134
+                    precursorMZ = scanHeaders$precursorMZ[scans],
135
+                    precursorIntensity = scanHeaders$precursorIntensity[scans],
136
+                    peaksCount=scanHeaders$peaksCount[scans],
137
+                    msLevel = scanHeaders$msLevel[scans],
138
+                    precursorCharge = scanHeaders$precursorCharge[scans],
139
+                    scanindex = sipeaks$scanindex, collisionEnergy = scanHeaders$collisionEnergy[scans],
140
+                    mz = sipeaks$mz,
141
+                    intensity =sipeaks$intensity);
142
+
143
+    return(retdata)
144
+}
0 145
new file mode 100644
... ...
@@ -0,0 +1,63 @@
1
+test.oldramp.mzXML <- function() {
2
+    library(msdata)
3
+    library(mzR)
4
+    cdfpath <- system.file("threonine", package = "msdata")
5
+    filename <- list.files(cdfpath, pattern="threonine_i2_e35_pH_tree.mzXML",
6
+                       full.names=TRUE, recursive = TRUE)
7
+
8
+    rampid <- mzR:::rampOpen(filename)
9
+    if (rampid < 0)
10
+       stop("Could not open mzXML/mzData file")
11
+
12
+    on.exit(mzR:::rampClose(rampid))
13
+    rawdata <- mzR:::rampRawData(rampid)
14
+    mzR:::rampClose(rampid)
15
+}
16
+
17
+test.oldramp.mzML <- function() {
18
+    library(msdata)
19
+    library(mzR)
20
+    cdfpath <- system.file("microtofq", package = "msdata")
21
+    filename <- list.files(cdfpath, pattern="MM14.mzML",
22
+                       full.names=TRUE, recursive = TRUE)
23
+
24
+    rampid <- mzR:::rampOpen(filename)
25
+    if (rampid < 0)
26
+       stop("Could not open mzXML/mzData file")
27
+
28
+    on.exit(mzR:::rampClose(rampid))
29
+    rawdata <- mzR:::rampRawData(rampid)
30
+    mzR:::rampClose(rampid)
31
+}
32
+
33
+test.oldramp.mzData <- function() {
34
+    library(msdata)
35
+    library(mzR)
36
+    cdfpath <- system.file("microtofq", package = "msdata")
37
+    filename <- list.files(cdfpath, pattern="MM14.mzdata$",
38
+                       full.names=TRUE, recursive = TRUE)
39
+
40
+    rampid <- mzR:::rampOpen(filename)
41
+    if (rampid < 0)
42
+       stop("Could not open mzXML/mzData file")
43
+
44
+    on.exit(mzR:::rampClose(rampid))
45
+    rawdata <- mzR:::rampRawData(rampid)
46
+    mzR:::rampClose(rampid)
47
+}
48
+
49
+test.oldramp.mzData.gz <- function() {
50
+    library(msdata)
51
+    library(mzR)
52
+    cdfpath <- system.file("microtofq", package = "msdata")
53
+    filename <- list.files(cdfpath, pattern="MM14.mzdata.gz",
54
+                       full.names=TRUE, recursive = TRUE)
55
+
56
+    rampid <- mzR:::rampOpen(filename)
57
+    if (rampid < 0)
58
+       stop("Could not open mzXML/mzData file")
59
+
60
+    on.exit(mzR:::rampClose(rampid))
61
+    rawdata <- mzR:::rampRawData(rampid)
62
+    mzR:::rampClose(rampid)
63
+}
... ...
@@ -73,7 +73,7 @@ PWIZOBJECTS=\
73 73
 
74 74
 MZROBJECTS=cramp.o  ramp_base64.o  ramp.o  RcppRamp.o RcppRampModule.o rnetCDF.o 
75 75
 
76
-OBJECTS= $(MZROBJECTS) $(PWIZOBJECTS) R_init_mzR.o
76
+OBJECTS= $(MZROBJECTS) $(PWIZOBJECTS) rampR.o R_init_mzR.o
77 77
 
78 78
 ##
79 79
 ## R complains about assert(), hence -D_NODEBUG
... ...
@@ -77,7 +77,7 @@ PWIZOBJECTS=\
77 77
 
78 78
 MZROBJECTS=cramp.o  ramp_base64.o  ramp.o  RcppRamp.o RcppRampModule.o rnetCDF.o 
79 79
 
80
-OBJECTS= $(MZROBJECTS) $(PWIZOBJECTS) R_init_mzR.o
80
+OBJECTS= $(MZROBJECTS) $(PWIZOBJECTS) rampR.o R_init_mzR.o
81 81
 
82 82
 PWIZ_CPPFLAGS=-I./boost_aux/ -I./boost -I. -DHAVE_PWIZ_MZML_LIB -fpermissive -DWINDOWS_NATIVE
83 83
 
84 84
new file mode 100644
... ...
@@ -0,0 +1,414 @@
1
+#include<Rcpp.h>
2
+
3
+#ifdef __MINGW32__
4
+#undef Realloc
5
+#undef Free
6
+//#include <winsock2.h>
7
+#endif
8
+
9
+
10
+#include <stdio.h>
11
+
12
+#include <R.h>
13
+#include <Rinternals.h>
14
+
15
+#include "ramp.h"
16
+
17
+#define MAX_RAMP_FILES 100
18
+#define NEW_LIST(n)		Rf_allocVector(VECSXP,n)
19
+#define NEW_CHARACTER(n)	Rf_allocVector(STRSXP,n)
20
+#define SET_NAMES(x, n)		Rf_setAttrib(x, R_NamesSymbol, n)
21
+#define SET_CLASS(x, n)		Rf_setAttrib(x, R_ClassSymbol, n)
22
+#define INTEGER_POINTER(x)	INTEGER(x)
23
+#define NEW_INTEGER(n)		Rf_allocVector(INTSXP,n)
24
+#define NEW_NUMERIC(n)		Rf_allocVector(REALSXP,n)
25
+#define NUMERIC_POINTER(x)	REAL(x)
26
+#define SET_LEVELS(x, l)       	Rf_setAttrib(x, R_LevelsSymbol, l)
27
+extern"C" {
28
+
29
+typedef struct {
30
+    RAMPFILE          *file;
31
+    ramp_fileoffset_t *index;
32
+    int               numscans;
33
+} RAMPSTRUCT;
34
+
35
+static int        rampInitalized = 0;
36
+static RAMPSTRUCT rampStructs[MAX_RAMP_FILES];
37
+
38
+void RampRInit(void) {
39
+
40
+    int    i;
41
+    
42
+    for (i = 0; i < MAX_RAMP_FILES; i++) {
43
+        rampStructs[i].file = NULL;
44
+        rampStructs[i].index = NULL;
45
+        rampStructs[i].numscans = 0;
46
+    }
47
+    
48
+    rampInitalized = 1;
49
+}
50
+
51
+void RampRPrintFiles(void) {
52
+
53
+    int i;
54
+    
55
+    if (!rampInitalized)
56
+        RampRInit();
57
+    
58
+    for (i = 0; i < MAX_RAMP_FILES; i++)
59
+        if (rampStructs[i].file)
60
+            Rprintf("File %i (%i scans)\n", i, rampStructs[i].numscans);
61
+}
62
+
63
+int RampRFreeHandle(void) {
64
+
65
+    int i;
66
+    
67
+    if (!rampInitalized)
68
+        RampRInit();
69
+    
70
+    for (i = 0; i < MAX_RAMP_FILES; i++)
71
+        if (!rampStructs[i].file)
72
+            return i;
73
+    
74
+    return -1;
75
+}
76
+
77
+void RampRIsFile(const char *fileName[], int *isfile) {
78
+
79
+    RAMPFILE *rampfile;
80
+    
81
+    *isfile = 0;
82
+    
83
+    rampfile = rampOpenFile(fileName[0]);
84
+    if (!rampfile)
85
+	    return;
86
+	
87
+	rampCloseFile(rampfile);
88
+	
89
+	*isfile = 1;
90
+}
91
+
92
+void RampROpen(const char *fileName[], int *rampid, int *status) {
93
+
94
+    ramp_fileoffset_t indexOffset;
95
+    int               numscans;
96
+    
97
+    if (!rampInitalized)
98
+        RampRInit();
99
+    
100
+    *status = -1;
101
+    
102
+    *rampid = RampRFreeHandle();
103
+    if (*rampid < 0) {
104
+        *status = *rampid;
105
+        return;
106
+    }
107
+    
108
+    rampStructs[*rampid].file = rampOpenFile(fileName[0]);
109
+	if (!rampStructs[*rampid].file)
110
+	    return;
111
+	
112
+	indexOffset = getIndexOffset(rampStructs[*rampid].file);
113
+	
114
+	rampStructs[*rampid].index = readIndex(rampStructs[*rampid].file,
115
+	                                       indexOffset, &numscans);
116
+	
117
+	if (!rampStructs[*rampid].index || numscans < 1) {
118
+	    rampStructs[*rampid].file = NULL;
119
+	    if (rampStructs[*rampid].index)
120
+	        free(rampStructs[*rampid].index);
121
+	    rampStructs[*rampid].index = NULL;
122
+	    return;
123
+	}
124
+	
125
+	rampStructs[*rampid].numscans = numscans;
126
+	*status = 0;
127
+}
128
+
129
+void RampRClose(const int *rampid) {
130
+
131
+    if (!rampInitalized)
132
+        return;
133
+
134
+    if (*rampid < 0 || *rampid >= MAX_RAMP_FILES)
135
+        return;
136
+
137
+    if (rampStructs[*rampid].file)
138
+        rampCloseFile(rampStructs[*rampid].file);
139
+    rampStructs[*rampid].file = NULL;
140
+    
141
+    if (rampStructs[*rampid].index)
142
+        free(rampStructs[*rampid].index);
143
+    rampStructs[*rampid].index = NULL;
144
+    
145
+    rampStructs[*rampid].numscans = 0;
146
+}
147
+
148
+void RampRCloseAll(void) {
149
+
150
+    int i;
151
+    
152
+    if (!rampInitalized)
153
+        return;
154
+    
155
+    for (i = 0; i < MAX_RAMP_FILES; i++)
156
+        if (rampStructs[i].file)
157
+            RampRClose(&i);
158
+}
159
+
160
+void RampRNumScans(const int *rampid, int *numscans, int *status) {
161
+
162
+    if (!rampInitalized)
163
+        RampRInit();
164
+
165
+    *status = -1;
166
+
167
+    if (*rampid < 0 || *rampid >= MAX_RAMP_FILES)
168
+        return;
169
+	
170
+	*numscans = rampStructs[*rampid].numscans;
171
+	
172
+	if (*numscans)
173
+	    *status = 0;
174
+}
175
+
176
+SEXP RampRScanHeaders(SEXP rampid) {
177
+
178
+    int               i, j, id, numscans, ncol = 18, ntypes = 0, stlen = 10;
179
+    SEXP              result = PROTECT(NEW_LIST(ncol));
180
+	SEXP              temp;
181
+	SEXP              names;
182
+    struct            ScanHeaderStruct scanHeader;
183
+    RAMPFILE          *file;
184
+    ramp_fileoffset_t *index;
185
+    char              rowname[20], *scanTypes;
186
+    int               *seqNum, *acquisitionNum, *msLevel, *peaksCount,
187
+      *precursorScanNum, *precursorCharge, *scanType;
188
+    double            *totIonCurrent, *retentionTime, *basePeakMZ, 
189
+                      *basePeakIntensity, *collisionEnergy, *ionisationEnergy,
190
+	              *lowMZ, *highMZ, *precursorMZ, *precursorIntensity;
191
+    
192
+    if (!rampInitalized)
193
+        RampRInit();
194
+    
195
+    if (Rf_length(rampid) != 1)
196
+        Rf_error("rampid must be of length 1");
197
+    
198
+    id = *INTEGER_POINTER(rampid);
199
+    if (id < 0 || id >= MAX_RAMP_FILES || !rampStructs[id].file)
200
+        Rf_error("invalid rampid");
201
+    
202
+    file = rampStructs[id].file;
203
+    index = rampStructs[id].index;
204
+    numscans = rampStructs[id].numscans;
205
+    
206
+    SET_NAMES(result, names = NEW_CHARACTER(ncol));
207
+    
208
+    Rf_setAttrib(result, Rf_install("row.names"), temp = NEW_CHARACTER(numscans));
209
+    for (i = 0; i < numscans; i++) {
210
+        sprintf(rowname, "%i", i+1);
211
+        SET_STRING_ELT(temp, i, Rf_mkChar(rowname));
212
+    }
213
+    
214
+    SET_CLASS(result, temp = NEW_CHARACTER(1));
215
+    SET_STRING_ELT(temp, 0, Rf_mkChar("data.frame"));
216
+    
217
+    SET_VECTOR_ELT(result, 0, temp = NEW_INTEGER(numscans));
218
+    seqNum = INTEGER_POINTER(temp);
219
+    SET_STRING_ELT(names, 0, Rf_mkChar("seqNum"));
220
+    
221
+    SET_VECTOR_ELT(result, 1, temp = NEW_INTEGER(numscans));
222
+    acquisitionNum = INTEGER_POINTER(temp);
223
+    SET_STRING_ELT(names, 1, Rf_mkChar("acquisitionNum"));
224
+    
225
+    SET_VECTOR_ELT(result, 2, temp = NEW_INTEGER(numscans));
226
+    msLevel = INTEGER_POINTER(temp);
227
+    SET_STRING_ELT(names, 2, Rf_mkChar("msLevel"));
228
+    
229
+    SET_VECTOR_ELT(result, 3, temp = NEW_INTEGER(numscans));
230
+    peaksCount = INTEGER_POINTER(temp);
231
+    SET_STRING_ELT(names, 3, Rf_mkChar("peaksCount"));
232
+    
233
+    SET_VECTOR_ELT(result, 4, temp = NEW_NUMERIC(numscans));
234
+    totIonCurrent = NUMERIC_POINTER(temp);
235
+    SET_STRING_ELT(names, 4, Rf_mkChar("totIonCurrent"));
236
+    
237
+    SET_VECTOR_ELT(result, 5, temp = NEW_NUMERIC(numscans));
238
+    retentionTime = NUMERIC_POINTER(temp);
239
+    SET_STRING_ELT(names, 5, Rf_mkChar("retentionTime"));
240
+    
241
+    SET_VECTOR_ELT(result, 6, temp = NEW_NUMERIC(numscans));
242
+    basePeakMZ = NUMERIC_POINTER(temp);
243
+    SET_STRING_ELT(names, 6, Rf_mkChar("basePeakMZ"));
244
+    
245
+    SET_VECTOR_ELT(result, 7, temp = NEW_NUMERIC(numscans));
246
+    basePeakIntensity = NUMERIC_POINTER(temp);
247
+    SET_STRING_ELT(names, 7, Rf_mkChar("basePeakIntensity"));
248
+    
249
+    SET_VECTOR_ELT(result, 8, temp = NEW_NUMERIC(numscans));
250
+    collisionEnergy = NUMERIC_POINTER(temp);
251
+    SET_STRING_ELT(names, 8, Rf_mkChar("collisionEnergy"));
252
+    
253
+    SET_VECTOR_ELT(result, 9, temp = NEW_NUMERIC(numscans));
254
+    ionisationEnergy = NUMERIC_POINTER(temp);
255
+    SET_STRING_ELT(names, 9, Rf_mkChar("ionisationEnergy"));
256
+    
257
+    SET_VECTOR_ELT(result, 10, temp = NEW_NUMERIC(numscans));
258
+    lowMZ = NUMERIC_POINTER(temp);
259
+    SET_STRING_ELT(names, 10, Rf_mkChar("lowMZ"));
260
+    
261
+    SET_VECTOR_ELT(result, 11, temp = NEW_NUMERIC(numscans));
262
+    highMZ = NUMERIC_POINTER(temp);
263
+    SET_STRING_ELT(names, 11, Rf_mkChar("highMZ"));
264
+    
265
+    SET_VECTOR_ELT(result, 12, temp = NEW_INTEGER(numscans));
266
+    precursorScanNum = INTEGER_POINTER(temp);
267
+    SET_STRING_ELT(names, 12, Rf_mkChar("precursorScanNum"));
268
+    
269
+    SET_VECTOR_ELT(result, 13, temp = NEW_NUMERIC(numscans));
270
+    precursorMZ = NUMERIC_POINTER(temp);
271
+    SET_STRING_ELT(names, 13, Rf_mkChar("precursorMZ"));
272
+    
273
+    SET_VECTOR_ELT(result, 14, temp = NEW_INTEGER(numscans));
274
+    precursorCharge = INTEGER_POINTER(temp);
275
+    SET_STRING_ELT(names, 14, Rf_mkChar("precursorCharge"));
276
+    
277
+    SET_VECTOR_ELT(result, 15, temp = NEW_INTEGER(numscans));
278
+    scanType = INTEGER_POINTER(temp);
279
+    SET_STRING_ELT(names, 15, Rf_mkChar("scanType"));
280
+
281
+    SET_VECTOR_ELT(result, 16, temp = NEW_NUMERIC(numscans));
282
+    precursorIntensity = NUMERIC_POINTER(temp);
283
+    SET_STRING_ELT(names, 16, Rf_mkChar("precursorIntensity"));
284
+
285
+    //SET_VECTOR_ELT(result, 17, temp = NEW_INTEGER(numscans));
286
+    //polarity = INTEGER_POINTER(temp);
287
+    //SET_STRING_ELT(names, 17, Rf_mkChar("polarity"));
288
+        
289
+    scanTypes = S_alloc(stlen*SCANTYPE_LENGTH, sizeof(char));
290
+    
291
+    for (i = 0; i < numscans; i++) {
292
+	    readHeader(file, index[i+1], &scanHeader);
293
+	    seqNum[i] = scanHeader.seqNum;
294
+	    acquisitionNum[i] = scanHeader.acquisitionNum;
295
+	    msLevel[i] = scanHeader.msLevel;
296
+	    peaksCount[i] = scanHeader.peaksCount;
297
+        totIonCurrent[i] = scanHeader.totIonCurrent;
298
+        retentionTime[i] = scanHeader.retentionTime;
299
+        basePeakMZ[i] = scanHeader.basePeakMZ;
300
+        basePeakIntensity[i] = scanHeader.basePeakIntensity;
301
+        collisionEnergy[i] = scanHeader.collisionEnergy;
302
+        ionisationEnergy[i] = scanHeader.ionisationEnergy;
303
+        lowMZ[i] = scanHeader.lowMZ;
304
+        highMZ[i] = scanHeader.highMZ;
305
+        precursorScanNum[i] = scanHeader.precursorScanNum;
306
+        precursorMZ[i] = scanHeader.precursorMZ;
307
+        precursorIntensity[i] = scanHeader.precursorIntensity;
308
+	//polarity[i] = scanHeader.polarity;
309
+        precursorCharge[i] = scanHeader.precursorCharge;
310
+        for (j = 0; j < ntypes; j++)
311
+            if (!strcmp(scanHeader.scanType, scanTypes+j*SCANTYPE_LENGTH)) {
312
+                scanType[i] = j+1;
313
+                break;
314
+            }
315
+        if (j == ntypes) {
316
+            if (j >= stlen) {
317
+                stlen *= 2;
318
+                scanTypes = S_realloc(scanTypes, stlen*SCANTYPE_LENGTH/2,
319
+                                      stlen*SCANTYPE_LENGTH, sizeof(char));
320
+            }
321
+            strcpy(scanTypes+j*SCANTYPE_LENGTH, scanHeader.scanType);
322
+            ntypes++;
323
+            scanType[i] = j+1;
324
+        }
325
+    }
326
+    
327
+    SET_LEVELS(VECTOR_ELT(result, 16), temp = NEW_CHARACTER(ntypes));
328
+    for (i = 0; i < ntypes; i++)
329
+        SET_STRING_ELT(temp, i, Rf_mkChar(scanTypes+i*SCANTYPE_LENGTH));
330
+    
331
+    SET_CLASS(VECTOR_ELT(result, 15), temp = NEW_CHARACTER(1));
332
+    SET_STRING_ELT(temp, 0, Rf_mkChar("factor"));
333
+    
334
+    UNPROTECT(1);
335
+    
336
+    return result;
337
+}
338
+
339
+SEXP RampRSIPeaks(SEXP rampid, SEXP seqNum, SEXP peaksCount) {
340
+
341
+    SEXP              result = PROTECT(NEW_LIST(3)), temp, names;
342
+    RAMPFILE          *file;
343
+    ramp_fileoffset_t *index;
344
+    int               i, j, id, numscans, numpeaks, *seqNumPtr, *peaksCountPtr,
345
+                      *scanindex;
346
+    double            *mz, *intensity;
347
+    RAMPREAL          *peaks, *peaksPtr;
348
+    
349
+    if (!rampInitalized)
350
+        RampRInit();
351
+    
352
+    if (Rf_length(rampid) != 1)
353
+        Rf_error("rampid must be of length 1");
354
+    
355
+    if (Rf_length(seqNum) != Rf_length(peaksCount))
356
+        Rf_error("seqNum and peaksCount must be the same length");
357
+    
358
+    id = *INTEGER_POINTER(rampid);
359
+    if (id < 0 || id >= MAX_RAMP_FILES || !rampStructs[id].file)
360
+        Rf_error("invalid rampid");
361
+    
362
+    file = rampStructs[id].file;
363
+    index = rampStructs[id].index;
364
+    
365
+    seqNumPtr = INTEGER_POINTER(seqNum);
366
+    peaksCountPtr = INTEGER_POINTER(peaksCount);
367
+    numscans = Rf_length(seqNum);
368
+    
369
+    SET_NAMES(result, names = NEW_CHARACTER(3));
370
+    
371
+    SET_VECTOR_ELT(result, 0, temp = NEW_INTEGER(numscans));
372
+    scanindex = INTEGER_POINTER(temp);
373
+    SET_STRING_ELT(names, 0, Rf_mkChar("scanindex"));
374
+    
375
+    numpeaks = 0;
376
+    for (i = 0; i < numscans; i++) {
377
+        if (seqNumPtr[i] > rampStructs[id].numscans)
378
+            Rf_error("invalid number in seqnum");
379
+        scanindex[i] = numpeaks;
380
+        numpeaks += peaksCountPtr[i];
381
+    }
382
+    
383
+    SET_VECTOR_ELT(result, 1, temp = NEW_NUMERIC(numpeaks));
384
+    mz = NUMERIC_POINTER(temp);
385
+    SET_STRING_ELT(names, 1, Rf_mkChar("mz"));
386
+    
387
+    SET_VECTOR_ELT(result, 2, temp = NEW_NUMERIC(numpeaks));
388
+    intensity = NUMERIC_POINTER(temp);
389
+    SET_STRING_ELT(names, 2, Rf_mkChar("intensity"));
390
+    
391
+    for (i = 0; i < numscans; i++) {
392
+        if (peaksCountPtr[i] != readPeaksCount(file, index[seqNumPtr[i]]))
393
+            Rf_error("invalid number in peaksCount");
394
+        if (peaksCountPtr[i] == 0)
395
+            continue;
396
+        peaks = readPeaks(file, index[seqNumPtr[i]]);
397
+        if (!peaks)
398
+            Rf_error("unknown problem while reading peaks");
399
+        peaksPtr = peaks;
400
+        for (j = 0; j < peaksCountPtr[i]; j++) {
401
+            if (*peaksPtr < 0)
402
+                Rf_error("unexpected end of peak list");
403
+            mz[j+scanindex[i]] = *(peaksPtr++);
404
+		    intensity[j+scanindex[i]] = *(peaksPtr++);
405
+        }
406
+        free(peaks);
407
+    }
408
+    
409
+    UNPROTECT(1);
410
+    
411
+    return result;
412
+}
413
+
414
+}