KFPInputData.h 12 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371
  1. //----------------------------------------------------------------------------
  2. // Structures with input data for KF Particle Finder
  3. // .
  4. // @author M.Zyzak
  5. // @version 1.0
  6. // @since 20.08.13
  7. //
  8. //
  9. // -= Copyright &copy ALICE HLT and CBM L1 Groups =-
  10. //____________________________________________________________________________
  11. #ifndef KFPINPUTDATA_H
  12. #define KFPINPUTDATA_H
  13. #include "KFPTrackVector.h"
  14. #include "KFParticle.h"
  15. #include <vector>
  16. #include <string>
  17. #include <fstream>
  18. /** @class KFPTrackIndex
  19. ** @brief Helper structure to sort tracks in the KFPTrackVector object.
  20. ** @author M.Zyzak, I.Kisel
  21. ** @date 05.02.2019
  22. ** @version 1.0
  23. **
  24. ** The structure is used in the KFParticleTopoReconstructor::SortTracks() function.
  25. ** Tracks are sorted according to their pdg hypothesis: electrons, muons, pions,
  26. ** tracks without pdg (-1), kaons, protons, deuterons, tritons, He3, He4.
  27. ** Teh structure contains pdg hypothesis of the track and its index in the
  28. ** KFPTrackVector object.
  29. **/
  30. struct KFPTrackIndex
  31. {
  32. int fIndex; ///< index of the track in the KFPTrackVector object.
  33. int fPdg; ///< PDG hypothesis of the track
  34. static bool Compare(const KFPTrackIndex& a, const KFPTrackIndex& b)
  35. {
  36. /** Static sorting function for comparison of the two input objects of class KFPTrackIndex.
  37. ** Objects are sorted according to the PDG hypothesis: electrons, muons, pions,
  38. ** tracks without pdg (-1), kaons, protons, deuterons, tritons, He3, He4.
  39. ** Return "true" if a.fPdg < b.fPdg, otherwise returns "false".
  40. ** \param[in] a - first object
  41. ** \param[in] b - second object
  42. **/
  43. int pdg1 = a.fPdg == -1 ? 250 : a.fPdg;
  44. int pdg2 = b.fPdg == -1 ? 250 : b.fPdg;
  45. return (abs(pdg1) < abs(pdg2));
  46. }
  47. };
  48. /** @class KFPInputData
  49. ** @brief Class with the input data for KF Particle Finder: tracks, primary vertex and magnetic field.
  50. ** @author M.Zyzak, I.Kisel
  51. ** @date 05.02.2019
  52. ** @version 1.0
  53. **
  54. ** The class is used to transfer the data between devices: CPU and Intel Xeon Phi. The memory is aligned
  55. ** with the size of the SIMD vectors.
  56. **/
  57. class KFPInputData
  58. {
  59. public:
  60. void *operator new(size_t size) { return _mm_malloc(size, sizeof(float_v)); } ///< new operator for allocation of the SIMD-alligned dynamic memory allocation
  61. void *operator new[](size_t size) { return _mm_malloc(size, sizeof(float_v)); } ///< new operator for allocation of the SIMD-alligned dynamic memory allocation
  62. void *operator new(size_t size, void *ptr) { return ::operator new(size, ptr);} ///< new operator for allocation of the SIMD-alligned dynamic memory allocation
  63. void *operator new[](size_t size, void *ptr) { return ::operator new(size, ptr);} ///< new operator for allocation of the SIMD-alligned dynamic memory allocation
  64. void operator delete(void *ptr, size_t) { _mm_free(ptr); } ///< delete operator for the SIMD-alligned dynamic memory release
  65. void operator delete[](void *ptr, size_t) { _mm_free(ptr); } ///< delete operator for the SIMD-alligned dynamic memory release
  66. KFPInputData():fPV(0),fBz(0.f) {};
  67. ~KFPInputData() {};
  68. bool ReadDataFromFile( std::string prefix )
  69. {
  70. /** Reads the input data from the input file with the name defined by "prefix".
  71. ** \param[in] prefix - string with the name of the input file
  72. **/
  73. std::ifstream ifile(prefix.data());
  74. if ( !ifile.is_open() ) return 0;
  75. int nSets;
  76. ifile >> fBz;
  77. ifile >> nSets;
  78. for(int iSet=0; iSet<nSets; iSet++)
  79. {
  80. int nTracks = 0;
  81. ifile >> nTracks;
  82. fTracks[iSet].Resize(nTracks);
  83. for(int iP=0; iP<6; iP++)
  84. {
  85. float value;
  86. for(int iTr=0; iTr<fTracks[iSet].Size(); iTr++)
  87. {
  88. ifile >> value;
  89. fTracks[iSet].SetParameter(value, iP, iTr);
  90. }
  91. }
  92. for(int iC=0; iC<21; iC++)
  93. {
  94. float value;
  95. for(int iTr=0; iTr<fTracks[iSet].Size(); iTr++)
  96. {
  97. ifile >> value;
  98. fTracks[iSet].SetCovariance(value, iC, iTr);
  99. }
  100. }
  101. int tmpInt;
  102. for(int iTr=0; iTr<fTracks[iSet].Size(); iTr++)
  103. {
  104. ifile >> tmpInt;
  105. fTracks[iSet].SetId(tmpInt, iTr);
  106. }
  107. for(int iTr=0; iTr<fTracks[iSet].Size(); iTr++)
  108. {
  109. ifile >> tmpInt;
  110. fTracks[iSet].SetPDG(tmpInt, iTr);
  111. }
  112. for(int iTr=0; iTr<fTracks[iSet].Size(); iTr++)
  113. {
  114. ifile >> tmpInt;
  115. fTracks[iSet].SetQ(tmpInt, iTr);
  116. }
  117. for(int iTr=0; iTr<fTracks[iSet].Size(); iTr++)
  118. {
  119. ifile >> tmpInt;
  120. fTracks[iSet].SetPVIndex(tmpInt, iTr);
  121. }
  122. ifile >> tmpInt;
  123. fTracks[iSet].SetLastElectron(tmpInt);
  124. ifile >> tmpInt;
  125. fTracks[iSet].SetLastMuon (tmpInt);
  126. ifile >> tmpInt;
  127. fTracks[iSet].SetLastPion (tmpInt);
  128. ifile >> tmpInt;
  129. fTracks[iSet].SetLastKaon (tmpInt);
  130. ifile >> tmpInt;
  131. fTracks[iSet].SetLastProton (tmpInt);
  132. }
  133. int nPV;
  134. ifile>>nPV;
  135. fPV.resize(nPV);
  136. for(unsigned int iPV=0; iPV < fPV.size(); iPV++)
  137. {
  138. for(int iP=0; iP<3; iP++)
  139. ifile >> fPV[iPV].Parameter(iP);
  140. for(int iC=0; iC<6; iC++)
  141. ifile >> fPV[iPV].Covariance(iC);
  142. }
  143. ifile.close();
  144. return 1;
  145. }
  146. void SetDataToVector(int* data, int& dataSize)
  147. {
  148. /** Stores information to the memory under pointer "data".
  149. ** \param[out] data - memory, where input information will be stored
  150. ** \param[out] dataSize - size of the stored memory in "int" (or bloks of 4 bytes, or 32 bits)
  151. **/
  152. dataSize = NInputSets + 1 + 1; //sizes of the track vectors and pv vector, and field
  153. for(int iSet=0; iSet<NInputSets; iSet++)
  154. dataSize += fTracks[iSet].DataSize();
  155. dataSize += fPV.size() * 9;
  156. for(int iSet=0; iSet<NInputSets; iSet++)
  157. data[iSet] = fTracks[iSet].Size();
  158. data[NInputSets] = fPV.size();
  159. float& field = reinterpret_cast<float&>(data[NInputSets+1]);
  160. field = fBz;
  161. int offset = NInputSets+2;
  162. for(int iSet=0; iSet<NInputSets; iSet++)
  163. fTracks[iSet].SetDataToVector(data, offset);
  164. for(int iP=0; iP<3; iP++)
  165. {
  166. for(unsigned int iPV=0; iPV<fPV.size(); iPV++)
  167. {
  168. float& tmpFloat = reinterpret_cast<float&>(data[offset + iPV]);
  169. tmpFloat = fPV[iPV].Parameter(iP);
  170. }
  171. offset += fPV.size();
  172. }
  173. for(int iC=0; iC<6; iC++)
  174. {
  175. for(unsigned int iPV=0; iPV<fPV.size(); iPV++)
  176. {
  177. float& tmpFloat = reinterpret_cast<float&>(data[offset + iPV]);
  178. tmpFloat = fPV[iPV].Covariance(iC);
  179. }
  180. offset += fPV.size();
  181. }
  182. }
  183. void ReadDataFromVector(int* data)
  184. {
  185. /** Reads input data from the given memory.
  186. ** \param[in] data - pointer to the memory with the input data
  187. **/
  188. int offset = NInputSets+2;
  189. for(int iSet=0; iSet<NInputSets; iSet++)
  190. {
  191. fTracks[iSet].Resize(data[iSet]);
  192. fTracks[iSet].ReadDataFromVector(data, offset);
  193. }
  194. float& field = reinterpret_cast<float&>(data[NInputSets+1]);
  195. fBz = field;
  196. fPV.resize(data[NInputSets]);
  197. for(int iP=0; iP<3; iP++)
  198. {
  199. for(unsigned int iPV=0; iPV<fPV.size(); iPV++)
  200. {
  201. float& tmpFloat = reinterpret_cast<float&>(data[offset + iPV]);
  202. fPV[iPV].Parameter(iP) = tmpFloat;
  203. }
  204. offset += fPV.size();
  205. }
  206. for(int iC=0; iC<6; iC++)
  207. {
  208. for(unsigned int iPV=0; iPV<fPV.size(); iPV++)
  209. {
  210. float& tmpFloat = reinterpret_cast<float&>(data[offset + iPV]);
  211. fPV[iPV].Covariance(iC) = tmpFloat;
  212. }
  213. offset += fPV.size();
  214. }
  215. }
  216. void Print()
  217. {
  218. /**Prints all fields of the current object.*/
  219. for(int iSet=0; iSet<NInputSets; iSet++)
  220. fTracks[iSet].Print();
  221. std::cout << "N PV: " << fPV.size() << std::endl;
  222. std::cout << "X: ";
  223. for(unsigned int iPV=0; iPV<fPV.size(); iPV++)
  224. std::cout << fPV[iPV].X() <<" ";
  225. std::cout << std::endl;
  226. std::cout << "Y: ";
  227. for(unsigned int iPV=0; iPV<fPV.size(); iPV++)
  228. std::cout << fPV[iPV].Y() <<" ";
  229. std::cout << std::endl;
  230. std::cout << "Z: ";
  231. for(unsigned int iPV=0; iPV<fPV.size(); iPV++)
  232. std::cout << fPV[iPV].Z() <<" ";
  233. std::cout << std::endl;
  234. std::cout << "Cov matrix: " << std::endl;
  235. for(int iC=0; iC<6; iC++)
  236. {
  237. std::cout << " iC " << iC << ": ";
  238. for(unsigned int iPV=0; iPV<fPV.size(); iPV++)
  239. std::cout << fPV[iPV].Covariance(iC) <<" ";
  240. std::cout << std::endl;
  241. }
  242. std::cout << "Field: " << fBz << std::endl;
  243. }
  244. KFPTrackVector* GetTracks() { return fTracks; } ///< Returns pointer to the array with track vectors.
  245. float GetBz() const { return fBz; } ///< Returns value of the constant field Bz.
  246. const std::vector<KFParticle>& GetPV() const { return fPV; } ///< Returns vector with primary vertices.
  247. const KFPInputData& operator = (const KFPInputData& data)
  248. {
  249. /** Copies input data from object "data" to the current object. Returns the current object. \param[in] data - input data*/
  250. for(int i=0; i<NInputSets; i++)
  251. fTracks[i] = data.fTracks[i];
  252. fPV = data.fPV;
  253. fBz = data.fBz;
  254. return *this;
  255. }
  256. KFPInputData(const KFPInputData& data):fPV(0),fBz(0.f)
  257. {
  258. /** Copies input data from object "data" to the current object. \param[in] data - input data */
  259. for(int i=0; i<NInputSets; i++)
  260. fTracks[i] = data.fTracks[i];
  261. fPV = data.fPV;
  262. fBz = data.fBz;
  263. }
  264. protected:
  265. /** Array of track vectors: \n
  266. ** 0 - positive secondary tracks stored at the first point; \n
  267. ** 1 - negative secondary tracks stored at the first point; \n
  268. ** 2 - positive primary tracks stored at the first point; \n
  269. ** 3 - positive primary tracks stored at the first point; \n
  270. ** 4 - positive secondary tracks stored at the last point; \n
  271. ** 5 - negative secondary tracks stored at the last point; \n
  272. ** 6 - positive primary tracks stored at the last point; \n
  273. ** 7 - positive primary tracks stored at the last point.
  274. ** \see KFPTrackVector for documentation.
  275. **/
  276. KFPTrackVector fTracks[NInputSets]__attribute__((aligned(sizeof(float_v))));
  277. std::vector<KFParticle> fPV; ///< Vector with primary vertices.
  278. float fBz; ///< Constant homogenious one-component magnetic field Bz.
  279. } __attribute__((aligned(sizeof(float_v))));
  280. /** @class KFPInputDataArray
  281. ** @brief Structure with the set of the input data for KF Particle Finder.
  282. ** @author M.Zyzak, I.Kisel
  283. ** @date 05.02.2019
  284. ** @version 1.0
  285. **
  286. ** The structure contains pointer to array of KFPInputData objects. Copying of the
  287. ** objects of this structure is disabled.
  288. **/
  289. struct KFPInputDataArray{
  290. KFPInputDataArray():fInput(0){};
  291. ~KFPInputDataArray() { if(fInput) delete [] fInput; }
  292. KFPInputData *fInput; ///< Pointer to the array of the input data objects.
  293. private:
  294. const KFPInputDataArray& operator = (const KFPInputDataArray&);
  295. KFPInputDataArray(const KFPInputDataArray&);
  296. };
  297. /** @class KFPLinkedList
  298. ** @brief Structure to creat a linked list of the input data.
  299. ** @author M.Zyzak, I.Kisel
  300. ** @date 05.02.2019
  301. ** @version 1.0
  302. **
  303. ** The structure contains pointer to array of KFPInputData objects. Copying of the
  304. ** objects of this structure is disabled. The list is used to create a queue for processing
  305. ** at the device side (Intel Xeon Phi).
  306. **/
  307. struct KFPLinkedList
  308. {
  309. void *operator new(size_t size) { return _mm_malloc(size, sizeof(float_v)); } ///< new operator for allocation of the SIMD-alligned dynamic memory allocation
  310. void *operator new[](size_t size) { return _mm_malloc(size, sizeof(float_v)); } ///< new operator for allocation of the SIMD-alligned dynamic memory allocation
  311. void *operator new(size_t size, void *ptr) { return ::operator new(size, ptr);} ///< new operator for allocation of the SIMD-alligned dynamic memory allocation
  312. void *operator new[](size_t size, void *ptr) { return ::operator new(size, ptr);} ///< new operator for allocation of the SIMD-alligned dynamic memory allocation
  313. void operator delete(void *ptr, size_t) { _mm_free(ptr); } ///< delete operator for the SIMD-alligned dynamic memory release
  314. void operator delete[](void *ptr, size_t) { _mm_free(ptr); } ///< delete operator for the SIMD-alligned dynamic memory release
  315. KFPInputData data __attribute__((aligned(sizeof(float_v)))); ///< Input data for KF Particle Finder \see KFPInputData.
  316. KFPLinkedList* next; ///< Link to the nex object in the linked list.
  317. } __attribute__((aligned(sizeof(float_v))));
  318. #endif