mzxml_parsing.cpp
来自「MS-Clustering is designed to rapidly clu」· C++ 代码 · 共 1,002 行 · 第 1/2 页
CPP
1,002 行
{
BytesRead = fread(Buffer + BufferEnd, sizeof(char), BytesToRead, MZXMLFile);
BufferEnd += BytesRead;
}
Buffer[BufferEnd] = '\0';
// Look for a new <scan tag opening:
ScanStr = strstr(Buffer, "<scan");
if (ScanStr)
{
Pos = ScanStr - Buffer;
}
else
{
Pos = 0;
}
if (!ScanStr || Pos > XML_BUFFER_HALF_SIZE)
{
// There's not a <scan tag in the first half of the buffer.
// If we're at EOF, then stop now:
if (BufferEnd < XML_BUFFER_HALF_SIZE)
{
break;
}
// Shunt the tail of the buffer to the front, and carry on:
memmove(Buffer, Buffer + XML_BUFFER_HALF_SIZE, BufferEnd - XML_BUFFER_HALF_SIZE);
BufferEnd -= XML_BUFFER_HALF_SIZE;
FilePos += XML_BUFFER_HALF_SIZE;
continue;
}
ScanNumberStr = strstr(ScanStr, "num=");
if (!ScanNumberStr)
{
printf("** Warning: mzXML parser encountered a scan with no scan number! File %s Pos %d\n",
mzxml_name.c_str(), FilePos + Pos);
ScanNumber = -1;
}
else
{
ScanNumber = ParseIntFromXML(ScanNumberStr);
}
retentionTimeStr = strstr(ScanStr,"retentionTime=\"PT");
if (! retentionTimeStr)
{
// printf("Error: mzXML parser encountered a scan with no retnetion time: File %s Pos %d\n",
// mzxml_name.c_str(), FilePos + Pos);
// exit(1);
retentionTime = -1;
}
else
{
retentionTime = ParseMassFromXML(retentionTimeStr);
}
char *PeakCountStr = strstr(ScanStr, "peaksCount=\"");
if (!PeakCountStr)
{
cout << "Warning: couldn't parse peaks from mzxml! " << endl;
cout << "Scan: " << ScanNumber << " Pos: " << Pos << endl;
continue;
}
int PeakCount = ParseIntFromXML(PeakCountStr);
MSLevelStr = strstr(ScanStr, "msLevel=");
if (!MSLevelStr)
{
printf("** Warning: mzXML parser encountered a scan with no MS level! File %s Scan %d Pos %d\n",
mzxml_name.c_str(), ScanNumber, FilePos + Pos);
MSLevel = -1;
continue;
}
else
{
MSLevel = ParseIntFromXML(MSLevelStr);
}
precursorIntensityStr = strstr(ScanStr,"precursorIntensity=");
if (! precursorIntensityStr)
{
//if (MSLevel>1)
//{
// printf("Error: mzXML parser encountered a scan with no precursor intenisty: File %s Pos %d\n",
// mzxml_name.c_str(), FilePos + Pos);
// exit(1);
//}
precursorIntensity = 0;
}
else
{
precursorIntensity = ParseMassFromXML(precursorIntensityStr);
}
PrecursorStr = strstr(ScanStr, "<precursorMz");
if (PrecursorStr)
{
PrecursorStr = strstr(PrecursorStr, ">");
precursorMZ = ParseMassFromXML(PrecursorStr);
}
if (!PrecursorStr && MSLevel > 1)
{
printf("Warning: mzXML parser encountered a scan with no m/z: File %s Pos %d\n",
mzxml_name.c_str(), FilePos + Pos);
continue;
}
if (MSLevel > 1 && ScanNumber >= 0 && PeakCount > 10)
{
spec_header.scan_number = ScanNumber;
spec_header.MS_level = MSLevel;
spec_header.precursor_intensity = precursorIntensity;
spec_header.retention_time = retentionTime;
spec_header.m_over_z = precursorMZ;
spec_header.num_peaks = PeakCount;
spec_header.charge = 0;
spec_header.type = MZXML;
spec_header.file_pos = FilePos + Pos;
single_spectra.push_back(spec_header);
}
}
free(Buffer);
fclose(MZXMLFile);
// cout << single_spectra.size() << " spectra..." << endl;
}
/******************************************************************************
This is a special function designed to overcome parsing problems I have
with mzXML in Linux enviornments. The function serially extracts spectra from
an mzXML file and stores the peak lists (floats of pairs (mass,intensity)
*******************************************************************************/
int MZXML_file::extract_peak_lists_from_mzXML(Config *config,
string& mzxml_name,
int file_idx,
mass_t min_m_over_z,
mass_t max_m_over_z)
{
static char* PeakBuffer = NULL;
static char* DecodedPeakBuffer = NULL;
static float* Peaks = NULL;
static float* FilteredPeaks = NULL;
static int PeakBufferSize = 0;
int BytesToRead;
char* Buffer;
int BufferStartPos = 0;
int Pos;
int BytesRead;
int BufferEnd = 0;
FILE* MZXMLFile;
int ParseState = 0;
int FilePos = 0;
int ByteOrderLittle = 1;
char* ScanStr;
char* ScanNumberStr;
char* PeakStr;
char* MSLevelStr;
char *retentionTimeStr;
char *precursorIntensityStr;
char *PrecursorStr;
char* ByteOrderStr;
int ScanNumber;
int MSLevel;
float retentionTime;
float precursorIntensity;
mass_t precursorMZ;
//
Buffer = (char*)calloc(XML_BUFFER_SIZE + 1, sizeof(char));
MZXMLFile = fopen(mzxml_name.c_str(), "rb");
if (!MZXMLFile)
{
cout << "Error: Can't open MZXML file " << mzxml_name << endl;
exit(1);
}
// cout << "Extracting peaks from: " << mzxml_name << endl;
// initialize the file_peak_buff. Initialize to file_size*0.08 floats
int file_size = getFileSize(mzxml_name.c_str());
int num_floats = (int)(file_size*0.08);
if (num_floats<20000)
num_floats=20000;
if (file_peak_buff.size()<(int)(num_floats*1.5))
file_peak_buff.resize((int)(num_floats*1.5));
int spec_counter=0;
char *scan_start_ptr = NULL;
while (1)
{
MZXML_single spec_header;
spec_header.file_idx = file_idx;
// check if the peak buff should be extended
if (file_peak_buff.size()-file_peak_buff_pos<5000)
{
file_peak_buff.resize((int)(file_peak_buff.size()*1.5));
}
// Read more data, to fill up the buffer:
if ( ! scan_start_ptr ||
( (Buffer + BufferEnd - scan_start_ptr) < XML_BUFFER_HALF_SIZE) )
{
// try shunt half of the buffer
if (scan_start_ptr)
{
if (BufferEnd - XML_BUFFER_HALF_SIZE>0)
{
memmove(Buffer, Buffer + XML_BUFFER_HALF_SIZE, BufferEnd - XML_BUFFER_HALF_SIZE);
BufferEnd -= XML_BUFFER_HALF_SIZE;
scan_start_ptr -= XML_BUFFER_HALF_SIZE;
// cout << "MOVED!" << endl;
}
}
else
scan_start_ptr = Buffer;
BytesToRead = XML_BUFFER_SIZE - BufferEnd;
BytesRead = fread(Buffer + BufferEnd, sizeof(char), BytesToRead, MZXMLFile);
if (BytesRead<5)
break;
BufferEnd += BytesRead;
Buffer[BufferEnd] = '\0';
FilePos += BytesRead;
}
// Look for a new <scan tag opening:
const char *last_pos = Buffer + BufferEnd - 5;
char *pos = scan_start_ptr;
while (++pos<last_pos)
{
if (*pos != '<')
continue;
if (*(pos+1)=='s' && *(pos+2)=='c' && *(pos+3)=='a' && *(pos+4)=='n')
break;
}
ScanStr = (pos<last_pos) ? pos : NULL;
if (ScanStr)
{
Pos = ScanStr - Buffer;
}
else
{
Pos = 0;
}
if (!ScanStr )
{
scan_start_ptr = Buffer + BufferEnd-5;
continue;
}
ScanNumberStr = strstr(ScanStr, "num=");
if (!ScanNumberStr)
{
// printf("** Warning: mzXML parser encountered a scan with no scan number! File %s Pos %d\n",
// mzxml_name.c_str(), FilePos + Pos);
ScanNumber = -1;
}
else
{
ScanNumber = ParseIntFromXML(ScanNumberStr);
}
retentionTimeStr = strstr(ScanStr,"retentionTime=\"PT");
if (! retentionTimeStr)
{
// printf("Error: mzXML parser encountered a scan with no retnetion time: File %s Pos %d\n",
// mzxml_name.c_str(), FilePos + Pos);
// exit(1);
retentionTime = -1;
}
else
{
retentionTime = ParseMassFromXML(retentionTimeStr);
}
char *PeakCountStr = strstr(ScanStr, "peaksCount=\"");
if (!PeakCountStr)
{
cout << "Warninig: bad parsing of peaks from mzxml! " << endl;
cout << "Scan: " << ScanNumber << " Pos: " << FilePos << endl;
scan_start_ptr += 50;
continue;
}
int PeakCount = ParseIntFromXML(PeakCountStr);
MSLevelStr = strstr(ScanStr, "msLevel=");
if (!MSLevelStr)
{
printf("** Warning: mzXML parser encountered a scan with no MS level! File %s Pos %d\n",
mzxml_name.c_str(), FilePos + Pos);
scan_start_ptr += 50;
continue;
}
else
{
MSLevel = ParseIntFromXML(MSLevelStr);
}
precursorIntensityStr = strstr(ScanStr,"precursorIntensity=");
if (! precursorIntensityStr)
{
scan_start_ptr += 50;
continue;
// if (MSLevel>1)
// {
// printf("Warning: mzXML parser encountered a scan with no precursor intenisty: File %s Pos %d\n",
// mzxml_name.c_str(), FilePos + Pos);
// exit(1);
// }
}
else
{
precursorIntensity = ParseMassFromXML(precursorIntensityStr);
}
PrecursorStr = strstr(ScanStr, "<precursorMz");
if (PrecursorStr)
{
PrecursorStr = strstr(PrecursorStr, ">");
precursorMZ = ParseMassFromXML(PrecursorStr);
}
if (!PrecursorStr && MSLevel > 1)
{
printf("Warning: mzXML parser encountered a scan with no m/z: File %s Pos %d\n",
mzxml_name.c_str(), FilePos + Pos);
scan_start_ptr += 50;
continue;
}
if (MSLevel > 1 && ScanNumber >= 0 && PeakCount > 7 &&
precursorMZ>=min_m_over_z && precursorMZ<=max_m_over_z)
{
// read peaks
PeakStr = strstr(PrecursorStr, "<peaks");
if (PeakStr)
{
// Get byte order:
ByteOrderStr = strstr(PeakStr, "byteOrder=\"");
if (ByteOrderStr)
{
ByteOrderStr += 11;
if (!strncmp(ByteOrderStr, "network", 7))
{
ByteOrderLittle = 0;
}
if (!strncmp(ByteOrderStr, "big", 3))
{
ByteOrderLittle = 0;
}
if (!strncmp(ByteOrderStr, "little", 6))
{
ByteOrderLittle = 1;
}
}
PeakStr = strstr(PeakStr, ">");
}
if (!PeakStr)
{
cout << "Warning: bad parsing of peaks from mzxml (scan " << ScanNumber << ") skipping..." << endl;
scan_start_ptr += 50;
continue;
}
PeakStr++;
PeakBuffer = PeakStr;
if (PeakBufferSize < PeakCount)
{
if (DecodedPeakBuffer)
{
char *dbf = DecodedPeakBuffer;
free(DecodedPeakBuffer);
DecodedPeakBuffer = NULL;
free(Peaks);
Peaks = NULL;
free(FilteredPeaks);
FilteredPeaks=NULL;
}
PeakBufferSize = (int)(PeakCount*1.5);
DecodedPeakBuffer = (char*)calloc(PeakBufferSize * 8 + 8, 1);
Peaks = (float*)calloc(PeakBufferSize * 2, sizeof(float));
FilteredPeaks = (float*)calloc(PeakBufferSize * 2, sizeof(float));
}
int Trail = (PeakCount % 3);
if (!(PeakCount % 3))
{
PeakBuffer[PeakCount * 32/3] = '\0';
}
else
{
PeakBuffer[(PeakCount * 32/3) + Trail + 1] = '\0';
}
b64_decode_mio( DecodedPeakBuffer, PeakBuffer);
int FloatIndex;
for (FloatIndex = 0; FloatIndex < (2 * PeakCount); FloatIndex++)
{
#ifdef BYTEORDER_LITTLE_ENDIAN
if (!ByteOrderLittle)
{
char ByteSwap = DecodedPeakBuffer[FloatIndex*4];
DecodedPeakBuffer[FloatIndex*4] = DecodedPeakBuffer[FloatIndex*4 + 3];
DecodedPeakBuffer[FloatIndex*4 + 3] = ByteSwap;
ByteSwap = DecodedPeakBuffer[FloatIndex*4 + 1];
DecodedPeakBuffer[FloatIndex*4 + 1] = DecodedPeakBuffer[FloatIndex*4 + 2];
DecodedPeakBuffer[FloatIndex*4 + 2] = ByteSwap;
}
memcpy(Peaks + FloatIndex, DecodedPeakBuffer + FloatIndex * 4, 4);
#else
if (ByteOrderLittle)
{
char ByteSwap = DecodedPeakBuffer[FloatIndex*4];
DecodedPeakBuffer[FloatIndex*4] = DecodedPeakBuffer[FloatIndex*4 + 3];
DecodedPeakBuffer[FloatIndex*4 + 3] = ByteSwap;
ByteSwap = DecodedPeakBuffer[FloatIndex*4 + 1];
DecodedPeakBuffer[FloatIndex*4 + 1] = DecodedPeakBuffer[FloatIndex*4 + 2];
DecodedPeakBuffer[FloatIndex*4 + 2] = ByteSwap;
}
memcpy(Peaks + FloatIndex, DecodedPeakBuffer + FloatIndex * 4, 4);
#endif
}
int num_new_peaks = join_and_filter_peak_list(config,precursorMZ,Peaks,
PeakCount, &file_peak_buff[file_peak_buff_pos]);
spec_header.scan_number = ScanNumber;
spec_header.MS_level = MSLevel;
spec_header.precursor_intensity = precursorIntensity;
spec_header.retention_time = retentionTime;
spec_header.m_over_z = precursorMZ;
spec_header.num_peaks = num_new_peaks;
spec_header.charge = 0;
spec_header.type = MZXML;
spec_header.file_pos = FilePos + Pos;
spec_header.peak_buff_start_idx = file_peak_buff_pos;
single_spectra.push_back(spec_header);
spec_counter++;
file_peak_buff_pos += num_new_peaks * 2;
scan_start_ptr = PeakStr + 8 * PeakCount;
}
else
scan_start_ptr = ScanStr +50;
}
free(Buffer);
fclose(MZXMLFile);
// cout << spec_counter << " spectra, " << file_peak_buff_pos << "/" << file_peak_buff.size() << endl;
return spec_counter;
}
⌨️ 快捷键说明
复制代码Ctrl + C
搜索代码Ctrl + F
全屏模式F11
增大字号Ctrl + =
减小字号Ctrl + -
显示快捷键?