Changeset 1716933 in git
- Timestamp:
- May 25, 2018, 8:56:33 AM (7 years ago)
- Branches:
- ImprovedOutputFile, Timing, dual_readout, llp, master
- Children:
- 7d0eb75
- Parents:
- 792092a
- Files:
-
- 4 added
- 8 edited
Legend:
- Unmodified
- Added
- Removed
-
Makefile
r792092a r1716933 15 15 PcmSuf = _rdict.pcm 16 16 17 CXXFLAGS += $(ROOTCFLAGS) -Wno-write-strings -D_FILE_OFFSET_BITS=64 -DDROP_CGAL -I. -Iexternal -Iexternal/tcl -I/usr/include/tirpc17 CXXFLAGS += $(ROOTCFLAGS) -Wno-write-strings -D_FILE_OFFSET_BITS=64 -DDROP_CGAL -I. -Iexternal -Iexternal/tcl 18 18 DELPHES_LIBS = $(shell $(RC) --libs) -lEG $(SYSLIBS) 19 19 DISPLAY_LIBS = $(shell $(RC) --evelibs) -lGuiHtml $(SYSLIBS) … … 481 481 tmp/classes/DelphesPileUpReader.$(ObjSuf): \ 482 482 classes/DelphesPileUpReader.$(SrcSuf) \ 483 classes/DelphesPileUpReader.h 483 classes/DelphesPileUpReader.h \ 484 classes/DelphesXDRReader.h 484 485 tmp/classes/DelphesPileUpWriter.$(ObjSuf): \ 485 486 classes/DelphesPileUpWriter.$(SrcSuf) \ 486 classes/DelphesPileUpWriter.h 487 classes/DelphesPileUpWriter.h \ 488 classes/DelphesXDRWriter.h 487 489 tmp/classes/DelphesSTDHEPReader.$(ObjSuf): \ 488 490 classes/DelphesSTDHEPReader.$(SrcSuf) \ … … 490 492 classes/DelphesClasses.h \ 491 493 classes/DelphesFactory.h \ 494 classes/DelphesXDRReader.h \ 492 495 external/ExRootAnalysis/ExRootTreeBranch.h 493 496 tmp/classes/DelphesStream.$(ObjSuf): \ … … 497 500 classes/DelphesTF2.$(SrcSuf) \ 498 501 classes/DelphesTF2.h 502 tmp/classes/DelphesXDRReader.$(ObjSuf): \ 503 classes/DelphesXDRReader.$(SrcSuf) \ 504 classes/DelphesXDRReader.h 505 tmp/classes/DelphesXDRWriter.$(ObjSuf): \ 506 classes/DelphesXDRWriter.$(SrcSuf) \ 507 classes/DelphesXDRWriter.h 499 508 tmp/external/ExRootAnalysis/ExRootConfReader.$(ObjSuf): \ 500 509 external/ExRootAnalysis/ExRootConfReader.$(SrcSuf) \ … … 1004 1013 tmp/classes/DelphesStream.$(ObjSuf) \ 1005 1014 tmp/classes/DelphesTF2.$(ObjSuf) \ 1015 tmp/classes/DelphesXDRReader.$(ObjSuf) \ 1016 tmp/classes/DelphesXDRWriter.$(ObjSuf) \ 1006 1017 tmp/external/ExRootAnalysis/ExRootConfReader.$(ObjSuf) \ 1007 1018 tmp/external/ExRootAnalysis/ExRootFilter.$(ObjSuf) \ … … 2018 2029 @touch $@ 2019 2030 2031 classes/DelphesSTDHEPReader.h: \ 2032 classes/DelphesXDRReader.h 2033 @touch $@ 2034 2020 2035 external/fastjet/plugins/CDFCones/fastjet/CDFMidPointPlugin.hh: \ 2021 2036 external/fastjet/JetDefinition.hh -
classes/DelphesPileUpReader.cc
r792092a r1716933 33 33 34 34 #include <stdio.h> 35 #include <rpc/types.h> 36 #include <rpc/xdr.h> 35 #include <stdint.h> 36 37 #include "classes/DelphesXDRReader.h" 37 38 38 39 using namespace std; … … 47 48 fEntries(0), fEntrySize(0), fCounter(0), 48 49 fPileUpFile(0), fIndex(0), fBuffer(0), 49 fInput XDR(0), fIndexXDR(0), fBufferXDR(0)50 fInputReader(0), fIndexReader(0), fBufferReader(0) 50 51 { 51 52 stringstream message; 52 53 53 fIndex = new char[kIndexSize*8]; 54 fBuffer = new char[kBufferSize*kRecordSize*4]; 55 fInputXDR = new XDR; 56 fIndexXDR = new XDR; 57 fBufferXDR = new XDR; 58 xdrmem_create(fIndexXDR, fIndex, kIndexSize*8, XDR_DECODE); 59 xdrmem_create(fBufferXDR, fBuffer, kBufferSize*kRecordSize*4, XDR_DECODE); 54 fIndex = new uint8_t[kIndexSize*8]; 55 fBuffer = new uint8_t[kBufferSize*kRecordSize*4]; 56 fInputReader = new DelphesXDRReader; 57 fIndexReader = new DelphesXDRReader; 58 fBufferReader = new DelphesXDRReader; 59 60 fIndexReader->SetBuffer(fIndex); 61 fBufferReader->SetBuffer(fBuffer); 60 62 61 63 fPileUpFile = fopen(fileName, "r"); … … 67 69 } 68 70 69 xdrstdio_create(fInputXDR, fPileUpFile, XDR_DECODE);71 fInputReader->SetFile(fPileUpFile); 70 72 71 73 // read number of events 72 74 fseeko(fPileUpFile, -8, SEEK_END); 73 xdr_hyper(fInputXDR, &fEntries);75 fInputReader->ReadValue(&fEntries, 8); 74 76 75 77 if(fEntries >= kIndexSize) … … 81 83 // read index of events 82 84 fseeko(fPileUpFile, -8 - 8*fEntries, SEEK_END); 83 xdr_opaque(fInputXDR,fIndex, fEntries*8);85 fInputReader->ReadRaw(fIndex, fEntries*8); 84 86 } 85 87 … … 88 90 DelphesPileUpReader::~DelphesPileUpReader() 89 91 { 90 xdr_destroy(fInputXDR);91 92 if(fPileUpFile) fclose(fPileUpFile); 92 xdr_destroy(fBufferXDR); 93 xdr_destroy(fIndexXDR); 94 if(fBufferXDR) delete fBufferXDR; 95 if(fIndexXDR) delete fIndexXDR; 96 if(fInputXDR) delete fInputXDR; 93 if(fBufferReader) delete fBufferReader; 94 if(fIndexReader) delete fIndexReader; 95 if(fInputReader) delete fInputReader; 97 96 if(fBuffer) delete[] fBuffer; 98 97 if(fIndex) delete[] fIndex; … … 101 100 //------------------------------------------------------------------------------ 102 101 103 bool DelphesPileUpReader::ReadParticle(int &pid,102 bool DelphesPileUpReader::ReadParticle(int32_t &pid, 104 103 float &x, float &y, float &z, float &t, 105 104 float &px, float &py, float &pz, float &e) … … 107 106 if(fCounter >= fEntrySize) return false; 108 107 109 xdr_int(fBufferXDR, &pid);110 xdr_float(fBufferXDR, &x);111 xdr_float(fBufferXDR, &y);112 xdr_float(fBufferXDR, &z);113 xdr_float(fBufferXDR, &t);114 xdr_float(fBufferXDR, &px);115 xdr_float(fBufferXDR, &py);116 xdr_float(fBufferXDR, &pz);117 xdr_float(fBufferXDR, &e);108 fBufferReader->ReadValue(&pid, 4); 109 fBufferReader->ReadValue(&x, 4); 110 fBufferReader->ReadValue(&y, 4); 111 fBufferReader->ReadValue(&z, 4); 112 fBufferReader->ReadValue(&t, 4); 113 fBufferReader->ReadValue(&px, 4); 114 fBufferReader->ReadValue(&py, 4); 115 fBufferReader->ReadValue(&pz, 4); 116 fBufferReader->ReadValue(&e, 4); 118 117 119 118 ++fCounter; … … 124 123 //------------------------------------------------------------------------------ 125 124 126 bool DelphesPileUpReader::ReadEntry( quad_t entry)125 bool DelphesPileUpReader::ReadEntry(int64_t entry) 127 126 { 128 quad_t offset;127 int64_t offset; 129 128 130 129 if(entry >= fEntries) return false; 131 130 132 131 // read event position 133 xdr_setpos(fIndexXDR,8*entry);134 xdr_hyper(fIndexXDR, &offset);132 fIndexReader->SetOffset(8*entry); 133 fIndexReader->ReadValue(&offset, 8); 135 134 136 135 // read event 137 136 fseeko(fPileUpFile, offset, SEEK_SET); 138 xdr_int(fInputXDR, &fEntrySize);137 fInputReader->ReadValue(&fEntrySize, 4); 139 138 140 139 if(fEntrySize >= kBufferSize) … … 143 142 } 144 143 145 xdr_opaque(fInputXDR,fBuffer, fEntrySize*kRecordSize*4);146 xdr_setpos(fBufferXDR,0);144 fInputReader->ReadRaw(fBuffer, fEntrySize*kRecordSize*4); 145 fBufferReader->SetOffset(0); 147 146 fCounter = 0; 148 147 -
classes/DelphesPileUpReader.h
r792092a r1716933 29 29 30 30 #include <stdio.h> 31 #include <rpc/types.h> 32 #include <rpc/xdr.h> 31 #include <stdint.h> 32 33 class DelphesXDRReader; 33 34 34 35 class DelphesPileUpReader … … 40 41 ~DelphesPileUpReader(); 41 42 42 bool ReadParticle(int &pid,43 bool ReadParticle(int32_t &pid, 43 44 float &x, float &y, float &z, float &t, 44 45 float &px, float &py, float &pz, float &e); 45 46 46 bool ReadEntry( quad_t entry);47 bool ReadEntry(int64_t entry); 47 48 48 quad_t GetEntries() const { return fEntries; }49 int64_t GetEntries() const { return fEntries; } 49 50 50 51 private: 51 52 52 quad_t fEntries;53 int64_t fEntries; 53 54 54 int fEntrySize;55 int fCounter;55 int32_t fEntrySize; 56 int32_t fCounter; 56 57 57 58 FILE *fPileUpFile; 58 char*fIndex;59 char*fBuffer;59 uint8_t *fIndex; 60 uint8_t *fBuffer; 60 61 61 XDR *fInputXDR;62 XDR *fIndexXDR;63 XDR *fBufferXDR;62 DelphesXDRReader *fInputReader; 63 DelphesXDRReader *fIndexReader; 64 DelphesXDRReader *fBufferReader; 64 65 }; 65 66 66 67 #endif // DelphesPileUpReader_h 67 68 -
classes/DelphesPileUpWriter.cc
r792092a r1716933 33 33 34 34 #include <stdio.h> 35 #include <rpc/types.h> 36 #include <rpc/xdr.h> 35 #include <stdint.h> 36 37 #include "classes/DelphesXDRWriter.h" 37 38 38 39 using namespace std; … … 47 48 fEntries(0), fEntrySize(0), fOffset(0), 48 49 fPileUpFile(0), fIndex(0), fBuffer(0), 49 fOutput XDR(0), fIndexXDR(0), fBufferXDR(0)50 fOutputWriter(0), fIndexWriter(0), fBufferWriter(0) 50 51 { 51 52 stringstream message; 52 53 53 fIndex = new char[kIndexSize*8]; 54 fBuffer = new char[kBufferSize*kRecordSize*4]; 55 fOutputXDR = new XDR; 56 fIndexXDR = new XDR; 57 fBufferXDR = new XDR; 58 xdrmem_create(fIndexXDR, fIndex, kIndexSize*8, XDR_ENCODE); 59 xdrmem_create(fBufferXDR, fBuffer, kBufferSize*kRecordSize*4, XDR_ENCODE); 54 fIndex = new uint8_t[kIndexSize*8]; 55 fBuffer = new uint8_t[kBufferSize*kRecordSize*4]; 56 fOutputWriter = new DelphesXDRWriter; 57 fIndexWriter = new DelphesXDRWriter; 58 fBufferWriter = new DelphesXDRWriter; 59 60 fIndexWriter->SetBuffer(fIndex); 61 fBufferWriter->SetBuffer(fBuffer); 60 62 61 63 fPileUpFile = fopen(fileName, "w+"); … … 67 69 } 68 70 69 xdrstdio_create(fOutputXDR, fPileUpFile, XDR_ENCODE);71 fOutputWriter->SetFile(fPileUpFile); 70 72 } 71 73 … … 74 76 DelphesPileUpWriter::~DelphesPileUpWriter() 75 77 { 76 xdr_destroy(fOutputXDR);77 78 if(fPileUpFile) fclose(fPileUpFile); 78 xdr_destroy(fBufferXDR); 79 xdr_destroy(fIndexXDR); 80 if(fBufferXDR) delete fBufferXDR; 81 if(fIndexXDR) delete fIndexXDR; 82 if(fOutputXDR) delete fOutputXDR; 79 if(fBufferWriter) delete fBufferWriter; 80 if(fIndexWriter) delete fIndexWriter; 81 if(fOutputWriter) delete fOutputWriter; 83 82 if(fBuffer) delete[] fBuffer; 84 83 if(fIndex) delete[] fIndex; … … 87 86 //------------------------------------------------------------------------------ 88 87 89 void DelphesPileUpWriter::WriteParticle(int pid,88 void DelphesPileUpWriter::WriteParticle(int32_t pid, 90 89 float x, float y, float z, float t, 91 90 float px, float py, float pz, float e) … … 96 95 } 97 96 98 xdr_int(fBufferXDR, &pid);99 xdr_float(fBufferXDR, &x);100 xdr_float(fBufferXDR, &y);101 xdr_float(fBufferXDR, &z);102 xdr_float(fBufferXDR, &t);103 xdr_float(fBufferXDR, &px);104 xdr_float(fBufferXDR, &py);105 xdr_float(fBufferXDR, &pz);106 xdr_float(fBufferXDR, &e);97 fBufferWriter->WriteValue(&pid, 4); 98 fBufferWriter->WriteValue(&x, 4); 99 fBufferWriter->WriteValue(&y, 4); 100 fBufferWriter->WriteValue(&z, 4); 101 fBufferWriter->WriteValue(&t, 4); 102 fBufferWriter->WriteValue(&px, 4); 103 fBufferWriter->WriteValue(&py, 4); 104 fBufferWriter->WriteValue(&pz, 4); 105 fBufferWriter->WriteValue(&e, 4); 107 106 108 107 ++fEntrySize; … … 118 117 } 119 118 120 xdr_int(fOutputXDR, &fEntrySize);121 xdr_opaque(fOutputXDR,fBuffer, fEntrySize*kRecordSize*4);119 fOutputWriter->WriteValue(&fEntrySize, 4); 120 fOutputWriter->WriteRaw(fBuffer, fEntrySize*kRecordSize*4); 122 121 123 xdr_hyper(fIndexXDR, &fOffset);122 fIndexWriter->WriteValue(&fOffset, 8); 124 123 fOffset += fEntrySize*kRecordSize*4 + 4; 125 124 126 xdr_setpos(fBufferXDR,0);125 fBufferWriter->SetOffset(0); 127 126 fEntrySize = 0; 128 127 129 128 ++fEntries; 130 129 } … … 134 133 void DelphesPileUpWriter::WriteIndex() 135 134 { 136 xdr_opaque(fOutputXDR,fIndex, fEntries*8);137 xdr_hyper(fOutputXDR, &fEntries);135 fOutputWriter->WriteRaw(fIndex, fEntries*8); 136 fOutputWriter->WriteValue(&fEntries, 8); 138 137 } 139 138 -
classes/DelphesPileUpWriter.h
r792092a r1716933 29 29 30 30 #include <stdio.h> 31 #include <rpc/types.h> 32 #include <rpc/xdr.h> 31 #include <stdint.h> 32 33 class DelphesXDRWriter; 33 34 34 35 class DelphesPileUpWriter … … 40 41 ~DelphesPileUpWriter(); 41 42 42 void WriteParticle(int pid,43 void WriteParticle(int32_t pid, 43 44 float x, float y, float z, float t, 44 45 float px, float py, float pz, float e); … … 50 51 private: 51 52 52 quad_t fEntries;53 int fEntrySize;54 quad_t fOffset;53 int64_t fEntries; 54 int32_t fEntrySize; 55 int64_t fOffset; 55 56 56 57 FILE *fPileUpFile; 57 char*fIndex;58 char*fBuffer;58 uint8_t *fIndex; 59 uint8_t *fBuffer; 59 60 60 XDR *fOutputXDR;61 XDR *fIndexXDR;62 XDR *fBufferXDR;61 DelphesXDRWriter *fOutputWriter; 62 DelphesXDRWriter *fIndexWriter; 63 DelphesXDRWriter *fBufferWriter; 63 64 }; 64 65 65 66 #endif // DelphesPileUpWriter_h 66 67 -
classes/DelphesSTDHEPReader.cc
r792092a r1716933 34 34 #include <stdio.h> 35 35 #include <errno.h> 36 #include < rpc/types.h>37 #include < rpc/xdr.h>36 #include <stdint.h> 37 #include <string.h> 38 38 39 39 #include "TObjArray.h" … … 45 45 #include "classes/DelphesClasses.h" 46 46 #include "classes/DelphesFactory.h" 47 #include "classes/DelphesXDRReader.h" 47 48 48 49 #include "ExRootAnalysis/ExRootTreeBranch.h" … … 55 56 56 57 DelphesSTDHEPReader::DelphesSTDHEPReader() : 57 fInputFile(0), fInputXDR(0), fBuffer(0), fPDG(0), fBlockType(-1) 58 { 59 fInputXDR = new XDR; 60 fBuffer = new char[kBufferSize*96 + 24]; 58 fInputFile(0), fBuffer(0), fPDG(0), fBlockType(-1) 59 { 60 fBuffer = new uint8_t[kBufferSize*96 + 24]; 61 61 62 62 fPDG = TDatabasePDG::Instance(); … … 68 68 { 69 69 if(fBuffer) delete fBuffer; 70 if(fInputXDR) delete fInputXDR;71 70 } 72 71 … … 76 75 { 77 76 fInputFile = inputFile; 78 xdrstdio_create(fInputXDR, inputFile, XDR_DECODE);77 fReader[0].SetFile(inputFile); 79 78 } 80 79 … … 100 99 TObjArray *partonOutputArray) 101 100 { 101 fReader[0].ReadValue(&fBlockType, 4); 102 102 103 if(feof(fInputFile)) return kFALSE; 103 104 xdr_int(fInputXDR, &fBlockType);105 104 106 105 SkipBytes(4); … … 146 145 //--------------------------------------------------------------------------- 147 146 148 void DelphesSTDHEPReader::SkipBytes( u_int size)147 void DelphesSTDHEPReader::SkipBytes(int size) 149 148 { 150 149 int rc; 151 u_int rndup;150 int rndup; 152 151 153 152 rndup = size % 4; … … 161 160 if(rc != 0 && errno == ESPIPE) 162 161 { 163 xdr_opaque(fInputXDR,fBuffer, size);164 } 165 } 166 167 //--------------------------------------------------------------------------- 168 169 void DelphesSTDHEPReader::SkipArray( u_int elsize)170 { 171 u _int size;172 xdr_u_int(fInputXDR, &size);162 fReader[0].ReadRaw(fBuffer, size); 163 } 164 } 165 166 //--------------------------------------------------------------------------- 167 168 void DelphesSTDHEPReader::SkipArray(int elsize) 169 { 170 uint32_t size; 171 fReader[0].ReadValue(&size, 4); 173 172 SkipBytes(size*elsize); 174 173 } … … 178 177 void DelphesSTDHEPReader::ReadFileHeader() 179 178 { 180 u _int i;179 uint32_t i; 181 180 enum STDHEPVersion {UNKNOWN, V1, V2, V21} version; 182 181 183 182 // version 184 xdr_string(fInputXDR, &fBuffer, 100);183 fReader[0].ReadString(fBuffer, 100); 185 184 if(fBuffer[0] == '\0' || fBuffer[1] == '\0') version = UNKNOWN; 186 185 else if(fBuffer[0] == '1') version = V1; 187 else if(strncmp( fBuffer, "2.01", 4) == 0) version = V21;186 else if(strncmp((char *)fBuffer, "2.01", 4) == 0) version = V21; 188 187 else if(fBuffer[0] == '2') version = V2; 189 188 else version = UNKNOWN; … … 207 206 208 207 // Number of events 209 xdr_u_int(fInputXDR, &fEntries);208 fReader[0].ReadValue(&fEntries, 4); 210 209 211 210 SkipBytes(8); 212 211 213 212 // Number of blocks 214 u _int nBlocks = 0;215 xdr_u_int(fInputXDR, &nBlocks);213 uint32_t nBlocks = 0; 214 fReader[0].ReadValue(&nBlocks, 4); 216 215 217 216 // Number of NTuples 218 u _int nNTuples = 0;217 uint32_t nNTuples = 0; 219 218 if(version != V1) 220 219 { 221 xdr_u_int(fInputXDR, &nNTuples);220 fReader[0].ReadValue(&nNTuples, 4); 222 221 } 223 222 … … 244 243 { 245 244 // version 246 xdr_string(fInputXDR, &fBuffer, 100);247 if(strncmp( fBuffer, "1.00", 4) == 0)245 fReader[0].ReadString(fBuffer, 100); 246 if(strncmp((char *)fBuffer, "1.00", 4) == 0) 248 247 { 249 248 SkipBytes(8); … … 255 254 SkipArray(4); 256 255 } 257 else if(strncmp( fBuffer, "2.00", 4) == 0)256 else if(strncmp((char *)fBuffer, "2.00", 4) == 0) 258 257 { 259 258 SkipBytes(12); … … 272 271 { 273 272 bool skipNTuples = false; 274 u_int skipSize = 4;273 int skipSize = 4; 275 274 276 275 // version 277 xdr_string(fInputXDR, &fBuffer, 100);278 if(strncmp( fBuffer, "2.00", 4) == 0)276 fReader[0].ReadString(fBuffer, 100); 277 if(strncmp((char *)fBuffer, "2.00", 4) == 0) 279 278 { 280 279 skipNTuples = true; 281 280 } 282 else if(strncmp( fBuffer, "3.00", 4) == 0)281 else if(strncmp((char *)fBuffer, "3.00", 4) == 0) 283 282 { 284 283 skipNTuples = true; … … 288 287 SkipBytes(20); 289 288 290 u _int dimBlocks = 0;291 xdr_u_int(fInputXDR, &dimBlocks);292 293 u _int dimNTuples = 0;289 uint32_t dimBlocks = 0; 290 fReader[0].ReadValue(&dimBlocks, 4); 291 292 uint32_t dimNTuples = 0; 294 293 if(skipNTuples) 295 294 { 296 295 SkipBytes(4); 297 xdr_u_int(fInputXDR, &dimNTuples);296 fReader[0].ReadValue(&dimNTuples, 4); 298 297 } 299 298 … … 318 317 { 319 318 // version 320 xdr_string(fInputXDR, &fBuffer, 100);319 fReader[0].ReadString(fBuffer, 100); 321 320 322 321 // skip 5*4 + 2*8 = 36 bytes 323 322 SkipBytes(36); 324 323 325 if((strncmp( fBuffer, "1.", 2) == 0) || (strncmp(fBuffer, "2.", 2) == 0) ||326 (strncmp( fBuffer, "3.", 2) == 0) || (strncmp(fBuffer, "4.", 2) == 0) ||327 (strncmp( fBuffer, "5.00", 4) == 0))324 if((strncmp((char *)fBuffer, "1.", 2) == 0) || (strncmp((char *)fBuffer, "2.", 2) == 0) || 325 (strncmp((char *)fBuffer, "3.", 2) == 0) || (strncmp((char *)fBuffer, "4.", 2) == 0) || 326 (strncmp((char *)fBuffer, "5.00", 4) == 0)) 328 327 { 329 328 return; … … 333 332 SkipArray(1); 334 333 335 if(strncmp( fBuffer, "5.01", 4) == 0)334 if(strncmp((char *)fBuffer, "5.01", 4) == 0) 336 335 { 337 336 return; … … 345 344 void DelphesSTDHEPReader::ReadSTDHEP() 346 345 { 347 u _int idhepSize, isthepSize, jmohepSize, jdahepSize, phepSize, vhepSize;346 uint32_t idhepSize, isthepSize, jmohepSize, jdahepSize, phepSize, vhepSize; 348 347 349 348 // version 350 xdr_string(fInputXDR, &fBuffer, 100);349 fReader[0].ReadString(fBuffer, 100); 351 350 352 351 // Extracting the event number 353 xdr_int(fInputXDR, &fEventNumber);352 fReader[0].ReadValue(&fEventNumber, 4); 354 353 355 354 // Extracting the number of particles 356 xdr_int(fInputXDR, &fEventSize);355 fReader[0].ReadValue(&fEventSize, 4); 357 356 358 357 if(fEventSize >= kBufferSize) … … 364 363 // 4 + 4 + 4 + 4 + 4 + 4 = 96*n + 24 365 364 366 xdr_opaque(fInputXDR, fBuffer, 96*fEventSize + 24); 367 368 idhepSize = ntohl(*(u_int*)(fBuffer)); 369 isthepSize = ntohl(*(u_int*)(fBuffer + 4*1 + 4*1*fEventSize)); 370 jmohepSize = ntohl(*(u_int*)(fBuffer + 4*2 + 4*2*fEventSize)); 371 jdahepSize = ntohl(*(u_int*)(fBuffer + 4*3 + 4*4*fEventSize)); 372 phepSize = ntohl(*(u_int*)(fBuffer + 4*4 + 4*6*fEventSize)); 373 vhepSize = ntohl(*(u_int*)(fBuffer + 4*5 + 4*16*fEventSize)); 365 fReader[0].ReadRaw(fBuffer, 96*fEventSize + 24); 366 367 fReader[1].SetBuffer(fBuffer); 368 fReader[2].SetBuffer(fBuffer + 4*1 + 4*1*fEventSize); 369 fReader[3].SetBuffer(fBuffer + 4*2 + 4*2*fEventSize); 370 fReader[4].SetBuffer(fBuffer + 4*3 + 4*4*fEventSize); 371 fReader[5].SetBuffer(fBuffer + 4*4 + 4*6*fEventSize); 372 fReader[6].SetBuffer(fBuffer + 4*5 + 4*16*fEventSize); 373 374 fReader[1].ReadValue(&idhepSize, 4); 375 fReader[2].ReadValue(&isthepSize, 4); 376 fReader[3].ReadValue(&jmohepSize, 4); 377 fReader[4].ReadValue(&jdahepSize, 4); 378 fReader[5].ReadValue(&phepSize, 4); 379 fReader[6].ReadValue(&vhepSize, 4); 374 380 375 381 if(fEventSize < 0 || … … 392 398 void DelphesSTDHEPReader::ReadSTDHEP4() 393 399 { 394 u _int number;400 uint32_t number; 395 401 396 402 // Extracting the event weight 397 xdr_double(fInputXDR, &fWeight);403 fReader[0].ReadValue(&fWeight, 8); 398 404 399 405 // Extracting alpha QED 400 xdr_double(fInputXDR, &fAlphaQED);406 fReader[0].ReadValue(&fAlphaQED, 8); 401 407 402 408 // Extracting alpha QCD 403 xdr_double(fInputXDR, &fAlphaQCD);409 fReader[0].ReadValue(&fAlphaQCD, 8); 404 410 405 411 // Extracting the event scale 406 xdr_u_int(fInputXDR, &fScaleSize);412 fReader[0].ReadValue(&fScaleSize, 4); 407 413 for(number = 0; number < fScaleSize; ++number) 408 414 { 409 xdr_double(fInputXDR, &fScale[number]);415 fReader[0].ReadValue(&fScale[number], 8); 410 416 } 411 417 … … 450 456 451 457 int number; 452 int pid, status, m1, m2, d1, d2;458 int32_t pid, status, m1, m2, d1, d2; 453 459 double px, py, pz, e, mass; 454 460 double x, y, z, t; 455 461 456 XDR bufferXDR[6];457 xdrmem_create(&bufferXDR[0], fBuffer + 4*1, 4*fEventSize, XDR_DECODE);458 xdrmem_create(&bufferXDR[1], fBuffer + 4*2 + 4*1*fEventSize, 4*fEventSize, XDR_DECODE);459 xdrmem_create(&bufferXDR[2], fBuffer + 4*3 + 4*2*fEventSize, 8*fEventSize, XDR_DECODE);460 xdrmem_create(&bufferXDR[3], fBuffer + 4*4 + 4*4*fEventSize, 8*fEventSize, XDR_DECODE);461 xdrmem_create(&bufferXDR[4], fBuffer + 4*5 + 4*6*fEventSize, 40*fEventSize, XDR_DECODE);462 xdrmem_create(&bufferXDR[5], fBuffer + 4*6 + 4*16*fEventSize, 32*fEventSize, XDR_DECODE);463 464 462 for(number = 0; number < fEventSize; ++number) 465 463 { 466 xdr_int(&bufferXDR[0], &status);467 xdr_int(&bufferXDR[1], &pid);468 xdr_int(&bufferXDR[2], &m1);469 xdr_int(&bufferXDR[2], &m2);470 xdr_int(&bufferXDR[3], &d1);471 xdr_int(&bufferXDR[3], &d2);472 473 xdr_double(&bufferXDR[4], &px);474 xdr_double(&bufferXDR[4], &py);475 xdr_double(&bufferXDR[4], &pz);476 xdr_double(&bufferXDR[4], &e);477 xdr_double(&bufferXDR[4], &mass);478 479 xdr_double(&bufferXDR[5], &x);480 xdr_double(&bufferXDR[5], &y);481 xdr_double(&bufferXDR[5], &z);482 xdr_double(&bufferXDR[5], &t);464 fReader[1].ReadValue(&status, 4); 465 fReader[2].ReadValue(&pid, 4); 466 fReader[3].ReadValue(&m1, 4); 467 fReader[3].ReadValue(&m2, 4); 468 fReader[4].ReadValue(&d1, 4); 469 fReader[4].ReadValue(&d2, 4); 470 471 fReader[5].ReadValue(&px, 8); 472 fReader[5].ReadValue(&py, 8); 473 fReader[5].ReadValue(&pz, 8); 474 fReader[5].ReadValue(&e, 8); 475 fReader[5].ReadValue(&mass, 8); 476 477 fReader[6].ReadValue(&x, 8); 478 fReader[6].ReadValue(&y, 8); 479 fReader[6].ReadValue(&z, 8); 480 fReader[6].ReadValue(&t, 8); 483 481 484 482 candidate = factory->NewCandidate(); -
classes/DelphesSTDHEPReader.h
r792092a r1716933 29 29 30 30 #include <stdio.h> 31 #include <rpc/types.h> 32 #include <rpc/xdr.h> 31 #include <stdint.h> 32 33 #include "classes/DelphesXDRReader.h" 33 34 34 35 class TObjArray; … … 37 38 class ExRootTreeBranch; 38 39 class DelphesFactory; 40 class DelphesXDRReader; 39 41 40 42 class DelphesSTDHEPReader … … 76 78 TObjArray *partonOutputArray); 77 79 78 void SkipBytes( u_int size);79 void SkipArray( u_int elsize);80 void SkipBytes(int size); 81 void SkipArray(int elsize); 80 82 81 83 void ReadFileHeader(); … … 88 90 FILE *fInputFile; 89 91 90 XDR *fInputXDR;92 DelphesXDRReader fReader[7]; 91 93 92 char*fBuffer;94 uint8_t *fBuffer; 93 95 94 96 TDatabasePDG *fPDG; 95 97 96 u _int fEntries;97 int fBlockType, fEventNumber, fEventSize;98 uint32_t fEntries; 99 int32_t fBlockType, fEventNumber, fEventSize; 98 100 double fWeight, fAlphaQCD, fAlphaQED; 99 101 100 u _int fScaleSize;102 uint32_t fScaleSize; 101 103 double fScale[10]; 102 104 }; 103 105 104 106 #endif // DelphesSTDHEPReader_h 105 106 -
doc/genMakefile.tcl
r792092a r1716933 208 208 PcmSuf = _rdict.pcm 209 209 210 CXXFLAGS += $(ROOTCFLAGS) -Wno-write-strings -D_FILE_OFFSET_BITS=64 -DDROP_CGAL -I. -Iexternal -Iexternal/tcl -I/usr/include/tirpc210 CXXFLAGS += $(ROOTCFLAGS) -Wno-write-strings -D_FILE_OFFSET_BITS=64 -DDROP_CGAL -I. -Iexternal -Iexternal/tcl 211 211 DELPHES_LIBS = $(shell $(RC) --libs) -lEG $(SYSLIBS) 212 212 DISPLAY_LIBS = $(shell $(RC) --evelibs) -lGuiHtml $(SYSLIBS)
Note:
See TracChangeset
for help on using the changeset viewer.