source: trunk/Mars/mcorsika/MCorsikaFormat.cc@ 10060

Last change on this file since 10060 was 10060, checked in by rohlfs, 11 years ago
two new command line arguments of readcorsika: -A=arrayNo and -T=telescopeNo. New design of program flow in MCorsikaRead: It is now determined by the order of the data blocks in the input file.
File size: 6.7 KB
Line 
1/* ======================================================================== *\
2!
3! *
4! * This file is part of MARS, the MAGIC Analysis and Reconstruction
5! * Software. It is distributed to you in the hope that it can be a useful
6! * and timesaving tool in analysing Data of imaging Cerenkov telescopes.
7! * It is distributed WITHOUT ANY WARRANTY.
8! *
9! * Permission to use, copy, modify and distribute this software and its
10! * documentation for any purpose is hereby granted without fee,
11! * provided that the above copyright notice appear in all copies and
12! * that both that copyright notice and this permission notice appear
13! * in supporting documentation. It is provided "as is" without express
14! * or implied warranty.
15! *
16!
17!
18! Author(s): Reiner Rohlfs 2010
19! Author(s): Thomas Bretz 2010 <mailto:thomas.bretz@epfl.ch>
20!
21! Copyright: Software Development, 2000-2010
22!
23!
24\* ======================================================================== */
25
26//////////////////////////////////////////////////////////////////////////////
27//
28// MCorsikaFormat
29//
30//////////////////////////////////////////////////////////////////////////////
31#include "MCorsikaFormat.h"
32
33#include <errno.h>
34#include <fstream>
35
36#include "MLogManip.h"
37
38
39using namespace std;
40
41
42const unsigned int MCorsikaFormat::kSyncMarker = 0xd41f8a37;
43
44// --------------------------------------------------------------------------
45//
46MCorsikaFormat *MCorsikaFormat::CorsikaFormatFactory(const char * fileName)
47{
48 ifstream * fileIn = new ifstream(fileName);
49
50 const Bool_t noexist = !(*fileIn);
51 if (noexist)
52 {
53 gLog << err << "Cannot open file " << fileName << ": ";
54 gLog << (errno!=0?strerror(errno):"Insufficient memory for decompression") << endl;
55 delete fileIn;
56 return NULL;
57 }
58
59 char *buffer = new char[5];
60 memset(buffer, 0, 5);
61 fileIn->read(buffer, 4);
62 fileIn->seekg(-4, ios::cur);
63
64 if (strcmp(buffer, "RUNH") == 0)
65 {
66 delete [] buffer;
67 return new MCorsikaFormatRaw(fileIn);
68 }
69
70 if (*reinterpret_cast<unsigned int*>(buffer) == kSyncMarker)
71 {
72 delete [] buffer;
73 return new MCorsikaFormatEventIO(fileIn);
74 }
75
76 gLog << err << "File " << fileName <<
77 " is neither a CORSIKA raw nor EventIO file" << endl;
78 delete fileIn;
79 delete [] buffer;
80
81 return NULL;
82}
83
84Bool_t MCorsikaFormat::Read(void *ptr, Int_t i) const
85{
86 fIn->read((char*)ptr, i);
87 return !fIn->fail();
88
89}
90// --------------------------------------------------------------------------
91//
92Bool_t MCorsikaFormat::Eof() const
93{
94 return fIn->eof();
95}
96
97// --------------------------------------------------------------------------
98//
99MCorsikaFormat::~MCorsikaFormat()
100{
101 delete fIn;
102}
103
104
105// --------------------------------------------------------------------------
106//
107// After a call to this function, the file pointer is located after the
108// header of the block. As the event block has no header it is located
109// at the beginning of this block, which is as the beginning of the data
110Bool_t MCorsikaFormatRaw::NextBlock(Bool_t subBlock,
111 Int_t & blockType,
112 Int_t & blockVersion,
113 Int_t & blockIdentifier,
114 Int_t & blockLength) const
115{
116 char blockHeader[5];
117 blockHeader[4] = 0;
118 fIn->read(blockHeader, 4);
119 if (fIn->eof())
120 return kFALSE;
121
122 blockVersion = 0;
123 blockIdentifier = 0;
124 blockLength = 272 * 4;
125
126 if (strcmp(blockHeader, "RUNH") == 0)
127 blockType = 1200;
128 else if (strcmp(blockHeader, "RUNE") == 0)
129 blockType = 1210;
130 else if (strcmp(blockHeader, "EVTH") == 0)
131 blockType = 1202;
132 else if (strcmp(blockHeader, "EVTE") == 0)
133 blockType = 1209;
134 else // the events, they don't have a specific header
135 {
136 blockType = 1105;
137 fIn->seekg(-4, ios::cur);
138 blockLength += 4;
139 }
140
141 return kTRUE;
142}
143// --------------------------------------------------------------------------
144//
145Bool_t MCorsikaFormatRaw::SeekEvtEnd()
146{
147 // Search subblockwise backward (Block: 5733*4 = 21*273*4)
148 for (int i=1; i<22; i++)
149 {
150 fIn->seekg(-i*273*4, ios::end);
151
152 char runh[5]="\0\0\0\0";
153 fIn->read(runh, 4);
154
155 if (!strcmp(runh, "RUNE"))
156 {
157// fIn->seekg(-4, ios::cur);
158 return kTRUE;
159 }
160 }
161
162 return kTRUE;
163}
164///////////////////////////////////////////////////////////////////////////////
165///////////////////////////////////////////////////////////////////////////////
166
167Bool_t MCorsikaFormatEventIO::NextBlock(Bool_t subBlock,
168 Int_t & blockType,
169 Int_t & blockVersion,
170 Int_t & blockIdentifier,
171 Int_t & blockLength) const
172{
173// we read - synchronisation markerif not subBlock
174// - type / version field
175// - identification field
176// - length
177// - unknown field
178// - id (first 4 bytes of data field)
179
180 if (subBlock)
181 {
182 int blockHeader[3];
183 fIn->read((char*)blockHeader, 3 * sizeof(int));
184
185 if (fIn->eof())
186 return kFALSE;
187
188
189 blockType = blockHeader[0] & 0xFFFF;
190 blockVersion = (blockHeader[0] & 0xFFF00000) >> 20;
191 blockIdentifier = blockHeader[1];
192 blockLength = blockHeader[2] & 0x3FFFFFFF;
193 }
194 else
195 {
196 int blockHeader[4];
197 fIn->read((char*)blockHeader, 4 * sizeof(int));
198
199 if (fIn->eof())
200 return kFALSE;
201
202
203 blockType = blockHeader[1] & 0xFFFF;
204 blockVersion = (blockHeader[1] & 0xFFF00000) >> 20;
205 blockIdentifier = blockHeader[2];
206 blockLength = blockHeader[3] & 0x3FFFFFFF;
207
208 if (blockType == 1200 || blockType == 1210 ||
209 blockType == 1202 || blockType == 1209 )
210 {
211 // read the "old" corsika header plus additional 4 unknown byte
212 char tmp[8];
213 fIn->read(tmp, 8);
214 blockLength -= 8;
215 }
216
217 }
218 return kTRUE;
219}
220
221// --------------------------------------------------------------------------
222//
223Bool_t MCorsikaFormatEventIO::SeekEvtEnd()
224{
225
226 // the RUNE block it at the very end of the file.
227 fIn->seekg(-32, ios::end);
228
229 unsigned int blockHeader[4];
230 fIn->read((char*)blockHeader, 4 * sizeof(int));
231
232 if ( blockHeader[0] == kSyncMarker &&
233 (blockHeader[1] & 0xffff) == 1210 &&
234 (blockHeader[3] & 0x3fffffff) == 16)
235 {
236 // this seams to be a RUNE (1210) block
237 fIn->seekg( 8, ios::cur);
238 return kTRUE;
239 }
240
241 return kFALSE;
242}
243
Note: See TracBrowser for help on using the repository browser.