musrfit  1.9.2
PRunDataHandler.cpp
Go to the documentation of this file.
1 /***************************************************************************
2 
3  PRunDataHandler.cpp
4 
5  Author: Andreas Suter
6  e-mail: andreas.suter@psi.ch
7 
8 ***************************************************************************/
9 
10 /***************************************************************************
11  * Copyright (C) 2007-2023 by Andreas Suter *
12  * andreas.suter@psi.ch *
13  * *
14  * This program is free software; you can redistribute it and/or modify *
15  * it under the terms of the GNU General Public License as published by *
16  * the Free Software Foundation; either version 2 of the License, or *
17  * (at your option) any later version. *
18  * *
19  * This program is distributed in the hope that it will be useful, *
20  * but WITHOUT ANY WARRANTY; without even the implied warranty of *
21  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
22  * GNU General Public License for more details. *
23  * *
24  * You should have received a copy of the GNU General Public License *
25  * along with this program; if not, write to the *
26  * Free Software Foundation, Inc., *
27  * 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. *
28  ***************************************************************************/
29 
30 #ifdef HAVE_CONFIG_H
31 #include "config.h"
32 #endif
33 
34 #include <sys/types.h>
35 #include <sys/stat.h>
36 #include <unistd.h>
37 #include <cstdio>
38 #include <ctime>
39 #include <iostream>
40 #include <iomanip>
41 #include <fstream>
42 #include <string>
43 #include <sstream>
44 #include <memory>
45 
46 #include <TROOT.h>
47 #include <TSystem.h>
48 #include <TString.h>
49 #include <TObjArray.h>
50 #include <TObjString.h>
51 #include <TFile.h>
52 #include <TFolder.h>
53 #include <TH1F.h>
54 #include <TDatime.h>
55 
56 #include "TMusrRunHeader.h"
57 #include "TLemRunHeader.h"
58 #include "MuSR_td_PSI_bin.h"
59 #include "mud.h"
60 
61 #ifdef PNEXUS_ENABLED
62 #include "PNeXus.h"
63 #endif
64 
65 #include "PRunDataHandler.h"
66 
67 #define PRH_MUSR_ROOT 0
68 #define PRH_LEM_ROOT 1
69 
70 #define PRH_NPP_OFFSET 0
71 #define PRH_PPC_OFFSET 20
72 
73 #define A2M_UNDEFINED 0
74 #define A2M_ROOT 1
75 #define A2M_MUSR_ROOT 2
76 #define A2M_PSIBIN 3
77 #define A2M_PSIMDU 4
78 #define A2M_MUD 5
79 #define A2M_NEXUS 6
80 #define A2M_WKM 7
81 #define A2M_ASCII 8
82 
83 #define PHR_INIT_ALL 0
84 #define PHR_INIT_MSR 1
85 #define PHR_INIT_ANY2MANY 2
86 
87 //--------------------------------------------------------------------------
88 // Constructor
89 //--------------------------------------------------------------------------
94 {
95  Init();
96 }
97 
98 //--------------------------------------------------------------------------
99 // Constructor
100 //--------------------------------------------------------------------------
107 PRunDataHandler::PRunDataHandler(TString fileName, const TString fileFormat) : fFileFormat(fileFormat)
108 {
109  Init();
110 
111  FileExistsCheck(fileName);
112 }
113 
114 //--------------------------------------------------------------------------
115 // Constructor
116 //--------------------------------------------------------------------------
124 PRunDataHandler::PRunDataHandler(TString fileName, const TString fileFormat, const PStringVector dataPath) : fDataPath(dataPath), fFileFormat(fileFormat)
125 {
126  Init();
127 
128  FileExistsCheck(fileName);
129 }
130 
131 //--------------------------------------------------------------------------
132 // Constructor
133 //--------------------------------------------------------------------------
142 PRunDataHandler::PRunDataHandler(TString fileName, const TString fileFormat, const TString dataPath, PRawRunData &runData)
143 {
144  Init();
145 }
146 
147 //--------------------------------------------------------------------------
148 // Constructor
149 //--------------------------------------------------------------------------
155 PRunDataHandler::PRunDataHandler(PAny2ManyInfo *any2ManyInfo) : fAny2ManyInfo(any2ManyInfo)
156 {
158 }
159 
160 //--------------------------------------------------------------------------
161 // Constructor
162 //--------------------------------------------------------------------------
170  fAny2ManyInfo(any2ManyInfo), fDataPath(dataPath)
171 {
173 }
174 
175 //--------------------------------------------------------------------------
176 // Constructor
177 //--------------------------------------------------------------------------
183 PRunDataHandler::PRunDataHandler(PMsrHandler *msrInfo) : fMsrInfo(msrInfo)
184 {
186 }
187 
188 //--------------------------------------------------------------------------
189 // Constructor
190 //--------------------------------------------------------------------------
199  fMsrInfo(msrInfo), fDataPath(dataPath)
200 {
202 }
203 
204 //--------------------------------------------------------------------------
205 // Destructor
206 //--------------------------------------------------------------------------
211 {
212  fDataPath.clear();
213  fData.clear();
214 }
215 
216 //--------------------------------------------------------------------------
217 // GetRunData
218 //--------------------------------------------------------------------------
229 PRawRunData* PRunDataHandler::GetRunData(const TString &runName)
230 {
231  UInt_t i;
232 
233  for (i=0; i<fData.size(); i++) {
234  if (!fData[i].GetRunName()->CompareTo(runName)) // run found
235  break;
236  }
237 
238  if (i == fData.size())
239  return nullptr;
240  else
241  return &fData[i];
242 }
243 
244 //--------------------------------------------------------------------------
245 // GetRunData
246 //--------------------------------------------------------------------------
257 {
258  if (idx >= fData.size())
259  return nullptr;
260  else
261  return &fData[idx];
262 }
263 
264 //--------------------------------------------------------------------------
265 // ReadData
266 //--------------------------------------------------------------------------
272 {
273  if (fMsrInfo) { // i.e. msr-file triggered
274  if (!ReadFilesMsr()) // couldn't read file
275  fAllDataAvailable = false;
276  else
277  fAllDataAvailable = true;
278  } else if (!fRunPathName.IsWhitespace()) { // i.e. file name triggered
279  if (fFileFormat == "") { // try to guess the file format if fromat hasn't been provided
280  TString ext=fRunPathName;
281  Ssiz_t pos=ext.Last('.');
282  ext.Remove(0, pos+1);
283  if (!ext.CompareTo("bin", TString::kIgnoreCase))
284  fFileFormat = "psibin";
285  else if (!ext.CompareTo("root", TString::kIgnoreCase))
286  fFileFormat = "musrroot";
287  else if (!ext.CompareTo("msr", TString::kIgnoreCase))
288  fFileFormat = "mud";
289  else if (!ext.CompareTo("nxs", TString::kIgnoreCase))
290  fFileFormat = "nexus";
291  else if (!ext.CompareTo("mdu", TString::kIgnoreCase))
292  fFileFormat = "psimdu";
293  }
294  if ((fFileFormat == "MusrRoot") || (fFileFormat == "musrroot")) {
296  } else if ((fFileFormat == "NeXus") || (fFileFormat == "nexus")) {
298  } else if ((fFileFormat == "PsiBin") || (fFileFormat == "psibin") ||
299  (fFileFormat == "PsiMdu") || (fFileFormat == "psimdu")) {
301  } else if ((fFileFormat == "Mud") || (fFileFormat == "mud")) {
303  } else if ((fFileFormat == "Wkm") || (fFileFormat == "wkm")) {
305  } else {
306  std::cerr << std::endl << ">> PRunDataHandler::ReadData(): **ERROR** unkown file format \"" << fFileFormat << "\" found." << std::endl;
307  fAllDataAvailable = false;
308  }
309  } else {
310  std::cerr << std::endl << ">> PRunDataHandler::ReadData(): **ERROR** Couldn't read files." << std::endl;
311  fAllDataAvailable = false;
312  }
313 }
314 
315 //--------------------------------------------------------------------------
316 // ConvertData
317 //--------------------------------------------------------------------------
322 {
323  if (!ReadWriteFilesList()) // couldn't read file
324  fAllDataAvailable = false;
325  else
326  fAllDataAvailable = true;
327 }
328 
329 //--------------------------------------------------------------------------
330 // SetRunData
331 //--------------------------------------------------------------------------
340 Bool_t PRunDataHandler::SetRunData(PRawRunData *data, UInt_t idx)
341 {
342  if ((idx == 0) && (fData.size() == 0)) {
343  fData.resize(1);
344  }
345  if (idx >= fData.size()) {
346  std::cerr << std::endl << ">>PRunDataHandler::SetRunData(): **ERROR** idx=" << idx << " is out-of-range (0.." << fData.size() << ")." << std::endl;
347  return false;
348  }
349 
350  fData[idx] = *data;
351 
352  return true;
353 }
354 
355 //--------------------------------------------------------------------------
356 // WriteData
357 //--------------------------------------------------------------------------
361 Bool_t PRunDataHandler::WriteData(TString fileName)
362 {
363  if ((fAny2ManyInfo == nullptr) && (fileName="")) {
364  std::cerr << std::endl << ">> PRunDataHandler::WriteData(): **ERROR** insufficient information: no fileName nor fAny2ManyInfo object.";
365  std::cerr << std::endl << " Cannot write data under this conditions." << std::endl;
366  return false;
367  }
368 
369  // get output file format tag, first try via fAny2ManyInfo
370  Int_t outTag = A2M_UNDEFINED;
371  if (fAny2ManyInfo != nullptr) {
372  if (!fAny2ManyInfo->outFormat.CompareTo("musrroot", TString::kIgnoreCase))
373  outTag = A2M_MUSR_ROOT;
374  else if (!fAny2ManyInfo->outFormat.CompareTo("psibin", TString::kIgnoreCase))
375  outTag = A2M_PSIBIN;
376  else if (!fAny2ManyInfo->outFormat.CompareTo("psimdu", TString::kIgnoreCase))
377  outTag = A2M_PSIMDU;
378  else if (!fAny2ManyInfo->outFormat.CompareTo("mud",TString::kIgnoreCase))
379  outTag = A2M_MUD;
380  else if (fAny2ManyInfo->outFormat.BeginsWith("nexus", TString::kIgnoreCase))
381  outTag = A2M_NEXUS;
382  else
383  outTag = A2M_UNDEFINED;
384  } else { // only fileName is given, try to guess from the extension
385  // STILL MISSING
386  }
387 
388  if (outTag == A2M_UNDEFINED) {
389  std::cerr << std::endl << ">> PRunDataHandler::WriteData(): **ERROR** no valid output data file format found: '" << fAny2ManyInfo->outFormat.Data() << "'" << std::endl;
390  return false;
391  }
392 
393  Bool_t success{true};
394  switch (outTag) {
395  case A2M_MUSR_ROOT:
396  if (fAny2ManyInfo->outFileName.Length() == 0)
397  success = WriteMusrRootFile(fileName);
398  else
400  break;
401  case A2M_PSIBIN:
402  case A2M_PSIMDU:
403  if (fAny2ManyInfo->outFileName.Length() == 0)
404  success = WritePsiBinFile(fileName);
405  else
407  break;
408  case A2M_MUD:
409  if (fAny2ManyInfo->outFileName.Length() == 0)
410  success = WriteMudFile(fileName);
411  else
413  break;
414  case A2M_NEXUS:
415  if (fAny2ManyInfo->outFileName.Length() == 0)
416  success = WriteNexusFile(fileName);
417  else
419  break;
420  default:
421  break;
422  }
423 
424  return success;
425 }
426 
427 //--------------------------------------------------------------------------
428 // Init (private)
429 //--------------------------------------------------------------------------
433 void PRunDataHandler::Init(const Int_t tag)
434 {
435  if ((tag==PHR_INIT_ALL) || (tag==PHR_INIT_ANY2MANY))
436  fMsrInfo = nullptr;
437  if ((tag==PHR_INIT_ALL) || (tag==PHR_INIT_MSR))
438  fAny2ManyInfo = nullptr;
439  fAllDataAvailable = false;
440  if (tag!=PHR_INIT_ALL)
441  fFileFormat = TString("");
442  fRunName = TString("");
443  fRunPathName = TString("");
444 }
445 
446 //--------------------------------------------------------------------------
447 // ReadFilesMsr (private)
448 //--------------------------------------------------------------------------
458 {
459  Bool_t success = true;
460 
461  // loop over the full RUN list to see what needs to be read
462  PMsrRunList *runList = nullptr;
463  runList = fMsrInfo->GetMsrRunList();
464  if (runList == nullptr) {
465  std::cerr << std::endl << ">> PRunDataHandler::ReadFilesMsr(): **ERROR** Couldn't obtain run list from PMsrHandler: something VERY fishy";
466  std::cerr << std::endl;
467  return false;
468  }
469 
470  char str[1024], *p_str=nullptr;
471  UInt_t year=0;
472  TString musrRoot("musr-root");
473 
474  for (UInt_t i=0; i<runList->size(); i++) {
475  for (UInt_t j=0; j<runList->at(i).GetRunNameSize(); j++) {
476  fRunName = *(runList->at(i).GetRunName(j));
477  // check if file actually exists
478  if (!FileExistsCheck(runList->at(i), j))
479  return false;
480 
481  // get year from string if LEM data file
482  strcpy(str, fRunName.Data());
483  p_str = strstr(str, "lem");
484  if (p_str != nullptr)
485  sscanf(p_str, "lem%d_his", &year);
486 
487  // check special case for ROOT-NPP/ROOT-PPC (LEM)
488  if (!runList->at(i).GetFileFormat(j)->CompareTo("root-npp")) { // not post pile up corrected histos
489  // if LEM file header is already TMusrRoot, change the data-file-format
490  if (year >= 12)
491  runList->at(i).SetFileFormat(musrRoot, j);
492  // check if forward/backward histoNo are within proper bounds, i.e. < PRH_PPC_OFFSET
493  for (UInt_t k=0; k<runList->at(i).GetForwardHistoNoSize(); k++) {
494  if (runList->at(i).GetForwardHistoNo(k) > PRH_PPC_OFFSET)
495  runList->at(i).SetForwardHistoNo(runList->at(i).GetForwardHistoNo(k)-PRH_PPC_OFFSET, k);
496  }
497  for (UInt_t k=0; k<runList->at(i).GetBackwardHistoNoSize(); k++) {
498  if (runList->at(i).GetBackwardHistoNo(k) > PRH_PPC_OFFSET)
499  runList->at(i).SetBackwardHistoNo(runList->at(i).GetBackwardHistoNo(k)-PRH_PPC_OFFSET, k);
500  }
501  } else if (!runList->at(i).GetFileFormat(j)->CompareTo("root-ppc")) { // post pile up corrected histos
502  // if LEM file header is already TMusrRoot, change the data-file-format
503  if (year >= 12)
504  runList->at(i).SetFileFormat(musrRoot, j);
505  // check if forward/backward histoNo are within proper bounds, i.e. > PRH_PPC_OFFSET
506  for (UInt_t k=0; k<runList->at(i).GetForwardHistoNoSize(); k++) {
507  if (runList->at(i).GetForwardHistoNo(k) < PRH_PPC_OFFSET)
508  runList->at(i).SetForwardHistoNo(runList->at(i).GetForwardHistoNo(k)+PRH_PPC_OFFSET, k);
509  }
510  for (UInt_t k=0; k<runList->at(i).GetBackwardHistoNoSize(); k++) {
511  if (runList->at(i).GetBackwardHistoNo(k) < PRH_PPC_OFFSET)
512  runList->at(i).SetBackwardHistoNo(runList->at(i).GetBackwardHistoNo(k)+PRH_PPC_OFFSET, k);
513  }
514  }
515 
516  // check is file is already read
517  if (FileAlreadyRead(*(runList->at(i).GetRunName(j))))
518  continue;
519  // everything looks fine, hence try to read the data file
520  if (!runList->at(i).GetFileFormat(j)->CompareTo("root-npp")) { // not post pile up corrected histos
521  success = ReadRootFile();
522  } else if (!runList->at(i).GetFileFormat(j)->CompareTo("root-ppc")) { // post pile up corrected histos
523  success = ReadRootFile();
524  } else if (!runList->at(i).GetFileFormat(j)->CompareTo("musr-root")) { // MusrRoot style file
525  success = ReadRootFile();
526  } else if (!runList->at(i).GetFileFormat(j)->CompareTo("nexus")) {
527  success = ReadNexusFile();
528  } else if (!runList->at(i).GetFileFormat(j)->CompareTo("psi-bin")) {
529  success = ReadPsiBinFile();
530  } else if (!runList->at(i).GetFileFormat(j)->CompareTo("psi-mdu")) {
531  success = ReadPsiBinFile();
532  } else if (!runList->at(i).GetFileFormat(j)->CompareTo("mud")) {
533  success = ReadMudFile();
534  } else if (!runList->at(i).GetFileFormat(j)->CompareTo("wkm")) {
535  success = ReadWkmFile();
536  } else if (!runList->at(i).GetFileFormat(j)->CompareTo("mdu-ascii")) {
537  success = ReadMduAsciiFile();
538  } else if (!runList->at(i).GetFileFormat(j)->CompareTo("ascii")) {
539  success = ReadAsciiFile();
540  } else if (!runList->at(i).GetFileFormat(j)->CompareTo("db")) {
541  success = ReadDBFile();
542  } else if (!runList->at(i).GetFileFormat(j)->CompareTo("dat")) {
543  success = ReadDatFile();
544  } else {
545  success = false;
546  }
547  }
548  }
549 
550  return success;
551 }
552 
553 //--------------------------------------------------------------------------
554 // ReadWriteFilesList (private)
555 //--------------------------------------------------------------------------
567 {
568  if ((fAny2ManyInfo->inFileName.size() == 0) && (fAny2ManyInfo->runList.size() == 0)) {
569  std::cerr << std::endl << ">> PRunDataHandler::ReadWriteFilesList(): **ERROR** Couldn't obtain run list from fAny2ManyInfo: something VERY fishy";
570  std::cerr << std::endl;
571  return false;
572  }
573 
574  // get input file format tag
575  Int_t inTag = A2M_UNDEFINED;
576  if (!fAny2ManyInfo->inFormat.CompareTo("root", TString::kIgnoreCase))
577  inTag = A2M_ROOT;
578  else if (!fAny2ManyInfo->inFormat.CompareTo("musrroot", TString::kIgnoreCase))
579  inTag = A2M_MUSR_ROOT;
580  else if (!fAny2ManyInfo->inFormat.CompareTo("psi-bin", TString::kIgnoreCase))
581  inTag = A2M_PSIBIN;
582  else if (!fAny2ManyInfo->inFormat.CompareTo("psi-mdu", TString::kIgnoreCase))
583  inTag = A2M_PSIMDU;
584  else if (!fAny2ManyInfo->inFormat.CompareTo("mud", TString::kIgnoreCase))
585  inTag = A2M_MUD;
586  else if (!fAny2ManyInfo->inFormat.CompareTo("nexus", TString::kIgnoreCase))
587  inTag = A2M_NEXUS;
588  else if (!fAny2ManyInfo->inFormat.CompareTo("wkm", TString::kIgnoreCase))
589  inTag = A2M_WKM;
590  else
591  inTag = A2M_UNDEFINED;
592 
593  if (inTag == A2M_UNDEFINED) {
594  std::cerr << std::endl << ">> PRunDataHandler::ReadWriteFilesList(): **ERROR** no valid input data file format found: '" << fAny2ManyInfo->inFormat.Data() << "'" << std::endl;
595  return false;
596  }
597 
598  // get output file format tag
599  Int_t outTag = A2M_UNDEFINED;
600  if (!fAny2ManyInfo->outFormat.CompareTo("root", TString::kIgnoreCase))
601  outTag = A2M_ROOT;
602  else if (!fAny2ManyInfo->outFormat.CompareTo("musrroot", TString::kIgnoreCase))
603  outTag = A2M_MUSR_ROOT;
604  else if (!fAny2ManyInfo->outFormat.CompareTo("psi-bin", TString::kIgnoreCase))
605  outTag = A2M_PSIBIN;
606  else if (!fAny2ManyInfo->outFormat.CompareTo("mud",TString::kIgnoreCase))
607  outTag = A2M_MUD;
608  else if (fAny2ManyInfo->outFormat.BeginsWith("nexus", TString::kIgnoreCase))
609  outTag = A2M_NEXUS;
610  else if (!fAny2ManyInfo->outFormat.CompareTo("wkm", TString::kIgnoreCase))
611  outTag = A2M_WKM;
612  else if (!fAny2ManyInfo->outFormat.CompareTo("ascii", TString::kIgnoreCase))
613  outTag = A2M_ASCII;
614  else
615  outTag = A2M_UNDEFINED;
616 
617  if (outTag == A2M_UNDEFINED) {
618  std::cerr << std::endl << ">> PRunDataHandler::ReadWriteFilesList(): **ERROR** no valid output data file format found: '" << fAny2ManyInfo->outFormat.Data() << "'" << std::endl;
619  return false;
620  }
621 
622  if (fAny2ManyInfo->inFileName.size() != 0) { // file name list given
623 
624  // loop over all runs of the run list
625  for (UInt_t i=0; i<fAny2ManyInfo->inFileName.size(); i++) {
626  if (!FileExistsCheck(true, i)) {
627  std::cerr << std::endl << ">> PRunDataHandler::ReadWriteFilesList(): **ERROR** Couldn't find file " << fAny2ManyInfo->inFileName[i].Data() << std::endl;
628  return false;
629  }
630 
631  // read input file
632  Bool_t success = false;
633  switch (inTag) {
634  case A2M_ROOT:
635  case A2M_MUSR_ROOT:
636  success = ReadRootFile();
637  break;
638  case A2M_PSIBIN:
639  case A2M_PSIMDU:
640  success = ReadPsiBinFile();
641  break;
642  case A2M_NEXUS:
643  success = ReadNexusFile();
644  break;
645  case A2M_MUD:
646  success = ReadMudFile();
647  break;
648  case A2M_WKM:
649  success = ReadWkmFile();
650  break;
651  default:
652  break;
653  }
654 
655  if (!success) {
656  std::cerr << std::endl << ">> PRunDataHandler::ReadWriteFilesList(): **ERROR** Couldn't read file " << fAny2ManyInfo->inFileName[i].Data() << std::endl;
657  return false;
658  }
659 
660  // write 'converted' output data file
661  success = false;
662  switch (outTag) {
663  case A2M_ROOT:
664  if (fAny2ManyInfo->outFileName.Length() == 0)
665  success = WriteRootFile();
666  else
668  break;
669  case A2M_MUSR_ROOT:
670  if (fAny2ManyInfo->outFileName.Length() == 0)
671  success = WriteMusrRootFile();
672  else
674  break;
675  case A2M_PSIBIN:
676  if (fAny2ManyInfo->outFileName.Length() == 0)
677  success = WritePsiBinFile();
678  else
680  break;
681  case A2M_PSIMDU:
682  if (fAny2ManyInfo->outFileName.Length() == 0)
683  success = WritePsiBinFile();
684  else
686  break;
687  case A2M_NEXUS:
688  if (fAny2ManyInfo->outFileName.Length() == 0)
689  success = WriteNexusFile();
690  else
692  break;
693  case A2M_MUD:
694  if (fAny2ManyInfo->outFileName.Length() == 0)
695  success = WriteMudFile();
696  else
698  break;
699  case A2M_WKM:
700  if (fAny2ManyInfo->outFileName.Length() == 0)
701  success = WriteWkmFile();
702  else
704  break;
705  case A2M_ASCII:
706  if (fAny2ManyInfo->outFileName.Length() == 0)
707  success = WriteAsciiFile();
708  else
710  break;
711  default:
712  break;
713  }
714 
715  if (success == false) {
716  std::cerr << std::endl << ">> PRunDataHandler::ReadWriteFilesList(): **ERROR** Couldn't write converted output file." << std::endl;
717  return false;
718  }
719 
720  // throw away the current data set
721  fData.clear();
722  }
723  }
724 
725  if (fAny2ManyInfo->runList.size() != 0) { // run list given
726 
727  // loop over all runs of the run list
728  for (UInt_t i=0; i<fAny2ManyInfo->runList.size(); i++) {
729  if (!FileExistsCheck(false, i)) {
730  std::cerr << std::endl << ">> PRunDataHandler::ReadWriteFilesList(): **ERROR** Couldn't find run " << fAny2ManyInfo->runList[i] << std::endl;
731  return false;
732  }
733 
734  // read input file
735  Bool_t success = false;
736  switch (inTag) {
737  case A2M_ROOT:
738  case A2M_MUSR_ROOT:
739  success = ReadRootFile();
740  break;
741  case A2M_PSIBIN:
742  case A2M_PSIMDU:
743  success = ReadPsiBinFile();
744  break;
745  case A2M_NEXUS:
746  success = ReadNexusFile();
747  break;
748  case A2M_MUD:
749  success = ReadMudFile();
750  break;
751  case A2M_WKM:
752  success = ReadWkmFile();
753  break;
754  default:
755  break;
756  }
757 
758  if (!success) {
759  std::cerr << std::endl << ">> PRunDataHandler::ReadWriteFilesList(): **ERROR** Couldn't read file " << fRunPathName.Data() << std::endl;
760  return false;
761  }
762 
763  TString year("");
764  TDatime dt;
765  year += dt.GetYear();
766  if (fAny2ManyInfo->year.Length() > 0)
767  year = fAny2ManyInfo->year;
768  Bool_t ok;
769  TString fln = FileNameFromTemplate(fAny2ManyInfo->outTemplate, fAny2ManyInfo->runList[i], year, ok);
770  if (!ok) {
771  std::cerr << std::endl << ">> PRunDataHandler::ReadWriteFilesList(): **ERROR** Couldn't create necessary output file name." << std::endl;
772  return false;
773  }
774 
775  // write 'converted' output data file
776  success = false;
777  switch (outTag) {
778  case A2M_ROOT:
779  success = WriteRootFile(fln);
780  break;
781  case A2M_MUSR_ROOT:
782  success = WriteMusrRootFile(fln);
783  break;
784  case A2M_PSIBIN:
785  success = WritePsiBinFile(fln);
786  break;
787  case A2M_PSIMDU:
788  success = WritePsiBinFile(fln);
789  break;
790  case A2M_NEXUS:
791  success = WriteNexusFile(fln);
792  break;
793  case A2M_MUD:
794  success = WriteMudFile(fln);
795  break;
796  case A2M_WKM:
797  success = WriteWkmFile(fln);
798  break;
799  case A2M_ASCII:
800  success = WriteAsciiFile(fln);
801  break;
802  default:
803  break;
804  }
805 
806  if (success == false) {
807  std::cerr << std::endl << ">> PRunDataHandler::ReadWriteFilesList(): **ERROR** Couldn't write converted output file." << std::endl;
808  return false;
809  }
810 
811  // throw away the current data set
812  fData.clear();
813  }
814 
815  }
816 
817  // check if compression is wished
818  if (fAny2ManyInfo->compressionTag > 0) {
820 
821  // currently system call is used, which means this is only running under Linux and Mac OS X but not under Windows
822  char cmd[256];
823  if (fAny2ManyInfo->outPathFileName.size() == 1) {
824  if (fAny2ManyInfo->compressionTag == 1) // gzip
825  fln += TString(".tar.gz");
826  else // bzip2
827  fln += TString(".tar.bz2");
828  if (fAny2ManyInfo->compressionTag == 1) // gzip
829  snprintf(cmd, sizeof(cmd), "tar -zcf %s %s", fln.Data(), fAny2ManyInfo->outPathFileName[0].Data());
830  else // bzip2
831  snprintf(cmd, sizeof(cmd), "tar -jcf %s %s", fln.Data(), fAny2ManyInfo->outPathFileName[0].Data());
832  if (system(cmd) == -1) {
833  std::cerr << "**ERROR** cmd: " << cmd << " failed." << std::endl;
834  }
835  } else {
836  fln += TString(".tar");
837  for (UInt_t i=0; i<fAny2ManyInfo->outPathFileName.size(); i++) {
838  if (i==0) {
839  snprintf(cmd, sizeof(cmd), "tar -cf %s %s", fln.Data(), fAny2ManyInfo->outPathFileName[i].Data());
840  } else {
841  snprintf(cmd, sizeof(cmd), "tar -rf %s %s", fln.Data(), fAny2ManyInfo->outPathFileName[i].Data());
842  }
843  if (system(cmd) == -1) {
844  std::cerr << "**ERROR** cmd: " << cmd << " failed." << std::endl;
845  }
846  }
847  if (fAny2ManyInfo->compressionTag == 1) { // gzip
848  snprintf(cmd, sizeof(cmd), "gzip %s", fln.Data());
849  fln += ".gz";
850  } else {
851  snprintf(cmd, sizeof(cmd), "bzip2 -z %s", fln.Data());
852  fln += ".bz2";
853  }
854  if (system(cmd) == -1) {
855  std::cerr << "**ERROR** cmd: " << cmd << " failed." << std::endl;
856  }
857  }
858 
859  // check if the compressed file shall be streamed to the stdout
861  // stream file to stdout
862  std::ifstream is;
863  int length=1024;
864  char *buffer;
865 
866  is.open(fln.Data(), std::ios::binary);
867  if (!is.is_open()) {
868  std::cerr << std::endl << "PRunDataHandler::ReadWriteFilesList(): **ERROR** Couldn't open the file for streaming." << std::endl;
869  remove(fln.Data());
870  return false;
871  }
872 
873  // get length of file
874  is.seekg(0, std::ios::end);
875  length = is.tellg();
876  is.seekg(0, std::ios::beg);
877 
878  if (length == -1) {
879  std::cerr << std::endl << "PRunDataHandler::ReadWriteFilesList(): **ERROR** Couldn't determine the file size." << std::endl;
880  remove(fln.Data());
881  return false;
882  }
883 
884  // allocate memory
885  buffer = new char [length];
886 
887  // read data as a block
888  while (!is.eof()) {
889  is.read(buffer, length);
890  std::cout.write(buffer, length);
891  }
892 
893  is.close();
894 
895  delete [] buffer;
896 
897  // delete temporary root file
898  remove(fln.Data());
899  }
900 
901  // remove all the converted files
902  for (UInt_t i=0; i<fAny2ManyInfo->outPathFileName.size(); i++)
903  remove(fAny2ManyInfo->outPathFileName[i].Data());
904  }
905 
906  return true;
907 }
908 
909 //--------------------------------------------------------------------------
910 // FileAlreadyRead (private)
911 //--------------------------------------------------------------------------
922 Bool_t PRunDataHandler::FileAlreadyRead(TString runName)
923 {
924  for (UInt_t i=0; i<fData.size(); i++) {
925  if (!fData[i].GetRunName()->CompareTo(runName)) { // run alread read
926  return true;
927  }
928  }
929 
930  return false;
931 }
932 
933 //--------------------------------------------------------------------------
934 // TestFileName (private)
935 //--------------------------------------------------------------------------
943 void PRunDataHandler::TestFileName(TString &runName, const TString &ext)
944 {
945  TString tmpStr(runName), tmpExt(ext);
946 
947  // check first if the file exists with the default extension which is not given with the run name
948  tmpStr += TString(".") + ext;
949  if (gSystem->AccessPathName(tmpStr.Data()) != true) { // found
950  runName = tmpStr;
951  return;
952  }
953 
954  // test if the file exists with the given run name but an only upper-case extension which is not included in the file name
955  tmpStr = runName;
956  tmpExt.ToUpper();
957  tmpStr += TString(".") + tmpExt;
958  if (gSystem->AccessPathName(tmpStr.Data()) != true) { // found
959  runName = tmpStr;
960  return;
961  }
962 
963  // test if the file exists with the given run name but an only lower-case extension which is not included in the file name
964  tmpStr = runName;
965  tmpExt.ToLower();
966  tmpStr += TString(".") + tmpExt;
967  if (gSystem->AccessPathName(tmpStr.Data()) != true) { // found
968  runName = tmpStr;
969  return;
970  }
971 
972  // test if the file exists with only upper-case letters
973  tmpStr = runName + TString(".") + tmpExt;
974  tmpStr.ToUpper();
975  if (gSystem->AccessPathName(tmpStr.Data()) != true) { // found
976  runName = tmpStr;
977  return;
978  }
979 
980  // test if the file exists with only lower-case letters
981  tmpStr.ToLower();
982  if (gSystem->AccessPathName(tmpStr.Data()) != true) { // found
983  runName = tmpStr;
984  return;
985  }
986 
987  tmpStr = runName;
988 
989  // the extension is already part of the file name, therefore, do not append it
990  if (tmpStr.EndsWith(ext.Data(), TString::kIgnoreCase)) {
991  if (gSystem->AccessPathName(tmpStr.Data()) != true) { // found
992  runName = tmpStr;
993  return;
994  }
995 
996  // assume some extension is part of the given file name but the real data file ends with an lower-case extension
997  tmpExt.ToLower();
998  tmpStr = tmpStr.Replace(static_cast<Ssiz_t>(tmpStr.Length()-tmpExt.Length()), tmpExt.Length(), tmpExt);
999  if (gSystem->AccessPathName(tmpStr.Data()) != true) { // found
1000  runName = tmpStr;
1001  return;
1002  }
1003 
1004  // assume some extension is part of the given file name but the real data file ends with an upper-case extension
1005  tmpExt.ToUpper();
1006  tmpStr = runName;
1007  tmpStr = tmpStr.Replace(static_cast<Ssiz_t>(tmpStr.Length()-tmpExt.Length()), tmpExt.Length(), tmpExt);
1008  if (gSystem->AccessPathName(tmpStr.Data()) != true) { // found
1009  runName = tmpStr;
1010  return;
1011  }
1012 
1013  // test if the file exists with only lower-case letters and the extension already included in the file name
1014  tmpStr = runName;
1015  tmpStr.ToLower();
1016  if (gSystem->AccessPathName(tmpStr.Data()) != true) { // found
1017  runName = tmpStr;
1018  return;
1019  }
1020 
1021  // test if the file exists with only upper-case letters and the extension already included in the file name
1022  tmpStr = runName;
1023  tmpStr.ToUpper();
1024  if (gSystem->AccessPathName(tmpStr.Data()) != true) { // found
1025  runName = tmpStr;
1026  return;
1027  }
1028  }
1029 
1030  // if the file has not been found, set the run name to be empty
1031  runName = "";
1032 
1033  return;
1034 }
1035 
1036 //--------------------------------------------------------------------------
1037 // FileExistsCheck (private)
1038 //--------------------------------------------------------------------------
1049 Bool_t PRunDataHandler::FileExistsCheck(PMsrRunBlock &runInfo, const UInt_t idx)
1050 {
1051  Bool_t success = true;
1052 
1053  // local init
1054  TROOT root("PRunBase", "PRunBase", nullptr);
1055  TString pathName = "???";
1056  TString str, *pstr;
1057  TString ext;
1058 
1059  pstr = runInfo.GetBeamline(idx);
1060  if (pstr == nullptr) {
1061  std::cerr << std::endl << ">> PRunDataHandler::FileExistsCheck: **ERROR** Couldn't obtain beamline data." << std::endl;
1062  assert(0);
1063  }
1064  pstr->ToLower();
1065  runInfo.SetBeamline(*pstr, idx);
1066  pstr = runInfo.GetInstitute(idx);
1067  if (pstr == nullptr) {
1068  std::cerr << std::endl << ">> PRunDataHandler::FileExistsCheck: **ERROR** Couldn't obtain institute data." << std::endl;
1069  assert(0);
1070  }
1071  pstr->ToLower();
1072  runInfo.SetInstitute(*pstr, idx);
1073  pstr = runInfo.GetFileFormat(idx);
1074  if (pstr == nullptr) {
1075  std::cerr << std::endl << ">> PRunDataHandler::FileExistsCheck: **ERROR** Couldn't obtain file format data." << std::endl;
1076  assert(0);
1077  }
1078  pstr->ToLower();
1079  runInfo.SetFileFormat(*pstr, idx);
1080 
1081  // file extensions for the various formats
1082  if (!runInfo.GetFileFormat(idx)->CompareTo("root-npp")) // not post pile up corrected histos
1083  ext = TString("root");
1084  else if (!runInfo.GetFileFormat(idx)->CompareTo("root-ppc")) // post pile up corrected histos
1085  ext = TString("root");
1086  else if (!runInfo.GetFileFormat(idx)->CompareTo("musr-root")) // post pile up corrected histos
1087  ext = TString("root");
1088  else if (!runInfo.GetFileFormat(idx)->CompareTo("nexus"))
1089  ext = TString("NXS");
1090  else if (!runInfo.GetFileFormat(idx)->CompareTo("psi-bin"))
1091  ext = TString("bin");
1092  else if (!runInfo.GetFileFormat(idx)->CompareTo("psi-mdu"))
1093  ext = TString("mdu");
1094  else if (!runInfo.GetFileFormat(idx)->CompareTo("mud"))
1095  ext = TString("msr");
1096  else if (!runInfo.GetFileFormat(idx)->CompareTo("wkm")) {
1097  if (!runInfo.GetBeamline(idx)->CompareTo("mue4"))
1098  ext = TString("nemu");
1099  else
1100  ext = *runInfo.GetBeamline(idx);
1101  }
1102  else if (!runInfo.GetFileFormat(idx)->CompareTo("mdu-ascii"))
1103  ext = TString("mdua");
1104  else if (!runInfo.GetFileFormat(idx)->CompareTo("ascii"))
1105  ext = TString("dat");
1106  else if (!runInfo.GetFileFormat(idx)->CompareTo("db"))
1107  ext = TString("db");
1108  else if (!runInfo.GetFileFormat(idx)->CompareTo("dat"))
1109  ext = TString("dat");
1110  else
1111  success = false;
1112 
1113  // unkown file format found
1114  if (!success) {
1115  pstr = runInfo.GetFileFormat(idx);
1116  if (pstr == nullptr) {
1117  std::cerr << std::endl << ">> PRunDataHandler::FileExistsCheck: **ERROR** Couldn't obtain file format data." << std::endl;
1118  assert(0);
1119  }
1120  pstr->ToUpper();
1121  std::cerr << std::endl << ">> PRunDataHandler::FileExistsCheck: **ERROR** File Format '" << pstr->Data() << "' unsupported.";
1122  std::cerr << std::endl << ">> support file formats are:";
1123  std::cerr << std::endl << ">> ROOT-NPP -> root not post pileup corrected for lem";
1124  std::cerr << std::endl << ">> ROOT-PPC -> root post pileup corrected for lem";
1125  std::cerr << std::endl << ">> MUSR-ROOT -> MusrRoot file format";
1126  std::cerr << std::endl << ">> NEXUS -> nexus file format, HDF4, HDF5, or XML";
1127  std::cerr << std::endl << ">> PSI-BIN -> psi bin file format";
1128  std::cerr << std::endl << ">> PSI-MDU -> psi mdu file format (see also MDU-ASCII)";
1129  std::cerr << std::endl << ">> MUD -> triumf mud file format";
1130  std::cerr << std::endl << ">> WKM -> wkm ascii file format";
1131  std::cerr << std::endl << ">> MDU-ASCII -> psi mdu ascii file format";
1132  std::cerr << std::endl << ">> ASCII -> column like file format";
1133  std::cerr << std::endl << ">> DB -> triumf db file \"format\"";
1134  std::cerr << std::endl << ">> DAT -> csv like file \"format\" as exported by msr2data.";
1135  std::cerr << std::endl;
1136  return success;
1137  }
1138 
1139  // check if the file is in the local directory
1140  str = *runInfo.GetRunName(idx);
1141  TestFileName(str, ext);
1142  if (!str.IsNull()) {
1143  pathName = str;
1144  }
1145 
1146  // check if the file is found in the <msr-file-directory>
1147  if (pathName.CompareTo("???") == 0) { // not found in local directory search
1148  str = *fMsrInfo->GetMsrFileDirectoryPath() + *runInfo.GetRunName(idx);
1149  TestFileName(str, ext);
1150  if (!str.IsNull()) {
1151  pathName = str;
1152  }
1153  }
1154 
1155  // check if the file is found in the directory given in the startup file
1156  if (pathName.CompareTo("???") == 0) { // not found in local directory search and not found in the <msr-file-directory> search
1157  for (UInt_t i=0; i<fDataPath.size(); i++) {
1158  str = fDataPath[i] + TString("/") + *runInfo.GetRunName(idx);
1159  TestFileName(str, ext);
1160  if (!str.IsNull()) {
1161  pathName = str;
1162  break;
1163  }
1164  }
1165  }
1166 
1167  // check if the file is found in the directories given by MUSRFULLDATAPATH
1168  const Char_t *musrpath = gSystem->Getenv("MUSRFULLDATAPATH");
1169  if (pathName.CompareTo("???") == 0) { // not found in local directory and xml path
1170  str = TString(musrpath);
1171  // MUSRFULLDATAPATH has the structure: path_1:path_2:...:path_n
1172  TObjArray *tokens = str.Tokenize(":");
1173  TObjString *ostr;
1174  for (Int_t i=0; i<tokens->GetEntries(); i++) {
1175  ostr = dynamic_cast<TObjString*>(tokens->At(i));
1176  str = ostr->GetString() + TString("/") + *runInfo.GetRunName(idx);
1177  TestFileName(str, ext);
1178  if (!str.IsNull()) {
1179  pathName = str;
1180  break;
1181  }
1182  }
1183  // cleanup
1184  if (tokens) {
1185  delete tokens;
1186  tokens = nullptr;
1187  }
1188  }
1189 
1190  // check if the file is found in the generated default path
1191  if (pathName.CompareTo("???") == 0) { // not found in MUSRFULLDATAPATH search
1192  str = TString(musrpath);
1193  // MUSRFULLDATAPATH has the structure: path_1:path_2:...:path_n
1194  TObjArray *tokens = str.Tokenize(":");
1195  TObjString *ostr;
1196  pstr = runInfo.GetInstitute(idx);
1197  if (pstr == nullptr) {
1198  std::cerr << std::endl << ">> PRunDataHandler::FileExistsCheck: **ERROR** Couldn't obtain institute data." << std::endl;
1199  assert(0);
1200  }
1201  pstr->ToUpper();
1202  runInfo.SetInstitute(*pstr, idx);
1203  pstr = runInfo.GetBeamline(idx);
1204  if (pstr == nullptr) {
1205  std::cerr << std::endl << ">> PRunDataHandler::FileExistsCheck: **ERROR** Couldn't obtain beamline data." << std::endl;
1206  assert(0);
1207  }
1208  pstr->ToUpper();
1209  runInfo.SetBeamline(*pstr, idx);
1210  for (Int_t i=0; i<tokens->GetEntries(); i++) {
1211  ostr = dynamic_cast<TObjString*>(tokens->At(i));
1212  str = ostr->GetString() + TString("/DATA/") +
1213  *runInfo.GetInstitute(idx) + TString("/") +
1214  *runInfo.GetBeamline(idx) + TString("/") +
1215  *runInfo.GetRunName(idx);
1216  TestFileName(str, ext);
1217  if (!str.IsNull()) { // found
1218  pathName = str;
1219  break;
1220  }
1221  }
1222  // clean up
1223  if (tokens) {
1224  delete tokens;
1225  tokens = nullptr;
1226  }
1227  }
1228 
1229  // no proper path name found
1230  if (pathName.CompareTo("???") == 0) {
1231  std::cerr << std::endl << ">> PRunDataHandler::FileExistsCheck(): **ERROR** Couldn't find '" << runInfo.GetRunName(idx)->Data() << "' in any standard path.";
1232  std::cerr << std::endl << ">> standard search pathes are:";
1233  std::cerr << std::endl << ">> 1. the local directory";
1234  std::cerr << std::endl << ">> 2. the data directory given in the startup XML file";
1235  std::cerr << std::endl << ">> 3. the directories listed in MUSRFULLDATAPATH";
1236  std::cerr << std::endl << ">> 4. default path construct which is described in the manual";
1237  return false;
1238  }
1239 
1240  fRunPathName = pathName;
1241 
1242  return true;
1243 }
1244 
1245 //--------------------------------------------------------------------------
1246 // FileExistsCheck (private)
1247 //--------------------------------------------------------------------------
1258 Bool_t PRunDataHandler::FileExistsCheck(const Bool_t fileName, const Int_t idx)
1259 {
1260  TString pathName("???");
1261  TString str("");
1262  TString fln("");
1263 
1264  if (fileName) { // single input file name
1265  if (idx >= static_cast<Int_t>(fAny2ManyInfo->inFileName.size())) {
1266  std::cerr << std::endl << ">> PRunDataHandler::FileExistsCheck(): **ERROR** idx=" << idx << " out of range. (inFileName.size()==" << fAny2ManyInfo->inFileName.size() << ")" << std::endl;
1267  return false;
1268  }
1269  fln = fAny2ManyInfo->inFileName[idx];
1270  } else { // run file list entry shall be handled
1271  if (idx >= static_cast<Int_t>(fAny2ManyInfo->runList.size())) {
1272  std::cerr << std::endl << ">> PRunDataHandler::FileExistsCheck(): **ERROR** idx=" << idx << " out of range. (inFileName.size()==" << fAny2ManyInfo->runList.size() << ")" << std::endl;
1273  return false;
1274  }
1275  // check for input/output templates
1276  if ((fAny2ManyInfo->inTemplate.Length() == 0) || (fAny2ManyInfo->outTemplate.Length() == 0)) {
1277  std::cerr << std::endl << ">> PRunDataHandler::FileExistsCheck(): **ERROR** when using run lists, input/output templates are needed as well." << std::endl;
1278  return false;
1279  }
1280  // make the input file name according to the input template
1281  TString year("");
1282  TDatime dt;
1283  year += dt.GetYear();
1284  if (fAny2ManyInfo->year.Length() > 0)
1285  year = fAny2ManyInfo->year;
1286  Bool_t ok;
1288  if (!ok)
1289  return false;
1290  }
1291 
1292  // check if the file is in the local directory
1293  if (gSystem->AccessPathName(fln) != true) { // found in the local dir
1294  pathName = fln;
1295  }
1296  // check if the file is found in the directory given in the startup file
1297  if (pathName.CompareTo("???") == 0) { // not found in local directory search
1298  for (UInt_t i=0; i<fDataPath.size(); i++) {
1299  str = fDataPath[i] + TString("/") + fln;
1300  if (gSystem->AccessPathName(str.Data())!=true) { // found
1301  pathName = str;
1302  break;
1303  }
1304  }
1305  }
1306  // check if the file is found in the directories given by MUSRFULLDATAPATH
1307  const Char_t *musrpath = gSystem->Getenv("MUSRFULLDATAPATH");
1308  if (pathName.CompareTo("???") == 0) { // not found in local directory and xml path
1309  str = TString(musrpath);
1310  // MUSRFULLDATAPATH has the structure: path_1:path_2:...:path_n
1311  TObjArray *tokens = str.Tokenize(":");
1312  TObjString *ostr;
1313  for (Int_t i=0; i<tokens->GetEntries(); i++) {
1314  ostr = dynamic_cast<TObjString*>(tokens->At(i));
1315  str = ostr->GetString() + TString("/") + fln;
1316  if (gSystem->AccessPathName(str.Data())!=true) { // found
1317  pathName = str;
1318  break;
1319  }
1320  }
1321  // cleanup
1322  if (tokens) {
1323  delete tokens;
1324  tokens = nullptr;
1325  }
1326  }
1327  // no proper path name found
1328  if (pathName.CompareTo("???") == 0) {
1329  std::cerr << std::endl << ">> PRunDataHandler::FileExistsCheck(): **ERROR** Couldn't find '" << fln.Data() << "' in any standard path.";
1330  std::cerr << std::endl << ">> standard search pathes are:";
1331  std::cerr << std::endl << ">> 1. the local directory";
1332  std::cerr << std::endl << ">> 2. the data directory given in the startup XML file";
1333  std::cerr << std::endl << ">> 3. the directories listed in MUSRFULLDATAPATH";
1334  return false;
1335  }
1336 
1337  fRunPathName = pathName;
1338 
1339  return true;
1340 }
1341 
1342 //--------------------------------------------------------------------------
1343 // FileExistsCheck (private)
1344 //--------------------------------------------------------------------------
1354 Bool_t PRunDataHandler::FileExistsCheck(const TString fileName)
1355 {
1356  TString pathName("???");
1357  TString str("");
1358 
1359  // check if the file is in the local directory
1360  if (gSystem->AccessPathName(fileName) != true) { // found in the local dir
1361  pathName = fileName;
1362  }
1363  // check if the file is found in the directory given in the startup file
1364  if (pathName.CompareTo("???") == 0) { // not found in local directory search
1365  for (UInt_t i=0; i<fDataPath.size(); i++) {
1366  str = fDataPath[i] + TString("/") + fileName;
1367  if (gSystem->AccessPathName(str.Data())!=true) { // found
1368  pathName = str;
1369  break;
1370  }
1371  }
1372  }
1373  // check if the file is found in the directories given by MUSRFULLDATAPATH
1374  const Char_t *musrpath = gSystem->Getenv("MUSRFULLDATAPATH");
1375  if (pathName.CompareTo("???") == 0) { // not found in local directory and xml path
1376  str = TString(musrpath);
1377  // MUSRFULLDATAPATH has the structure: path_1:path_2:...:path_n
1378  TObjArray *tokens = str.Tokenize(":");
1379  TObjString *ostr;
1380  for (Int_t i=0; i<tokens->GetEntries(); i++) {
1381  ostr = dynamic_cast<TObjString*>(tokens->At(i));
1382  str = ostr->GetString() + TString("/") + fileName;
1383  if (gSystem->AccessPathName(str.Data())!=true) { // found
1384  pathName = str;
1385  break;
1386  }
1387  }
1388  // cleanup
1389  if (tokens) {
1390  delete tokens;
1391  tokens = nullptr;
1392  }
1393  }
1394  // no proper path name found
1395  if (pathName.CompareTo("???") == 0) {
1396  std::cerr << std::endl << ">> PRunDataHandler::FileExistsCheck(): **ERROR** Couldn't find '" << fileName.Data() << "' in any standard path.";
1397  std::cerr << std::endl << ">> standard search pathes are:";
1398  std::cerr << std::endl << ">> 1. the local directory";
1399  std::cerr << std::endl << ">> 2. the data directory given in the startup XML file";
1400  std::cerr << std::endl << ">> 3. the directories listed in MUSRFULLDATAPATH";
1401  return false;
1402  }
1403 
1404  fRunPathName = pathName;
1405 
1406  return true;
1407 }
1408 
1409 //--------------------------------------------------------------------------
1410 // ReadRootFile (private)
1411 //--------------------------------------------------------------------------
1421 {
1422  PDoubleVector histoData;
1423  PRawRunData runData;
1424  PRawRunDataSet dataSet;
1425 
1426  TFile f(fRunPathName.Data());
1427  if (f.IsZombie()) {
1428  return false;
1429  }
1430 
1431  UInt_t fileType = PRH_MUSR_ROOT;
1432 
1433  TFolder *folder;
1434  f.GetObject("RunInfo", folder); // try first LEM-ROOT style file (used until 2011).
1435  if (!folder) { // either something is wrong, or it is a MusrRoot file
1436  f.GetObject("RunHeader", folder);
1437  if (!folder) { // something is wrong!!
1438  std::cerr << std::endl << ">> PRunDataHandler::ReadRootFile: **ERROR** Couldn't neither obtain RunInfo (LEM),";
1439  std::cerr << std::endl << " nor RunHeader (MusrRoot) from " << fRunPathName.Data() << std::endl;
1440  f.Close();
1441  return false;
1442  } else {
1443  fileType = PRH_MUSR_ROOT;
1444  }
1445  } else {
1446  fileType = PRH_LEM_ROOT;
1447  }
1448 
1449  if (fileType == PRH_LEM_ROOT) {
1450  // read header and check if some missing run info need to be fed
1451  TLemRunHeader *runHeader = dynamic_cast<TLemRunHeader*>(folder->FindObjectAny("TLemRunHeader"));
1452 
1453  // check if run header is valid
1454  if (!runHeader) {
1455  std::cerr << std::endl << ">> PRunDataHandler::ReadRootFile: **ERROR** Couldn't obtain run header info from ROOT file " << fRunPathName.Data() << std::endl;
1456  f.Close();
1457  return false;
1458  }
1459 
1460  // set laboratory / beamline / instrument
1461  runData.SetLaboratory("PSI");
1462  runData.SetBeamline("muE4");
1463  runData.SetInstrument("LEM");
1464  runData.SetMuonSource("low energy muon source");
1465  runData.SetMuonSpecies("positive muons");
1466 
1467  // get run title
1468  TObjString ostr = runHeader->GetRunTitle();
1469  runData.SetRunTitle(ostr.GetString());
1470 
1471  // get run number
1472  runData.SetRunNumber(static_cast<Int_t>(runHeader->GetRunNumber()));
1473 
1474  // get temperature
1475  runData.ClearTemperature();
1476  runData.SetTemperature(0, runHeader->GetSampleTemperature(), runHeader->GetSampleTemperatureError());
1477 
1478  // get field
1479  runData.SetField(runHeader->GetSampleBField());
1480 
1481  // get implantation energy
1482  runData.SetEnergy(runHeader->GetImpEnergy());
1483 
1484  // get moderator HV
1485  runData.SetTransport(runHeader->GetModeratorHV());
1486 
1487  // get setup
1488  runData.SetSetup(runHeader->GetLemSetup().GetString());
1489 
1490  // get start time/date
1491  // start date
1492  time_t idt = static_cast<time_t>(runHeader->GetStartTime());
1493  runData.SetStartDateTime(idt);
1494  struct tm *dt = localtime(&idt);
1495  char str[128];
1496  strftime(str, sizeof(str), "%F", dt);
1497  TString stime(str);
1498  runData.SetStartDate(stime);
1499  // start time
1500  memset(str, 0, sizeof(str));
1501  strftime(str, sizeof(str), "%T", dt);
1502  stime = str;
1503  runData.SetStartTime(stime);
1504 
1505  // get stop time/date
1506  // stop date
1507  idt = static_cast<time_t>(runHeader->GetStopTime());
1508  runData.SetStopDateTime(idt);
1509  dt = localtime(&idt);
1510  memset(str, 0, sizeof(str));
1511  strftime(str, sizeof(str), "%F", dt);
1512  stime = str;
1513  runData.SetStopDate(stime);
1514  // stop time
1515  memset(str, 0, sizeof(str));
1516  strftime(str, sizeof(str), "%T", dt);
1517  stime = str;
1518  runData.SetStopTime(stime);
1519 
1520  // get time resolution
1521  runData.SetTimeResolution(runHeader->GetTimeResolution());
1522 
1523  // get number of histogramms
1524  Int_t noOfHistos = runHeader->GetNHist();
1525 
1526  // get t0's there will be handled together with the data
1527  Double_t *t0 = runHeader->GetTimeZero();
1528 
1529  // read run summary to obtain ring anode HV values
1530  TObjArray *runSummary = dynamic_cast<TObjArray*>(folder->FindObjectAny("RunSummary"));
1531 
1532  // check if run summary is valid
1533  if (!runSummary) {
1534  std::cout << std::endl << "**INFO** Couldn't obtain run summary info from ROOT file " << fRunPathName.Data() << std::endl;
1535  // this is not fatal... only RA-HV values are not available
1536  } else { // it follows a (at least) little bit strange extraction of the RA values from Thomas' TObjArray...
1537  //streaming of a ASCII-file would be more easy
1538  TString s;
1539  TObjArrayIter summIter(runSummary);
1540  TObjString *os(dynamic_cast<TObjString*>(summIter.Next()));
1541  TObjArray *oa(nullptr);
1542  TObjString *objTok(nullptr);
1543  while (os != nullptr) {
1544  s = os->GetString();
1545  // will put four parallel if's since it may be that more than one RA-values are on one line
1546  if (s.Contains("RA-L")) {
1547  oa = s.Tokenize(" ");
1548  TObjArrayIter lineIter(oa);
1549  objTok = dynamic_cast<TObjString*>(lineIter.Next());
1550  while (objTok != nullptr) {
1551  if (!objTok->GetString().CompareTo("RA-L")) {
1552  objTok = dynamic_cast<TObjString*>(lineIter.Next()); // "="
1553  if ((objTok != nullptr) && !objTok->GetString().CompareTo("=")) {
1554  objTok = dynamic_cast<TObjString*>(lineIter.Next()); // HV value
1555  runData.SetRingAnode(0, objTok->GetString().Atof()); // fill RA-R value into the runData structure
1556  break;
1557  }
1558  }
1559  objTok = dynamic_cast<TObjString*>(lineIter.Next()); // next token...
1560  }
1561  // clean up
1562  if (oa) {
1563  delete oa;
1564  oa = nullptr;
1565  }
1566  }
1567 
1568  if (s.Contains("RA-R")) {
1569  oa = s.Tokenize(" ");
1570  TObjArrayIter lineIter(oa);
1571  objTok = dynamic_cast<TObjString*>(lineIter.Next());
1572  while (objTok != nullptr){
1573  if (!objTok->GetString().CompareTo("RA-R")) {
1574  objTok = dynamic_cast<TObjString*>(lineIter.Next()); // "="
1575  if (objTok != nullptr && !objTok->GetString().CompareTo("=")) {
1576  objTok = dynamic_cast<TObjString*>(lineIter.Next()); // HV value
1577  runData.SetRingAnode(1, objTok->GetString().Atof()); // fill RA-R value into the runData structure
1578  break;
1579  }
1580  }
1581  objTok = dynamic_cast<TObjString*>(lineIter.Next()); // next token...
1582  }
1583  // clean up
1584  if (oa) {
1585  delete oa;
1586  oa = nullptr;
1587  }
1588  }
1589 
1590  if (s.Contains("RA-T")) {
1591  oa = s.Tokenize(" ");
1592  TObjArrayIter lineIter(oa);
1593  objTok = dynamic_cast<TObjString*>(lineIter.Next());
1594  while (objTok != nullptr){
1595  if (!objTok->GetString().CompareTo("RA-T")) {
1596  objTok = dynamic_cast<TObjString*>(lineIter.Next()); // "="
1597  if ((objTok != nullptr) && !objTok->GetString().CompareTo("=")) {
1598  objTok = dynamic_cast<TObjString*>(lineIter.Next()); // HV value
1599  runData.SetRingAnode(2, objTok->GetString().Atof()); // fill RA-T value into the runData structure
1600  break;
1601  }
1602  }
1603  objTok = dynamic_cast<TObjString*>(lineIter.Next()); // next token...
1604  }
1605  // clean up
1606  if (oa) {
1607  delete oa;
1608  oa = nullptr;
1609  }
1610  }
1611 
1612  if (s.Contains("RA-B")) {
1613  oa = s.Tokenize(" ");
1614  TObjArrayIter lineIter(oa);
1615  objTok = dynamic_cast<TObjString*>(lineIter.Next());
1616  while (objTok != nullptr){
1617  if (!objTok->GetString().CompareTo("RA-B")) {
1618  objTok = dynamic_cast<TObjString*>(lineIter.Next()); // "="
1619  if ((objTok != nullptr) && !objTok->GetString().CompareTo("=")) {
1620  objTok = dynamic_cast<TObjString*>(lineIter.Next()); // HV value
1621  runData.SetRingAnode(3, objTok->GetString().Atof()); // fill RA-B value into the runData structure
1622  break;
1623  }
1624  }
1625  objTok = dynamic_cast<TObjString*>(lineIter.Next()); // next token...
1626  }
1627  // clean up
1628  if (oa) {
1629  delete oa;
1630  oa = nullptr;
1631  }
1632  }
1633 
1634  os = dynamic_cast<TObjString*>(summIter.Next()); // next summary line...
1635  }
1636  }
1637 
1638  // read data ---------------------------------------------------------
1639 
1640  // check if histos folder is found
1641  f.GetObject("histos", folder);
1642  if (!folder) {
1643  std::cerr << std::endl << ">> PRunDataHandler::ReadRootFile: **ERROR** Couldn't obtain histos from " << fRunPathName.Data() << std::endl;
1644  f.Close();
1645  return false;
1646  }
1647 
1648  // get all the data
1649  Char_t histoName[32];
1650  for (Int_t i=0; i<noOfHistos; i++) {
1651  snprintf(histoName, sizeof(histoName), "hDecay%02d", i);
1652  TH1F *histo = dynamic_cast<TH1F*>(folder->FindObjectAny(histoName));
1653  if (!histo) {
1654  std::cerr << std::endl << ">> PRunDataHandler::ReadRootFile: **ERROR** Couldn't get histo " << histoName;
1655  std::cerr << std::endl;
1656  f.Close();
1657  return false;
1658  }
1659 
1660  // fill data
1661  for (Int_t j=1; j<=histo->GetNbinsX(); j++) {
1662  histoData.push_back(histo->GetBinContent(j));
1663  }
1664 
1665  // fill the data set
1666  dataSet.Clear();
1667  dataSet.SetName(histo->GetTitle());
1668  dataSet.SetHistoNo(i+1);
1669  dataSet.SetTimeZeroBin(t0[i]);
1670  dataSet.SetTimeZeroBinEstimated(histo->GetMaximumBin());
1671  dataSet.SetFirstGoodBin(static_cast<Int_t>(t0[i]));
1672  dataSet.SetLastGoodBin(histo->GetNbinsX()-1);
1673  dataSet.SetData(histoData);
1674  runData.SetDataSet(dataSet);
1675 
1676  // clear histoData for the next histo
1677  histoData.clear();
1678  }
1679  // check if any post pileup histos are present at all (this is not the case for LEM data 2006 and earlier)
1680  snprintf(histoName, sizeof(histoName), "hDecay%02d", POST_PILEUP_HISTO_OFFSET);
1681  if (!folder->FindObjectAny(histoName)) {
1682  std::cerr << std::endl << ">> PRunDataHandler::ReadRootFile: **WARNING** Couldn't get histo " << histoName;
1683  std::cerr << std::endl << ">> most probably this is an old (2006 or earlier) LEM file without post pileup histos.";
1684  std::cerr << std::endl;
1685  } else {
1686  for (Int_t i=0; i<noOfHistos; i++) {
1687  snprintf(histoName, sizeof(histoName), "hDecay%02d", i+POST_PILEUP_HISTO_OFFSET);
1688  TH1F *histo = dynamic_cast<TH1F*>(folder->FindObjectAny(histoName));
1689  if (!histo) {
1690  std::cerr << std::endl << ">> PRunDataHandler::ReadRootFile: **ERROR** Couldn't get histo " << histoName;
1691  std::cerr << std::endl;
1692  f.Close();
1693  return false;
1694  }
1695 
1696  // fill data
1697  for (Int_t j=1; j<=histo->GetNbinsX(); j++)
1698  histoData.push_back(histo->GetBinContent(j));
1699 
1700  // fill the data set
1701  dataSet.Clear();
1702  dataSet.SetName(histo->GetTitle());
1703  dataSet.SetHistoNo(i+1+POST_PILEUP_HISTO_OFFSET);
1704  dataSet.SetTimeZeroBin(t0[i]);
1705  dataSet.SetTimeZeroBinEstimated(histo->GetMaximumBin());
1706  dataSet.SetFirstGoodBin(static_cast<Int_t>(t0[i]));
1707  dataSet.SetLastGoodBin(histo->GetNbinsX()-1);
1708  dataSet.SetData(histoData);
1709  runData.SetDataSet(dataSet);
1710 
1711  // clear histoData for the next histo
1712  histoData.clear();
1713  }
1714  }
1715  } else { // MusrRoot file
1716  // invoke the MusrRoot header object
1717  std::unique_ptr<TMusrRunHeader> header = std::make_unique<TMusrRunHeader>(true); // read quiet
1718  if (header == nullptr) {
1719  std::cerr << std::endl << ">> PRunDataHandler::ReadRootFile: **ERROR** Couldn't invoke MusrRoot RunHeader in file:" << fRunPathName;
1720  std::cerr << std::endl;
1721  f.Close();
1722  return false;
1723  }
1724 
1725  // try to populate the MusrRoot header object
1726  if (!header->ExtractAll(folder)) {
1727  std::cerr << std::endl << ">> PRunDataHandler::ReadRootFile: **ERROR** Couldn't invoke MusrRoot RunHeader in file:" << fRunPathName;
1728  std::cerr << std::endl;
1729  f.Close();
1730  return false;
1731  }
1732 
1733  // get all the header information
1734  Bool_t ok;
1735  TString str, path, pathName;
1736  Int_t ival, noOfHistos=0;
1737  Double_t dval;
1739  PIntVector ivec, redGreenOffsets;
1740 
1741  header->Get("RunInfo/Version", str, ok);
1742  if (ok)
1743  runData.SetVersion(str);
1744 
1745  header->Get("RunInfo/Generic Validator URL", str, ok);
1746  if (ok)
1747  runData.SetGenericValidatorUrl(str);
1748 
1749  header->Get("RunInfo/Specific Validator URL", str, ok);
1750  if (ok)
1751  runData.SetSpecificValidatorUrl(str);
1752 
1753  header->Get("RunInfo/Generator", str, ok);
1754  if (ok)
1755  runData.SetGenerator(str);
1756 
1757  header->Get("RunInfo/File Name", str, ok);
1758  if (ok)
1759  runData.SetFileName(str);
1760 
1761  header->Get("RunInfo/Run Title", str, ok);
1762  if (ok)
1763  runData.SetRunTitle(str);
1764 
1765  header->Get("RunInfo/Run Number", ival, ok);
1766  if (ok)
1767  runData.SetRunNumber(ival);
1768 
1769  header->Get("RunInfo/Run Start Time", str, ok);
1770  if (ok) {
1771  Ssiz_t pos = str.Index(' ');
1772  TString substr = str;
1773  substr.Remove(pos, str.Length());
1774  runData.SetStartDate(substr);
1775  substr = str;
1776  substr.Remove(0, pos+1);
1777  runData.SetStartTime(substr);
1778  }
1779 
1780  header->Get("RunInfo/Run Stop Time", str, ok);
1781  if (ok) {
1782  Ssiz_t pos = str.Index(' ');
1783  TString substr = str;
1784  substr.Remove(pos, str.Length());
1785  runData.SetStopDate(substr);
1786  substr = str;
1787  substr.Remove(0, pos+1);
1788  runData.SetStopTime(substr);
1789  }
1790 
1791  header->Get("RunInfo/Laboratory", str, ok);
1792  if (ok)
1793  runData.SetLaboratory(str);
1794 
1795  header->Get("RunInfo/Instrument", str, ok);
1796  if (ok)
1797  runData.SetInstrument(str);
1798 
1799  header->Get("RunInfo/Muon Beam Momentum", prop, ok);
1800  if (ok) {
1801  if (!prop.GetUnit().CompareTo("MeV/c"))
1802  runData.SetMuonBeamMomentum(prop.GetValue());
1803  }
1804 
1805  header->Get("RunInfo/Muon Species", str, ok);
1806  if (ok)
1807  runData.SetMuonSpecies(str);
1808 
1809  header->Get("RunInfo/Muon Source", str, ok);
1810  if (ok)
1811  runData.SetMuonSource(str);
1812 
1813  header->Get("RunInfo/Muon Spin Angle", prop, ok);
1814  if (ok) {
1815  runData.SetMuonSpinAngle(prop.GetValue());
1816  }
1817 
1818  header->Get("RunInfo/Setup", str, ok);
1819  if (ok)
1820  runData.SetSetup(str);
1821 
1822  header->Get("RunInfo/Comment", str, ok);
1823  if (ok)
1824  runData.SetComment(str);
1825 
1826  header->Get("RunInfo/Sample Name", str, ok);
1827  if (ok)
1828  runData.SetSample(str);
1829 
1830  header->Get("RunInfo/Sample Temperature", prop, ok);
1831  if (ok)
1832  runData.SetTemperature(0, prop.GetValue(), prop.GetError());
1833 
1834  header->Get("RunInfo/Sample Magnetic Field", prop, ok);
1835  if (ok) {
1836  dval = MRH_UNDEFINED;
1837  if (!prop.GetUnit().CompareTo("G") || !prop.GetUnit().CompareTo("Gauss"))
1838  dval = prop.GetValue();
1839  else if (!prop.GetUnit().CompareTo("T") || !prop.GetUnit().CompareTo("Tesla"))
1840  dval = prop.GetValue() * 1.0e4;
1841  runData.SetField(dval);
1842  }
1843 
1844  header->Get("RunInfo/No of Histos", ival, ok);
1845  if (ok) {
1846  noOfHistos = ival;
1847  }
1848 
1849  header->Get("RunInfo/Time Resolution", prop, ok);
1850  if (ok) {
1851  dval = -1.0;
1852  if (!prop.GetUnit().CompareTo("ps") || !prop.GetUnit().CompareTo("picosec"))
1853  dval = prop.GetValue()/1.0e3;
1854  else if (!prop.GetUnit().CompareTo("ns") || !prop.GetUnit().CompareTo("nanosec"))
1855  dval = prop.GetValue();
1856  else if (!prop.GetUnit().CompareTo("us") || !prop.GetUnit().CompareTo("microsec"))
1857  dval = prop.GetValue()*1.0e3;
1858  else
1859  std::cerr << std::endl << ">> PRunDataHandler::ReadRootFile: **ERROR** Found unrecognized Time Resolution unit: " << prop.GetUnit() << std::endl;
1860  runData.SetTimeResolution(dval);
1861  }
1862 
1863  header->Get("RunInfo/RedGreen Offsets", ivec, ok);
1864  if (ok) {
1865  // check if any2many is used and a group histo list is defined, if NOT, only take the 0-offset data!
1866  if (fAny2ManyInfo) { // i.e. any2many is called
1867  if (fAny2ManyInfo->groupHistoList.size() == 0) { // NO group list defined -> use only the 0-offset data
1868  redGreenOffsets.push_back(0);
1869  } else { // group list defined
1870  // make sure that the group list elements is a subset of present RedGreen offsets
1871  Bool_t found = false;
1872  for (UInt_t i=0; i<fAny2ManyInfo->groupHistoList.size(); i++) {
1873  found = false;
1874  for (UInt_t j=0; j<ivec.size(); j++) {
1875  if (fAny2ManyInfo->groupHistoList[i] == ivec[j])
1876  found = true;
1877  }
1878  if (!found) {
1879  std::cerr << std::endl << ">> PRunDataHandler::ReadRootFile: **ERROR** requested histo group " << fAny2ManyInfo->groupHistoList[i];
1880  std::cerr << std::endl << ">> which is NOT present in the data file." << std::endl;
1881  return false;
1882  }
1883  }
1884  // found all requested histo groups, hence stuff it to the right places
1885  redGreenOffsets = fAny2ManyInfo->groupHistoList;
1887  }
1888  } else { // not any2many, i.e. musrfit, musrview, ...
1889  redGreenOffsets = ivec;
1890  runData.SetRedGreenOffset(ivec);
1891  }
1892  }
1893 
1894  // check further for LEM specific stuff in RunInfo
1895 
1896  header->Get("RunInfo/Moderator HV", prop, ok);
1897  if (ok)
1898  runData.SetTransport(prop.GetValue());
1899 
1900  header->Get("RunInfo/Implantation Energy", prop, ok);
1901  if (ok)
1902  runData.SetEnergy(prop.GetValue());
1903 
1904  // read the SampleEnvironmentInfo
1905 
1906  header->Get("SampleEnvironmentInfo/Cryo", str, ok);
1907  if (ok)
1908  runData.SetCryoName(str);
1909 
1910  // read the MagneticFieldEnvironmentInfo
1911 
1912  header->Get("MagneticFieldEnvironmentInfo/Magnet Name", str, ok);
1913  if (ok)
1914  runData.SetMagnetName(str);
1915 
1916  // read the BeamlineInfo
1917 
1918  header->Get("BeamlineInfo/Name", str, ok);
1919  if (ok)
1920  runData.SetBeamline(str);
1921 
1922  // read run summary to obtain ring anode HV values
1923  TObjArray *runSummary = dynamic_cast<TObjArray*>(folder->FindObjectAny("RunSummary"));
1924 
1925  // check if run summary is valid
1926  if (!runSummary) {
1927  std::cout << std::endl << "**INFO** Couldn't obtain run summary info from ROOT file " << fRunPathName.Data() << std::endl;
1928  // this is not fatal... only RA-HV values are not available
1929  } else { // it follows a (at least) little bit strange extraction of the RA values from Thomas' TObjArray...
1930  //streaming of a ASCII-file would be more easy
1931  TString s;
1932  TObjArrayIter summIter(runSummary);
1933  TObjString *os(dynamic_cast<TObjString*>(summIter.Next()));
1934  TObjArray *oa(nullptr);
1935  TObjString *objTok(nullptr);
1936  while (os != nullptr) {
1937  s = os->GetString();
1938  // will put four parallel if's since it may be that more than one RA-values are on one line
1939  if (s.Contains("RA-L")) {
1940  oa = s.Tokenize(" ");
1941  TObjArrayIter lineIter(oa);
1942  objTok = dynamic_cast<TObjString*>(lineIter.Next());
1943  while (objTok != nullptr) {
1944  if (!objTok->GetString().CompareTo("RA-L")) {
1945  objTok = dynamic_cast<TObjString*>(lineIter.Next()); // "="
1946  if ((objTok != nullptr) && !objTok->GetString().CompareTo("=")) {
1947  objTok = dynamic_cast<TObjString*>(lineIter.Next()); // HV value
1948  runData.SetRingAnode(0, objTok->GetString().Atof()); // fill RA-R value into the runData structure
1949  break;
1950  }
1951  }
1952  objTok = dynamic_cast<TObjString*>(lineIter.Next()); // next token...
1953  }
1954  // clean up
1955  if (oa) {
1956  delete oa;
1957  oa = nullptr;
1958  }
1959  }
1960 
1961  if (s.Contains("RA-R")) {
1962  oa = s.Tokenize(" ");
1963  TObjArrayIter lineIter(oa);
1964  objTok = dynamic_cast<TObjString*>(lineIter.Next());
1965  while (objTok != nullptr){
1966  if (!objTok->GetString().CompareTo("RA-R")) {
1967  objTok = dynamic_cast<TObjString*>(lineIter.Next()); // "="
1968  if ((objTok != nullptr) && !objTok->GetString().CompareTo("=")) {
1969  objTok = dynamic_cast<TObjString*>(lineIter.Next()); // HV value
1970  runData.SetRingAnode(1, objTok->GetString().Atof()); // fill RA-R value into the runData structure
1971  break;
1972  }
1973  }
1974  objTok = dynamic_cast<TObjString*>(lineIter.Next()); // next token...
1975  }
1976  // clean up
1977  if (oa) {
1978  delete oa;
1979  oa = nullptr;
1980  }
1981  }
1982 
1983  if (s.Contains("RA-T")) {
1984  oa = s.Tokenize(" ");
1985  TObjArrayIter lineIter(oa);
1986  objTok = dynamic_cast<TObjString*>(lineIter.Next());
1987  while (objTok != nullptr){
1988  if (!objTok->GetString().CompareTo("RA-T")) {
1989  objTok = dynamic_cast<TObjString*>(lineIter.Next()); // "="
1990  if ((objTok != nullptr) && !objTok->GetString().CompareTo("=")) {
1991  objTok = dynamic_cast<TObjString*>(lineIter.Next()); // HV value
1992  runData.SetRingAnode(2, objTok->GetString().Atof()); // fill RA-T value into the runData structure
1993  break;
1994  }
1995  }
1996  objTok = dynamic_cast<TObjString*>(lineIter.Next()); // next token...
1997  }
1998  // clean up
1999  if (oa) {
2000  delete oa;
2001  oa = nullptr;
2002  }
2003  }
2004 
2005  if (s.Contains("RA-B")) {
2006  oa = s.Tokenize(" ");
2007  TObjArrayIter lineIter(oa);
2008  objTok = dynamic_cast<TObjString*>(lineIter.Next());
2009  while (objTok != nullptr){
2010  if (!objTok->GetString().CompareTo("RA-B")) {
2011  objTok = dynamic_cast<TObjString*>(lineIter.Next()); // "="
2012  if ((objTok != nullptr) && !objTok->GetString().CompareTo("=")) {
2013  objTok = dynamic_cast<TObjString*>(lineIter.Next()); // HV value
2014  runData.SetRingAnode(3, objTok->GetString().Atof()); // fill RA-B value into the runData structure
2015  break;
2016  }
2017  }
2018  objTok = dynamic_cast<TObjString*>(lineIter.Next()); // next token...
2019  }
2020  // clean up
2021  if (oa) {
2022  delete oa;
2023  oa = nullptr;
2024  }
2025  }
2026 
2027  os = dynamic_cast<TObjString*>(summIter.Next()); // next summary line...
2028  }
2029  }
2030 
2031  // read data ---------------------------------------------------------
2032 
2033  // check if histos folder is found
2034  f.GetObject("histos", folder);
2035  if (!folder) {
2036  std::cerr << std::endl << ">> PRunDataHandler::ReadRootFile: **ERROR** Couldn't obtain histos from " << fRunPathName.Data() << std::endl;
2037  f.Close();
2038  return false;
2039  }
2040 
2041  // get all the data
2042  for (UInt_t i=0; i<redGreenOffsets.size(); i++) {
2043  for (Int_t j=0; j<noOfHistos; j++) {
2044  str.Form("hDecay%03d", redGreenOffsets[i]+j+1);
2045  TH1F *histo = dynamic_cast<TH1F*>(folder->FindObjectAny(str.Data()));
2046  if (!histo) {
2047  std::cerr << std::endl << ">> PRunDataHandler::ReadRootFile: **ERROR** Couldn't get histo " << str;
2048  std::cerr << std::endl;
2049  f.Close();
2050  return false;
2051  }
2052 
2053  dataSet.Clear();
2054  dataSet.SetName(histo->GetTitle());
2055  dataSet.SetHistoNo(redGreenOffsets[i]+j+1);
2056 
2057  // get detector info
2058  path.Form("DetectorInfo/Detector%03d/", redGreenOffsets[i]+j+1);
2059  pathName = path + "Time Zero Bin";
2060  header->Get(pathName, dval, ok);
2061  if (ok)
2062  dataSet.SetTimeZeroBin(dval);
2063  pathName = path + "First Good Bin";
2064  header->Get(pathName, ival, ok);
2065  if (ok)
2066  dataSet.SetFirstGoodBin(ival);
2067  pathName = path + "Last Good Bin";
2068  header->Get(pathName, ival, ok);
2069  if (ok)
2070  dataSet.SetLastGoodBin(ival);
2071  dataSet.SetTimeZeroBinEstimated(histo->GetMaximumBin());
2072 
2073  // fill data
2074  for (Int_t j=1; j<=histo->GetNbinsX(); j++) {
2075  histoData.push_back(histo->GetBinContent(j));
2076  }
2077  dataSet.SetData(histoData);
2078  runData.SetDataSet(dataSet);
2079 
2080  // clear histoData for the next histo
2081  histoData.clear();
2082  }
2083  }
2084  }
2085 
2086  f.Close();
2087 
2088  // keep run name
2089  runData.SetRunName(fRunName);
2090 
2091  // add run to the run list
2092  fData.push_back(runData);
2093 
2094  return true;
2095 }
2096 
2097 //--------------------------------------------------------------------------
2098 // ReadNexusFile (private)
2099 //--------------------------------------------------------------------------
2108 {
2109 #ifdef PNEXUS_ENABLED
2110  std::cout << std::endl << ">> PRunDataHandler::ReadNexusFile(): Will read nexus file " << fRunPathName.Data() << " ...";
2111 
2112  PDoubleVector histoData;
2113  PRawRunData runData;
2114  PRawRunDataSet dataSet;
2115  TString str;
2116  std::string sstr;
2117  Double_t dval;
2118  bool ok;
2119 
2120  std::unique_ptr<PNeXus> nxs_file = std::make_unique<PNeXus>(fRunPathName.Data());
2121  if (!nxs_file->IsValid()) {
2122  std::cerr << std::endl << ">> PRunDataHandler::ReadNexusFile(): Not a valid NeXus file.";
2123  std::cerr << std::endl << ">> Error Message: " << nxs_file->GetErrorMsg().data() << std::endl;
2124  return false;
2125  }
2126 
2127  if (nxs_file->GetIdfVersion() == 1) {
2128  if (!nxs_file->IsValid()) {
2129  std::cout << std::endl << "**ERROR** invalid NeXus IDF 2 version file found." << std::endl;
2130  return false;
2131  }
2132 
2133  // get header information
2134 
2135  // get/set laboratory
2136  str = TString(nxs_file->GetEntryIdf1()->GetLaboratory());
2137  runData.SetLaboratory(str);
2138 
2139  // get/set beamline
2140  str = TString(nxs_file->GetEntryIdf1()->GetBeamline());
2141  runData.SetBeamline(str);
2142 
2143  // get/set instrument
2144  str = TString(nxs_file->GetEntryIdf1()->GetInstrument()->GetName());
2145  runData.SetInstrument(str);
2146 
2147  // get/set run title
2148  str = TString(nxs_file->GetEntryIdf1()->GetTitle());
2149  runData.SetRunTitle(str);
2150 
2151  // get/set run number
2152  runData.SetRunNumber(nxs_file->GetEntryIdf1()->GetRunNumber());
2153 
2154  // get/set temperature
2155  dval = nxs_file->GetEntryIdf1()->GetSample()->GetPhysPropValue("temperature", ok);
2156  if (ok)
2157  runData.SetTemperature(0, dval, 0.0);
2158 
2159  // get/set field
2160  dval = nxs_file->GetEntryIdf1()->GetSample()->GetPhysPropValue("magnetic_field", ok);
2161  nxs_file->GetEntryIdf1()->GetSample()->GetPhysPropUnit("magnetic_field", sstr, ok);
2162  str = sstr;
2163  // since field has to be given in Gauss, check the units
2164  Double_t factor=1.0;
2165  if (!str.CompareTo("gauss", TString::kIgnoreCase))
2166  factor=1.0;
2167  else if (!str.CompareTo("tesla", TString::kIgnoreCase))
2168  factor=1.0e4;
2169  else
2170  factor=1.0;
2171  runData.SetField(factor*dval);
2172 
2173  // get/set implantation energy
2174  runData.SetEnergy(PMUSR_UNDEFINED);
2175 
2176  // get/set moderator HV
2177  runData.SetTransport(PMUSR_UNDEFINED);
2178 
2179  // get/set RA HV's (LEM specific)
2180  for (UInt_t i=0; i<4; i++)
2181  runData.SetRingAnode(i, PMUSR_UNDEFINED);
2182 
2183  // get/set setup
2184  runData.SetSetup(nxs_file->GetEntryIdf1()->GetNotes());
2185 
2186  // get/set sample
2187  runData.SetSample(nxs_file->GetEntryIdf1()->GetSample()->GetName());
2188 
2189  // get/set orientation
2190  runData.SetOrientation("??");
2191 
2192  // get/set time resolution (ns)
2193  runData.SetTimeResolution(nxs_file->GetEntryIdf1()->GetData()->GetTimeResolution("ns"));
2194 
2195  // get/set start/stop time
2196  sstr = nxs_file->GetEntryIdf1()->GetStartTime();
2197  str = sstr;
2198  TString date, time;
2199  SplitTimeDate(str, time, date, ok);
2200  if (ok) {
2201  runData.SetStartTime(time);
2202  runData.SetStartDate(date);
2203  }
2204 
2205  sstr = nxs_file->GetEntryIdf1()->GetStopTime();
2206  str = sstr;
2207  SplitTimeDate(str, time, date, ok);
2208  if (ok) {
2209  runData.SetStopTime(time);
2210  runData.SetStopDate(date);
2211  }
2212 
2213  // get/set t0, firstGoodBin, lastGoodBin
2214  std::vector<unsigned int> *t0 = nxs_file->GetEntryIdf1()->GetData()->GetT0s();
2215  std::vector<unsigned int> *fgb = nxs_file->GetEntryIdf1()->GetData()->GetFirstGoodBins();
2216  std::vector<unsigned int> *lgb = nxs_file->GetEntryIdf1()->GetData()->GetLastGoodBins();
2217 
2218  // get/set data
2219  std::vector<unsigned int> *pdata;
2220  unsigned int max=0, binMax=0;
2221  PDoubleVector data;
2222  for (UInt_t i=0; i<nxs_file->GetEntryIdf1()->GetData()->GetNoOfHistos(); i++) {
2223  pdata = nxs_file->GetEntryIdf1()->GetData()->GetHisto(i);
2224  for (UInt_t j=0; j<pdata->size(); j++) {
2225  data.push_back(pdata->at(j));
2226  if (pdata->at(j) > max) {
2227  max = pdata->at(j);
2228  binMax = j;
2229  }
2230  }
2231 
2232  // fill data set
2233  dataSet.Clear();
2234  dataSet.SetHistoNo(i+1); // i.e. histo numbers start with 1
2235  // set time zero bin
2236  if (i<t0->size())
2237  dataSet.SetTimeZeroBin(t0->at(i));
2238  else
2239  dataSet.SetTimeZeroBin(t0->at(0));
2240  // set time zero bin estimate
2241  dataSet.SetTimeZeroBinEstimated(binMax);
2242  // set first good bin
2243  if (i<fgb->size())
2244  dataSet.SetFirstGoodBin(fgb->at(i));
2245  else
2246  dataSet.SetFirstGoodBin(fgb->at(0));
2247  // set last good bin
2248  if (i<lgb->size())
2249  dataSet.SetFirstGoodBin(lgb->at(i));
2250  else
2251  dataSet.SetFirstGoodBin(lgb->at(0));
2252  dataSet.SetData(data);
2253 
2254  runData.SetDataSet(dataSet);
2255  data.clear();
2256  }
2257 
2258  // keep run name from the msr-file
2259  runData.SetRunName(fRunName);
2260 
2261  // keep the information
2262  fData.push_back(runData);
2263  } else if (nxs_file->GetIdfVersion() == 2) {
2264  if (!nxs_file->IsValid()) {
2265  std::cout << std::endl << "**ERROR** invalid NeXus IDF 2 version file found." << std::endl;
2266  return false;
2267  }
2268  // get header information
2269 
2270  // get/set laboratory
2271  str = TString(nxs_file->GetEntryIdf2()->GetInstrument()->GetSource()->GetName());
2272  runData.SetLaboratory(str);
2273 
2274  // get/set beamline
2275  str = TString(nxs_file->GetEntryIdf2()->GetInstrument()->GetName());
2276  runData.SetBeamline(str);
2277 
2278  // get/set instrument
2279  str = TString(nxs_file->GetEntryIdf2()->GetInstrument()->GetName());
2280  runData.SetInstrument(str);
2281 
2282  // get/set muon source
2283  str = TString(nxs_file->GetEntryIdf2()->GetInstrument()->GetSource()->GetType());
2284  runData.SetMuonSource(str);
2285 
2286  // get/set muon species
2287  str = TString(nxs_file->GetEntryIdf2()->GetInstrument()->GetSource()->GetProbe());
2288  runData.SetMuonSpecies(str);
2289 
2290  // get/set run title
2291  str = TString(nxs_file->GetEntryIdf2()->GetTitle());
2292  runData.SetRunTitle(str);
2293 
2294  // get/set run number
2295  runData.SetRunNumber(nxs_file->GetEntryIdf2()->GetRunNumber());
2296 
2297  // get/set temperature
2298  dval = nxs_file->GetEntryIdf2()->GetSample()->GetPhysPropValue("temperature", ok);
2299  if (ok)
2300  runData.SetTemperature(0, dval, 0.0);
2301 
2302  // get/set field
2303  dval = nxs_file->GetEntryIdf2()->GetSample()->GetPhysPropValue("magnetic_field", ok);
2304  nxs_file->GetEntryIdf2()->GetSample()->GetPhysPropUnit("magnetic_field", sstr, ok);
2305  str = sstr;
2306  // since field has to be given in Gauss, check the units
2307  Double_t factor=1.0;
2308  if (!str.CompareTo("gauss", TString::kIgnoreCase))
2309  factor=1.0;
2310  else if (!str.CompareTo("tesla", TString::kIgnoreCase))
2311  factor=1.0e4;
2312  else
2313  factor=1.0;
2314  runData.SetField(factor*dval);
2315 
2316  // get/set implantation energy
2317  runData.SetEnergy(PMUSR_UNDEFINED);
2318 
2319  // get/set moderator HV
2320  runData.SetTransport(PMUSR_UNDEFINED);
2321 
2322  // get/set RA HV's (LEM specific)
2323  for (UInt_t i=0; i<4; i++)
2324  runData.SetRingAnode(i, PMUSR_UNDEFINED);
2325 
2326  // get/set setup take NXsample/temperature_1_env and NXsample/magnetic_field_1_env
2327  sstr = nxs_file->GetEntryIdf2()->GetSample()->GetEnvironmentTemp() + std::string("/");
2328  sstr += nxs_file->GetEntryIdf2()->GetSample()->GetEnvironmentField();
2329  str = sstr;
2330  runData.SetSetup(str);
2331 
2332  // get/set sample
2333  runData.SetSample(nxs_file->GetEntryIdf2()->GetSample()->GetName());
2334 
2335  // get/set orientation
2336  runData.SetOrientation("??");
2337 
2338  // get/set time resolution (ns)
2339  runData.SetTimeResolution(nxs_file->GetEntryIdf2()->GetInstrument()->GetDetector()->GetTimeResolution("ns"));
2340 
2341  // get/set start/stop time
2342  sstr = nxs_file->GetEntryIdf2()->GetStartTime();
2343  str = sstr;
2344  TString date, time;
2345  SplitTimeDate(str, time, date, ok);
2346  if (ok) {
2347  runData.SetStartTime(time);
2348  runData.SetStartDate(date);
2349  }
2350 
2351  sstr = nxs_file->GetEntryIdf2()->GetStopTime();
2352  str = sstr;
2353  SplitTimeDate(str, time, date, ok);
2354  if (ok) {
2355  runData.SetStopTime(time);
2356  runData.SetStopDate(date);
2357  }
2358 
2359  // get/set data, t0, fgb, lgb
2360  PDoubleVector data;
2361  PRawRunDataSet dataSet;
2362  UInt_t histoNo = 0;
2363  Int_t ival;
2364  int *histos = nxs_file->GetEntryIdf2()->GetInstrument()->GetDetector()->GetHistos();
2365  if (nxs_file->GetEntryIdf2()->GetInstrument()->GetDetector()->GetNoOfPeriods() > 0) { // counts[][][]
2366  for (int i=0; i<nxs_file->GetEntryIdf2()->GetInstrument()->GetDetector()->GetNoOfPeriods(); i++) {
2367  for (int j=0; j<nxs_file->GetEntryIdf2()->GetInstrument()->GetDetector()->GetNoOfSpectra(); j++) {
2368  for (int k=0; k<nxs_file->GetEntryIdf2()->GetInstrument()->GetDetector()->GetNoOfBins(); k++) {
2369  data.push_back(*(histos+i*nxs_file->GetEntryIdf2()->GetInstrument()->GetDetector()->GetNoOfSpectra()+j*nxs_file->GetEntryIdf2()->GetInstrument()->GetDetector()->GetNoOfBins()+k));
2370  }
2371  dataSet.Clear();
2372  dataSet.SetHistoNo(++histoNo); // i.e. histo numbers start with 1
2373  // get t0
2374  ival = nxs_file->GetEntryIdf2()->GetInstrument()->GetDetector()->GetT0(i,j);
2375  if (ival == -1) // i.e. single value only
2376  ival = nxs_file->GetEntryIdf2()->GetInstrument()->GetDetector()->GetT0();
2377  dataSet.SetTimeZeroBin(ival);
2378  // get first good bin
2379  ival = nxs_file->GetEntryIdf2()->GetInstrument()->GetDetector()->GetFirstGoodBin(i,j);
2380  if (ival == -1) // i.e. single value only
2381  ival = nxs_file->GetEntryIdf2()->GetInstrument()->GetDetector()->GetFirstGoodBin();
2382  dataSet.SetFirstGoodBin(ival);
2383  // get last good bin
2384  ival = nxs_file->GetEntryIdf2()->GetInstrument()->GetDetector()->GetLastGoodBin(i,j);
2385  if (ival == -1) // i.e. single value only
2386  ival = nxs_file->GetEntryIdf2()->GetInstrument()->GetDetector()->GetLastGoodBin();
2387  dataSet.SetLastGoodBin(ival);
2388 
2389  dataSet.SetData(data);
2390  runData.SetDataSet(dataSet);
2391  data.clear();
2392  }
2393  }
2394  } else {
2395  if (nxs_file->GetEntryIdf2()->GetInstrument()->GetDetector()->GetNoOfSpectra() > 0) { // counts[][]
2396  for (int i=0; i<nxs_file->GetEntryIdf2()->GetInstrument()->GetDetector()->GetNoOfSpectra(); i++) {
2397  for (int j=0; j<nxs_file->GetEntryIdf2()->GetInstrument()->GetDetector()->GetNoOfBins(); j++) {
2398  data.push_back(*(histos+i*nxs_file->GetEntryIdf2()->GetInstrument()->GetDetector()->GetNoOfBins()+j));
2399  }
2400  dataSet.Clear();
2401  dataSet.SetHistoNo(++histoNo); // i.e. histo numbers start with 1
2402  // get t0
2403  ival = nxs_file->GetEntryIdf2()->GetInstrument()->GetDetector()->GetT0(-1,i);
2404  if (ival == -1) // i.e. single value only
2405  ival = nxs_file->GetEntryIdf2()->GetInstrument()->GetDetector()->GetT0();
2406  dataSet.SetTimeZeroBin(ival);
2407  // get first good bin
2408  ival = nxs_file->GetEntryIdf2()->GetInstrument()->GetDetector()->GetFirstGoodBin(-1,i);
2409  if (ival == -1) // i.e. single value only
2410  ival = nxs_file->GetEntryIdf2()->GetInstrument()->GetDetector()->GetFirstGoodBin();
2411  dataSet.SetFirstGoodBin(ival);
2412  // get last good bin
2413  ival = nxs_file->GetEntryIdf2()->GetInstrument()->GetDetector()->GetLastGoodBin(-1,i);
2414  if (ival == -1) // i.e. single value only
2415  ival = nxs_file->GetEntryIdf2()->GetInstrument()->GetDetector()->GetLastGoodBin();
2416  dataSet.SetLastGoodBin(ival);
2417 
2418  dataSet.SetData(data);
2419  runData.SetDataSet(dataSet);
2420  data.clear();
2421  }
2422  } else { // counts[]
2423  for (int i=0; i<nxs_file->GetEntryIdf2()->GetInstrument()->GetDetector()->GetNoOfBins(); i++) {
2424  data.push_back(*(histos+i));
2425  }
2426  dataSet.Clear();
2427  dataSet.SetHistoNo(++histoNo); // i.e. histo numbers start with 1
2428  // get t0
2429  ival = nxs_file->GetEntryIdf2()->GetInstrument()->GetDetector()->GetT0();
2430  dataSet.SetTimeZeroBin(ival);
2431  // get first good bin
2432  ival = nxs_file->GetEntryIdf2()->GetInstrument()->GetDetector()->GetFirstGoodBin();
2433  dataSet.SetFirstGoodBin(ival);
2434  // get last good bin
2435  ival = nxs_file->GetEntryIdf2()->GetInstrument()->GetDetector()->GetLastGoodBin();
2436  dataSet.SetLastGoodBin(ival);
2437 
2438  dataSet.SetData(data);
2439  runData.SetDataSet(dataSet);
2440  data.clear();
2441  }
2442  }
2443 
2444  // keep run name from the msr-file
2445  runData.SetRunName(fRunName);
2446 
2447  // keep the information
2448  fData.push_back(runData);
2449  } else {
2450  std::cout << std::endl << ">> PRunDataHandler::ReadNexusFile(): IDF version " << nxs_file->GetIdfVersion() << ", not implemented." << std::endl;
2451  }
2452 #else
2453  std::cout << std::endl << ">> PRunDataHandler::ReadNexusFile(): Sorry, not enabled at configuration level, i.e. --enable-NeXus when executing configure" << std::endl << std::endl;
2454 #endif
2455  return true;
2456 }
2457 
2458 //--------------------------------------------------------------------------
2459 // ReadWkmFile (private)
2460 //--------------------------------------------------------------------------
2470 {
2471  PDoubleVector histoData;
2472  PRawRunData runData;
2473 
2474  // open file
2475  std::ifstream f;
2476 
2477  // open wkm-file
2478  f.open(fRunPathName.Data(), std::ifstream::in);
2479  if (!f.is_open()) {
2480  std::cerr << std::endl << ">> PRunDataHandler::ReadWkmFile: **ERROR** Couldn't open run data (" << fRunPathName.Data() << ") file for reading, sorry ...";
2481  std::cerr << std::endl;
2482  return false;
2483  }
2484 
2485  // read header
2486  Bool_t headerInfo = true;
2487  Char_t instr[512];
2488  TString line, linecp;
2489  Double_t dval;
2490  Int_t ival;
2491  Bool_t ok;
2492  Int_t groups = 0, channels = 0;
2493 
2494  // skip leading empty lines
2495  do {
2496  f.getline(instr, sizeof(instr));
2497  line = TString(instr);
2498  if (!line.IsWhitespace())
2499  break;
2500  } while (!f.eof());
2501 
2502  // real header data should start here
2503  Ssiz_t idx;
2504  do {
2505  line = TString(instr);
2506  if (line.IsDigit()) { // end of header reached
2507  headerInfo = false;
2508  } else { // real stuff, hence filter data
2509  if (line.Contains("Title") || line.Contains("Titel")) {
2510  idx = line.Index(":");
2511  line.Replace(0, idx+1, nullptr, 0); // remove 'Title:'
2512  StripWhitespace(line);
2513  runData.SetRunTitle(line);
2514  } else if (line.Contains("Run:")) {
2515  idx = line.Index(":");
2516  line.Replace(0, idx+1, nullptr, 0); // remove 'Run:'
2517  StripWhitespace(line);
2518  ival = ToInt(line, ok);
2519  if (ok)
2520  runData.SetRunNumber(ival);
2521  } else if (line.Contains("Field")) {
2522  idx = line.Index(":");
2523  line.Replace(0, idx+1, nullptr, 0); // remove 'Field:'
2524  StripWhitespace(line);
2525  idx = line.Index("G"); // check if unit is given
2526  if (idx > 0) // unit is indeed given
2527  line.Resize(idx);
2528  dval = ToDouble(line, ok);
2529  if (ok)
2530  runData.SetField(dval);
2531  } else if (line.Contains("Setup")) {
2532  idx = line.Index(":");
2533  line.Replace(0, idx+1, nullptr, 0); // remove 'Setup:'
2534  StripWhitespace(line);
2535  runData.SetSetup(line);
2536  } else if (line.Contains("Temp:") || line.Contains("Temp(meas1):")) {
2537  linecp = line;
2538  idx = line.Index(":");
2539  line.Replace(0, idx+1, nullptr, 0); // remove 'Temp:'
2540  StripWhitespace(line);
2541  idx = line.Index("+/-"); // remove "+/- ..." part
2542  if (idx > 0)
2543  line.Resize(idx);
2544  idx = line.Index("K"); // remove "K ..." part
2545  if (idx > 0)
2546  line.Resize(idx);
2547  dval = ToDouble(line, ok);
2548  if (ok) {
2549  runData.SetTemperature(0, dval, 0.0);
2550  }
2551  idx = linecp.Index("+/-"); // get the error
2552  linecp.Replace(0, idx+3, nullptr, 0);
2553  StripWhitespace(linecp);
2554  dval = ToDouble(linecp, ok);
2555  if (ok) {
2556  runData.SetTempError(0, dval);
2557  }
2558  } else if (line.Contains("Temp(meas2):")) {
2559  linecp = line;
2560  idx = line.Index(":");
2561  line.Replace(0, idx+1, nullptr, 0); // remove 'Temp(meas2):'
2562  StripWhitespace(line);
2563  idx = line.Index("+/-"); // remove "+/- ..." part
2564  if (idx > 0)
2565  line.Resize(idx);
2566  idx = line.Index("K"); // remove "K ..." part
2567  if (idx > 0)
2568  line.Resize(idx);
2569  dval = ToDouble(line, ok);
2570  if (ok) {
2571  runData.SetTemperature(1, dval, 0.0);
2572  }
2573  idx = linecp.Index("+/-"); // get the error
2574  linecp.Replace(0, idx+3, nullptr, 0);
2575  StripWhitespace(linecp);
2576  dval = ToDouble(linecp, ok);
2577  if (ok) {
2578  runData.SetTempError(1, dval);
2579  }
2580  } else if (line.Contains("Groups")) {
2581  idx = line.Index(":");
2582  line.Replace(0, idx+1, nullptr, 0); // remove 'Groups:'
2583  StripWhitespace(line);
2584  ival = ToInt(line, ok);
2585  if (ok)
2586  groups = ival;
2587  } else if (line.Contains("Channels")) {
2588  idx = line.Index(":");
2589  line.Replace(0, idx+1, nullptr, 0); // remove 'Channels:'
2590  StripWhitespace(line);
2591  ival = ToInt(line, ok);
2592  if (ok)
2593  channels = ival;
2594  } else if (line.Contains("Resolution")) {
2595  idx = line.Index(":");
2596  line.Replace(0, idx+1, nullptr, 0); // remove 'Resolution:'
2597  StripWhitespace(line);
2598  dval = ToDouble(line, ok);
2599  if (ok)
2600  runData.SetTimeResolution(dval*1.0e3); // us -> ns
2601  }
2602  }
2603 
2604  if (headerInfo)
2605  f.getline(instr, sizeof(instr));
2606  } while (headerInfo && !f.eof());
2607 
2608  if ((groups == 0) || (channels == 0) || runData.GetTimeResolution() == 0.0) {
2609  std::cerr << std::endl << ">> PRunDataHandler::ReadWkmFile(): **ERROR** essential header informations are missing!";
2610  std::cerr << std::endl << ">> groups = " << groups;
2611  std::cerr << std::endl << ">> channels = " << channels;
2612  std::cerr << std::endl << ">> time resolution = " << runData.GetTimeResolution();
2613  std::cerr << std::endl;
2614  f.close();
2615  return false;
2616  }
2617 
2618  // read data ---------------------------------------------------------
2619  UInt_t group_counter = 0;
2620  Int_t val;
2621  TObjArray *tokens;
2622  TObjString *ostr;
2623  TString str;
2624  UInt_t histoNo = 0;
2625  PRawRunDataSet dataSet;
2626  do {
2627  // check if empty line, i.e. new group
2628  if (IsWhitespace(instr)) {
2629  dataSet.Clear();
2630  dataSet.SetHistoNo(++histoNo);
2631  dataSet.SetData(histoData);
2632 
2633  // get a T0 estimate
2634  Double_t maxVal = 0.0;
2635  Int_t maxBin = 0;
2636  for (UInt_t i=0; i<histoData.size(); i++) {
2637  if (histoData[i] > maxVal) {
2638  maxVal = histoData[i];
2639  maxBin = i;
2640  }
2641  }
2642  dataSet.SetTimeZeroBinEstimated(maxBin);
2643 
2644  runData.SetDataSet(dataSet);
2645 
2646  histoData.clear();
2647  group_counter++;
2648  } else {
2649  // extract values
2650  line = TString(instr);
2651  // check if line starts with character. Needed for RAL WKM format
2652  if (!line.IsDigit()) {
2653  f.getline(instr, sizeof(instr));
2654  continue;
2655  }
2656  tokens = line.Tokenize(" ");
2657 
2658  if (!tokens) { // no tokens found
2659  std::cerr << std::endl << ">> PRunDataHandler::ReadWkmFile(): **ERROR** while reading data: coulnd't tokenize run data.";
2660  return false;
2661  }
2662  for (Int_t i=0; i<tokens->GetEntries(); i++) {
2663  ostr = dynamic_cast<TObjString*>(tokens->At(i));
2664  str = ostr->GetString();
2665  val = ToInt(str, ok);
2666  if (ok) {
2667  histoData.push_back(val);
2668  } else {
2669  std::cerr << std::endl << ">> PRunDataHandler::ReadWkmFile(): **ERROR** while reading data: data line contains non-integer values.";
2670  // clean up
2671  delete tokens;
2672  return false;
2673  }
2674  }
2675  // clean up
2676  if (tokens) {
2677  delete tokens;
2678  tokens = nullptr;
2679  }
2680  }
2681 
2682  f.getline(instr, sizeof(instr));
2683 
2684  } while (!f.eof());
2685 
2686  // handle last line if present
2687  if (strlen(instr) != 0) {
2688  // extract values
2689  line = TString(instr);
2690  tokens = line.Tokenize(" ");
2691  if (!tokens) { // no tokens found
2692  std::cerr << std::endl << ">> PRunDataHandler::ReadWkmFile(): **ERROR** while reading data: coulnd't tokenize run data.";
2693  return false;
2694  }
2695  for (Int_t i=0; i<tokens->GetEntries(); i++) {
2696  ostr = dynamic_cast<TObjString*>(tokens->At(i));
2697  str = ostr->GetString();
2698  val = ToInt(str, ok);
2699  if (ok) {
2700  histoData.push_back(val);
2701  } else {
2702  std::cerr << std::endl << ">> PRunDataHandler::ReadWkmFile(): **ERROR** while reading data: data line contains non-integer values.";
2703  // clean up
2704  delete tokens;
2705  return false;
2706  }
2707  }
2708  // clean up
2709  if (tokens) {
2710  delete tokens;
2711  tokens = nullptr;
2712  }
2713  }
2714 
2715  // save the last histo if not empty
2716  if (histoData.size() > 0) {
2717  dataSet.Clear();
2718  dataSet.SetHistoNo(++histoNo);
2719  dataSet.SetData(histoData);
2720 
2721  // get a T0 estimate
2722  Double_t maxVal = 0.0;
2723  Int_t maxBin = 0;
2724  for (UInt_t i=0; i<histoData.size(); i++) {
2725  if (histoData[i] > maxVal) {
2726  maxVal = histoData[i];
2727  maxBin = i;
2728  }
2729  }
2730  dataSet.SetTimeZeroBinEstimated(maxBin);
2731 
2732  runData.SetDataSet(dataSet);
2733  histoData.clear();
2734  }
2735 
2736  // close file
2737  f.close();
2738 
2739  // check if all groups are found
2740  if (static_cast<Int_t>(runData.GetNoOfHistos()) != groups) {
2741  std::cerr << std::endl << ">> PRunDataHandler::ReadWkmFile(): **ERROR**";
2742  std::cerr << std::endl << ">> expected " << groups << " histos, but found " << runData.GetNoOfHistos();
2743  return false;
2744  }
2745 
2746  // check if all groups have enough channels
2747  for (UInt_t i=0; i<runData.GetNoOfHistos(); i++) {
2748  if (static_cast<Int_t>(runData.GetDataBin(i+1)->size()) != channels) {
2749  std::cerr << std::endl << ">> PRunDataHandler::ReadWkmFile(): **ERROR**";
2750  std::cerr << std::endl << ">> expected " << channels << " bins in histo " << i+1 << ", but found " << runData.GetDataBin(i)->size();
2751  return false;
2752  }
2753  }
2754 
2755  // keep run name
2756  runData.SetRunName(fRunName);
2757 
2758  // add run to the run list
2759  fData.push_back(runData);
2760 
2761  return true;
2762 }
2763 
2764 //--------------------------------------------------------------------------
2765 // ReadPsiBinFile (private)
2766 //--------------------------------------------------------------------------
2776 {
2777  MuSR_td_PSI_bin psiBin;
2778  Int_t status;
2779  Bool_t success;
2780 
2781  // read psi bin file
2782  status = psiBin.Read(fRunPathName.Data());
2783  switch (status) {
2784  case 0: // everything perfect
2785  success = true;
2786  break;
2787  case 1: // couldn't open file, or failed while reading the header
2788  std::cerr << std::endl << ">> **ERROR** couldn't open psi-bin file, or failed while reading the header";
2789  std::cerr << std::endl;
2790  success = false;
2791  break;
2792  case 2: // unsupported version of the data
2793  std::cerr << std::endl << ">> **ERROR** psi-bin file: unsupported version of the data";
2794  std::cerr << std::endl;
2795  success = false;
2796  break;
2797  case 3: // error when allocating data buffer
2798  std::cerr << std::endl << ">> **ERROR** psi-bin file: error when allocating data buffer";
2799  std::cerr << std::endl;
2800  success = false;
2801  break;
2802  case 4: // number of histograms/record not equals 1
2803  std::cerr << std::endl << ">> **ERROR** psi-bin file: number of histograms/record not equals 1";
2804  std::cerr << std::endl;
2805  success = false;
2806  break;
2807  default: // you never should have reached this point
2808  success = false;
2809  break;
2810  }
2811 
2812  // if any reading error happend, get out of here
2813  if (!success)
2814  return success;
2815 
2816  // fill necessary header informations
2817  PRawRunData runData;
2818  Double_t dval;
2819 
2820  // set file name
2821  Ssiz_t pos = fRunPathName.Last('/');
2822  TString fln(fRunPathName);
2823  fln.Remove(0, pos+1);
2824  runData.SetFileName(fln);
2825 
2826  // set laboratory
2827  runData.SetLaboratory("PSI");
2828 
2829  // filter from the file name the instrument
2830  TString instrument("n/a"), beamline("n/a");
2831  TString muonSource("n/a"), muonSpecies("n/a");
2832  if (fRunPathName.Contains("_gps_", TString::kIgnoreCase)) {
2833  instrument = "GPS";
2834  beamline = "piM3.2";
2835  muonSource = "continuous surface muon source";
2836  muonSpecies = "positive muons";
2837  } else if (fRunPathName.Contains("_ltf_", TString::kIgnoreCase)) {
2838  instrument = "LTF";
2839  beamline = "piM3.3";
2840  muonSource = "continuous surface muon source";
2841  muonSpecies = "positive muons";
2842  } else if (fRunPathName.Contains("_gpd_", TString::kIgnoreCase)) {
2843  instrument = "GPD";
2844  beamline = "muE1";
2845  muonSource = "continuous decay channel muon source";
2846  muonSpecies = "positive muons";
2847  } else if (fRunPathName.Contains("_dolly_", TString::kIgnoreCase)) {
2848  instrument = "DOLLY";
2849  beamline = "piE1";
2850  muonSource = "continuous surface muon source";
2851  muonSpecies = "positive muons";
2852  } else if (fRunPathName.Contains("_alc_", TString::kIgnoreCase)) {
2853  instrument = "ALC";
2854  beamline = "piE3";
2855  muonSource = "continuous surface muon source";
2856  muonSpecies = "positive muons";
2857  } else if (fRunPathName.Contains("_hifi_", TString::kIgnoreCase)) {
2858  instrument = "HIFI";
2859  beamline = "piE3";
2860  muonSource = "continuous surface muon source";
2861  muonSpecies = "positive muons";
2862  }
2863  runData.SetInstrument(instrument);
2864  runData.SetBeamline(beamline);
2865  runData.SetMuonSource(muonSource);
2866  runData.SetMuonSpecies(muonSpecies);
2867 
2868  // keep run name
2869  runData.SetRunName(fRunName);
2870  // get run title
2871  runData.SetRunTitle(TString(psiBin.GetComment().c_str())); // run title
2872  // get run number
2873  runData.SetRunNumber(psiBin.GetRunNumberInt());
2874  // get setup
2875  runData.SetSetup(TString(psiBin.GetComment().c_str()));
2876  // get sample
2877  runData.SetSample(TString(psiBin.GetSample().c_str()));
2878  // get orientation
2879  runData.SetOrientation(TString(psiBin.GetOrient().c_str()));
2880  // get comment
2881  runData.SetComment(TString(psiBin.GetComment().c_str()));
2882  // set LEM specific information to default value since it is not in the file and not used...
2883  runData.SetEnergy(PMUSR_UNDEFINED);
2884  runData.SetTransport(PMUSR_UNDEFINED);
2885  // get field
2886  Double_t scale = 0.0;
2887  if (psiBin.GetField().rfind("G") != std::string::npos)
2888  scale = 1.0;
2889  if (psiBin.GetField().rfind("T") != std::string::npos)
2890  scale = 1.0e4;
2891  status = sscanf(psiBin.GetField().c_str(), "%lf", &dval);
2892  if (status == 1)
2893  runData.SetField(scale*dval);
2894  // get temperature
2895  PDoubleVector tempVec(psiBin.GetTemperaturesVector());
2896  PDoubleVector tempDevVec(psiBin.GetDevTemperaturesVector());
2897  if ((tempVec.size() > 1) && (tempDevVec.size() > 1) && tempVec[0] && tempVec[1]) {
2898  // take only the first two values for now...
2899  //maybe that's not enough - e.g. in older GPD data I saw the "correct values in the second and third entry..."
2900  for (UInt_t i(0); i<2; i++) {
2901  runData.SetTemperature(i, tempVec[i], tempDevVec[i]);
2902  }
2903  tempVec.clear();
2904  tempDevVec.clear();
2905  } else {
2906  status = sscanf(psiBin.GetTemp().c_str(), "%lfK", &dval);
2907  if (status == 1)
2908  runData.SetTemperature(0, dval, 0.0);
2909  }
2910 
2911  // get time resolution (ns)
2912  runData.SetTimeResolution(psiBin.GetBinWidthNanoSec());
2913 
2914  // get start/stop time
2915  std::vector<std::string> sDateTime = psiBin.GetTimeStartVector();
2916  if (sDateTime.size() < 2) {
2917  std::cerr << std::endl << ">> **WARNING** psi-bin file: couldn't obtain run start date/time" << std::endl;
2918  }
2919  std::string date("");
2920  if (DateToISO8601(sDateTime[0], date)) {
2921  runData.SetStartDate(date);
2922  } else {
2923  std::cerr << std::endl << ">> **WARNING** failed to convert start date: " << sDateTime[0] << " into ISO 8601 date." << std::endl;
2924  runData.SetStartDate(sDateTime[0]);
2925  }
2926  runData.SetStartTime(sDateTime[1]);
2927  sDateTime.clear();
2928 
2929  sDateTime = psiBin.GetTimeStopVector();
2930  if (sDateTime.size() < 2) {
2931  std::cerr << std::endl << ">> **WARNING** psi-bin file: couldn't obtain run stop date/time" << std::endl;
2932  }
2933  date = std::string("");
2934  if (DateToISO8601(sDateTime[0], date)) {
2935  runData.SetStopDate(date);
2936  } else {
2937  std::cerr << std::endl << ">> **WARNING** failed to convert stop date: " << sDateTime[0] << " into ISO 8601 date." << std::endl;
2938  runData.SetStopDate(sDateTime[0]);
2939  }
2940  runData.SetStopTime(sDateTime[1]);
2941  sDateTime.clear();
2942 
2943  // get t0's
2944  PIntVector t0 = psiBin.GetT0Vector();
2945 
2946  if (t0.empty()) {
2947  std::cerr << std::endl << ">> **ERROR** psi-bin file: couldn't obtain any t0's";
2948  std::cerr << std::endl;
2949  return false;
2950  }
2951 
2952  // get first good bin
2953  PIntVector fgb = psiBin.GetFirstGoodVector();
2954  if (fgb.empty()) {
2955  std::cerr << std::endl << ">> **ERROR** psi-bin file: couldn't obtain any fgb's";
2956  std::cerr << std::endl;
2957  return false;
2958  }
2959 
2960  // get last good bin
2961  PIntVector lgb = psiBin.GetLastGoodVector();
2962  if (lgb.empty()) {
2963  std::cerr << std::endl << ">> **ERROR** psi-bin file: couldn't obtain any lgb's";
2964  std::cerr << std::endl;
2965  return false;
2966  }
2967 
2968  // fill raw data
2969  PRawRunDataSet dataSet;
2970  PDoubleVector histoData;
2971  std::vector<Int_t> histo;
2972  for (Int_t i=0; i<psiBin.GetNumberHistoInt(); i++) {
2973  histo = psiBin.GetHistoArrayInt(i);
2974  for (Int_t j=0; j<psiBin.GetHistoLengthBin(); j++) {
2975  histoData.push_back(histo[j]);
2976  }
2977 
2978  // estimate T0 from maximum of the data
2979  Double_t maxVal = 0.0;
2980  Int_t maxBin = 0;
2981  for (UInt_t j=0; j<histoData.size(); j++) {
2982  if (histoData[j] > maxVal) {
2983  maxVal = histoData[j];
2984  maxBin = j;
2985  }
2986  }
2987 
2988  dataSet.Clear();
2989  dataSet.SetName(psiBin.GetNameHisto(i).c_str());
2990  dataSet.SetHistoNo(i+1); // i.e. hist numbering starts at 1
2991  if (i < static_cast<Int_t>(t0.size()))
2992  dataSet.SetTimeZeroBin(t0[i]);
2993  dataSet.SetTimeZeroBinEstimated(maxBin);
2994  if (i < static_cast<Int_t>(fgb.size()))
2995  dataSet.SetFirstGoodBin(fgb[i]);
2996  if (i < static_cast<Int_t>(lgb.size()))
2997  dataSet.SetLastGoodBin(lgb[i]);
2998  dataSet.SetData(histoData);
2999 
3000  runData.SetDataSet(dataSet);
3001 
3002  histoData.clear();
3003  }
3004 
3005  // add run to the run list
3006  fData.push_back(runData);
3007 
3008  return success;
3009 }
3010 
3011 //--------------------------------------------------------------------------
3012 // ReadMudFile (private)
3013 //--------------------------------------------------------------------------
3022 {
3023  Int_t fh;
3024  UINT32 type, val;
3025  Int_t success;
3026  Char_t str[1024];
3027  Double_t dval;
3028 
3029  PRawRunData runData;
3030 
3031  fh = MUD_openRead((char*)fRunPathName.Data(), &type);
3032  if (fh == -1) {
3033  std::cerr << std::endl << ">> **ERROR** Couldn't open mud-file " << fRunPathName.Data() << ", sorry.";
3034  std::cerr << std::endl;
3035  return false;
3036  }
3037 
3038  // read necessary header information
3039 
3040  // keep run name
3041  runData.SetRunName(fRunName);
3042 
3043  // get/set the lab
3044  success = MUD_getLab( fh, str, sizeof(str) );
3045  if ( !success ) {
3046  std::cerr << std::endl << ">> **WARNING** Couldn't obtain the laboratory name of run " << fRunName.Data();
3047  std::cerr << std::endl;
3048  strcpy(str, "n/a");
3049  }
3050  runData.SetLaboratory(TString(str));
3051 
3052  // get/set the beamline
3053  success = MUD_getArea( fh, str, sizeof(str) );
3054  if ( !success ) {
3055  std::cerr << std::endl << ">> **WARNING** Couldn't obtain the beamline of run " << fRunName.Data();
3056  std::cerr << std::endl;
3057  strcpy(str, "n/a");
3058  }
3059  runData.SetBeamline(TString(str));
3060 
3061  // get/set the instrument
3062  success = MUD_getApparatus( fh, str, sizeof(str) );
3063  if ( !success ) {
3064  std::cerr << std::endl << ">> **WARNING** Couldn't obtain the instrument name of run " << fRunName.Data();
3065  std::cerr << std::endl;
3066  strcpy(str, "n/a");
3067  }
3068  runData.SetInstrument(TString(str));
3069 
3070  // get run title
3071  success = MUD_getTitle( fh, str, sizeof(str) );
3072  if ( !success ) {
3073  std::cerr << std::endl << ">> **WARNING** Couldn't obtain the run title of run " << fRunName.Data();
3074  std::cerr << std::endl;
3075  }
3076  runData.SetRunTitle(TString(str));
3077 
3078  // get run number
3079  success = MUD_getRunNumber( fh, &val );
3080  if (success) {
3081  runData.SetRunNumber(static_cast<Int_t>(val));
3082  }
3083 
3084  // get start/stop time of the run
3085  time_t tval;
3086  struct tm *dt;
3087  TString stime("");
3088  success = MUD_getTimeBegin( fh, (UINT32*)&tval );
3089  if (success) {
3090  runData.SetStartDateTime(static_cast<const time_t>(tval));
3091  dt = localtime((const time_t*)&tval);
3092 
3093  if (dt) {
3094  // start date
3095  strftime(str, sizeof(str), "%F", dt);
3096  stime = str;
3097  runData.SetStartDate(stime);
3098  // start time
3099  memset(str, 0, sizeof(str));
3100  strftime(str, sizeof(str), "%T", dt);
3101  stime = str;
3102  runData.SetStartTime(stime);
3103  } else {
3104  std::cerr << "PRunDataHandler::ReadMudFile: **WARNING** run start time readback wrong, will set it to 1900-01-01, 00:00:00" << std::endl;
3105  stime = "1900-01-01";
3106  runData.SetStartDate(stime);
3107  stime = "00:00:00";
3108  runData.SetStartTime(stime);
3109  }
3110  }
3111 
3112  stime = TString("");
3113  success = MUD_getTimeEnd( fh, (UINT32*)&tval );
3114  if (success) {
3115  runData.SetStopDateTime((const time_t)tval);
3116  dt = localtime((const time_t*)&tval);
3117 
3118  if (dt) {
3119  // stop date
3120  strftime(str, sizeof(str), "%F", dt);
3121  stime = str;
3122  runData.SetStopDate(stime);
3123  // stop time
3124  memset(str, 0, sizeof(str));
3125  strftime(str, sizeof(str), "%T", dt);
3126  stime = str;
3127  runData.SetStopTime(stime);
3128  } else {
3129  std::cerr << "PRunDataHandler::ReadMudFile: **WARNING** run stop time readback wrong, will set it to 1900-01-01, 00:00:00" << std::endl;
3130  stime = "1900-01-01";
3131  runData.SetStopDate(stime);
3132  stime = "00:00:00";
3133  runData.SetStopTime(stime);
3134  }
3135  }
3136 
3137  // get setup
3138  TString setup;
3139  REAL64 timeResMultiplier = 1.0e9; // Multiplier time resolution
3140  success = MUD_getLab( fh, str, sizeof(str) );
3141  if (success) {
3142  setup = TString(str) + TString("/");
3143  }
3144  success = MUD_getArea( fh, str, sizeof(str) );
3145  if (success) {
3146  setup += TString(str) + TString("/");
3147  if (TString(str) == "BNQR" || TString(str) == "BNMR") {
3148  std::cerr << "PRunDataHandler::ReadMudFile: **INFORMATION** this run was performed on " << str << std::endl;
3149  std::cerr << "PRunDataHandler::ReadMudFile: **INFORMATION** apply correction to time resolution" << std::endl;
3150  // identified BNMR/BNQR, correct time resolution.
3151  timeResMultiplier = 1.0e9;
3152  }
3153  }
3154  success = MUD_getApparatus( fh, str, sizeof(str) );
3155  if (success) {
3156  setup += TString(str) + TString("/");
3157  }
3158  success = MUD_getSample( fh, str, sizeof(str) );
3159  if (success) {
3160  setup += TString(str);
3161  runData.SetSample(str);
3162  }
3163  runData.SetSetup(setup);
3164 
3165  // set LEM specific information to default value since it is not in the file and not used...
3166  runData.SetEnergy(PMUSR_UNDEFINED);
3167  runData.SetTransport(PMUSR_UNDEFINED);
3168 
3169  // get field
3170  success = MUD_getField( fh, str, sizeof(str) );
3171  if (success) {
3172  success = sscanf(str, "%lf G", &dval);
3173  if (success == 1) {
3174  runData.SetField(dval);
3175  } else {
3176  runData.SetField(PMUSR_UNDEFINED);
3177  }
3178  } else {
3179  runData.SetField(PMUSR_UNDEFINED);
3180  }
3181 
3182  // get temperature
3183  success = MUD_getTemperature( fh, str, sizeof(str) );
3184  if (success) {
3185  success = sscanf(str, "%lf K", &dval);
3186  if (success == 1) {
3187  runData.SetTemperature(0, dval, 0.0);
3188  } else {
3189  runData.SetTemperature(0, PMUSR_UNDEFINED, 0.0);
3190  }
3191  } else {
3192  runData.SetTemperature(0, PMUSR_UNDEFINED, 0.0);
3193  }
3194 
3195  // get number of histogramms
3196  success = MUD_getHists(fh, &type, &val);
3197  if ( !success ) {
3198  std::cerr << std::endl << ">> **ERROR** Couldn't obtain the number of histograms of run " << fRunName.Data();
3199  std::cerr << std::endl;
3200  MUD_closeRead(fh);
3201  return false;
3202  }
3203  Int_t noOfHistos = static_cast<Int_t>(val);
3204 
3205  // get time resolution (ns)
3206  // check that time resolution is identical for all histograms
3207  // >> currently it is not forseen to handle histos with different time resolutions <<
3208  // >> perhaps this needs to be reconsidered later on <<
3209  REAL64 timeResolution = 0.0; // in seconds!!
3210  REAL64 lrval = 0.0;
3211  for (Int_t i=1; i<=noOfHistos; i++) {
3212  success = MUD_getHistSecondsPerBin( fh, i, &lrval );
3213  if (!success) {
3214  std::cerr << std::endl << ">> **ERROR** Couldn't obtain the time resolution of run " << fRunName.Data();
3215  std::cerr << std::endl << ">> which is fatal, sorry.";
3216  std::cerr << std::endl;
3217  MUD_closeRead(fh);
3218  return false;
3219  }
3220  if (i==1) {
3221  timeResolution = lrval;
3222  } else {
3223  if (lrval != timeResolution) {
3224  std::cerr << std::endl << ">> **ERROR** various time resolutions found in run " << fRunName.Data();
3225  std::cerr << std::endl << ">> this is currently not supported, sorry.";
3226  std::cerr << std::endl;
3227  MUD_closeRead(fh);
3228  return false;
3229  }
3230  }
3231  }
3232 
3233  runData.SetTimeResolution(static_cast<Double_t>(timeResolution) * timeResMultiplier); // s -> ns or s -> ms for bNMR
3234  // Other possibility:
3235  // Check if it is a bNMR run and fix it or check if "timeres" line
3236  // was introduced in the msr file
3237 
3238 
3239  // read histograms
3240  UINT32 *pData; // histo memory
3241  pData = nullptr;
3242  PDoubleVector histoData;
3243  PRawRunDataSet dataSet;
3244  UInt_t noOfBins;
3245 
3246  for (Int_t i=1; i<=noOfHistos; i++) {
3247  dataSet.Clear();
3248 
3249  dataSet.SetHistoNo(i);
3250 
3251  // get t0's
3252  success = MUD_getHistT0_Bin( fh, i, &val );
3253  if ( !success ) {
3254  std::cerr << std::endl << ">> **WARNING** Couldn't get t0 of histo " << i << " of run " << fRunName.Data();
3255  std::cerr << std::endl;
3256  }
3257  dataSet.SetTimeZeroBin(static_cast<Double_t>(val));
3258 
3259  // get bkg bins
3260  success = MUD_getHistBkgd1( fh, i, &val );
3261  if ( !success ) {
3262  std::cerr << std::endl << ">> **WARNING** Couldn't get bkg bin 1 of histo " << i << " of run " << fRunName.Data();
3263  std::cerr << std::endl;
3264  val = 0;
3265  }
3266  dataSet.SetFirstBkgBin(static_cast<Int_t>(val));
3267 
3268  success = MUD_getHistBkgd2( fh, i, &val );
3269  if ( !success ) {
3270  std::cerr << std::endl << ">> **WARNING** Couldn't get bkg bin 2 of histo " << i << " of run " << fRunName.Data();
3271  std::cerr << std::endl;
3272  val = 0;
3273  }
3274  dataSet.SetLastBkgBin(static_cast<Int_t>(val));
3275 
3276  // get good data bins
3277  success = MUD_getHistGoodBin1( fh, i, &val );
3278  if ( !success ) {
3279  std::cerr << std::endl << ">> **WARNING** Couldn't get good bin 1 of histo " << i << " of run " << fRunName.Data();
3280  std::cerr << std::endl;
3281  val = 0;
3282  }
3283  dataSet.SetFirstGoodBin(static_cast<Int_t>(val));
3284 
3285  success = MUD_getHistGoodBin2( fh, i, &val );
3286  if ( !success ) {
3287  std::cerr << std::endl << ">> **WARNING** Couldn't get good bin 2 of histo " << i << " of run " << fRunName.Data();
3288  std::cerr << std::endl;
3289  val = 0;
3290  }
3291  dataSet.SetLastGoodBin(static_cast<Int_t>(val));
3292 
3293  // get number of bins
3294  success = MUD_getHistNumBins( fh, i, &val );
3295  if ( !success ) {
3296  std::cerr << std::endl << ">> **ERROR** Couldn't get the number of bins of histo " << i << ".";
3297  std::cerr << std::endl << ">> This is fatal, sorry.";
3298  std::cerr << std::endl;
3299  MUD_closeRead( fh );
3300  return false;
3301  }
3302  noOfBins = static_cast<UInt_t>(val);
3303 
3304  pData = (UINT32*)malloc(noOfBins*sizeof(pData));
3305  if (pData == nullptr) {
3306  std::cerr << std::endl << ">> **ERROR** Couldn't allocate memory for data.";
3307  std::cerr << std::endl << ">> This is fatal, sorry.";
3308  std::cerr << std::endl;
3309  MUD_closeRead( fh );
3310  return false;
3311  }
3312 
3313  // get histogram
3314  success = MUD_getHistData( fh, i, pData );
3315  if ( !success ) {
3316  std::cerr << std::endl << ">> **ERROR** Couldn't get histo no " << i << ".";
3317  std::cerr << std::endl << ">> This is fatal, sorry.";
3318  std::cerr << std::endl;
3319  MUD_closeRead( fh );
3320  return false;
3321  }
3322 
3323  for (UInt_t j=0; j<noOfBins; j++) {
3324  histoData.push_back(pData[j]);
3325  }
3326  dataSet.SetData(histoData);
3327 
3328  // estimate T0 from maximum of the data
3329  Double_t maxVal = 0.0;
3330  Int_t maxBin = 0;
3331  for (UInt_t j=0; j<histoData.size(); j++) {
3332  if (histoData[j] > maxVal) {
3333  maxVal = histoData[j];
3334  maxBin = j;
3335  }
3336  }
3337  dataSet.SetTimeZeroBinEstimated(maxBin);
3338 
3339  runData.SetDataSet(dataSet);
3340 
3341  histoData.clear();
3342 
3343  free(pData);
3344  }
3345 
3346  MUD_closeRead(fh);
3347 
3348  // add run to the run list
3349  fData.push_back(runData);
3350 
3351  return true;
3352 }
3353 
3354 //--------------------------------------------------------------------------
3355 // ReadMduAsciiFile (private)
3356 //--------------------------------------------------------------------------
3384 {
3385  Bool_t success = true;
3386 
3387  // open file
3388  std::ifstream f;
3389 
3390  // open data-file
3391  f.open(fRunPathName.Data(), std::ifstream::in);
3392  if (!f.is_open()) {
3393  std::cerr << std::endl << ">> PRunDataHandler::ReadMduAsciiFile **ERROR** Couldn't open data file (" << fRunPathName.Data() << ") for reading, sorry ...";
3394  std::cerr << std::endl;
3395  return false;
3396  }
3397 
3398  PRawRunData runData;
3399 
3400  // keep run name
3401  runData.SetRunName(fRunName);
3402 
3403  Int_t lineNo = 0;
3404  Char_t instr[512];
3405  TString line, workStr;
3406  Bool_t headerTag = false;
3407  Bool_t dataTag = false;
3408  Int_t dataLineCounter = 0;
3409  TObjString *ostr;
3410  TObjArray *tokens = nullptr;
3411  TString str;
3412  Int_t groups = 0;
3413  Int_t channels = 0;
3414  Double_t dval = 0.0, unitScaling = 0.0;
3415  std::vector<PDoubleVector> data;
3416 
3417  while (!f.eof()) {
3418  f.getline(instr, sizeof(instr));
3419  line = TString(instr);
3420  lineNo++;
3421 
3422  // ignore comment lines
3423  if (line.BeginsWith("#") || line.BeginsWith("%"))
3424  continue;
3425 
3426  // ignore empty lines
3427  if (line.IsWhitespace())
3428  continue;
3429 
3430  // check if header tag
3431  workStr = line;
3432  workStr.Remove(TString::kLeading, ' '); // remove spaces from the begining
3433  if (workStr.BeginsWith("header", TString::kIgnoreCase)) {
3434  headerTag = true;
3435  dataTag = false;
3436  continue;
3437  }
3438 
3439  // check if data tag
3440  workStr = line;
3441  workStr.Remove(TString::kLeading, ' '); // remove spaces from the beining
3442  if (workStr.BeginsWith("data", TString::kIgnoreCase)) {
3443  headerTag = false;
3444  dataTag = true;
3445  continue;
3446  }
3447 
3448  if (headerTag) {
3449  workStr = line;
3450  workStr.Remove(TString::kLeading, ' '); // remove spaces from the beining
3451  if (workStr.BeginsWith("title:", TString::kIgnoreCase)) {
3452  runData.SetRunTitle(TString(workStr.Data()+workStr.First(":")+2));
3453  } else if (workStr.BeginsWith("field:", TString::kIgnoreCase)) {
3454  tokens = workStr.Tokenize(":("); // field: val (units)
3455  // check if expected number of tokens present
3456  if (tokens->GetEntries() != 3) {
3457  std::cerr << std::endl << ">> PRunDataHandler::ReadMduAsciiFile **ERROR** line no " << lineNo << ", invalid field entry in header.";
3458  std::cerr << std::endl << ">> " << line.Data();
3459  std::cerr << std::endl;
3460  success = false;
3461  break;
3462  }
3463  // check if field value is a number
3464  ostr = dynamic_cast<TObjString*>(tokens->At(1));
3465  if (ostr->GetString().IsFloat()) {
3466  dval = ostr->GetString().Atof();
3467  } else {
3468  std::cerr << std::endl << ">> PRunDataHandler::ReadMduAsciiFile **ERROR** line no " << lineNo << ", field value is not float/doulbe.";
3469  std::cerr << std::endl << ">> " << line.Data();
3470  std::cerr << std::endl;
3471  success = false;
3472  break;
3473  }
3474  // check units, accept (G), (T)
3475  ostr = dynamic_cast<TObjString*>(tokens->At(2));
3476  if (ostr->GetString().Contains("G"))
3477  unitScaling = 1.0;
3478  else if (ostr->GetString().Contains("T"))
3479  unitScaling = 1.0e4;
3480  else {
3481  std::cerr << std::endl << ">> PRunDataHandler::ReadMduAsciiFile **ERROR** line no " << lineNo << ", unkown field units.";
3482  std::cerr << std::endl << ">> " << line.Data();
3483  std::cerr << std::endl;
3484  success = false;
3485  break;
3486  }
3487  runData.SetField(dval*unitScaling);
3488 
3489  // clean up tokens
3490  if (tokens) {
3491  delete tokens;
3492  tokens = nullptr;
3493  }
3494  } else if (workStr.BeginsWith("temp:", TString::kIgnoreCase)) {
3495  tokens = workStr.Tokenize(":("); // temp: val (units)
3496  // check if expected number of tokens present
3497  if (tokens->GetEntries() != 3) {
3498  std::cerr << std::endl << ">> PRunDataHandler::ReadMduAsciiFile **ERROR** line no " << lineNo << ", invalid temperatue entry in header.";
3499  std::cerr << std::endl << ">> " << line.Data();
3500  std::cerr << std::endl;
3501  success = false;
3502  break;
3503  }
3504  // check if field value is a number
3505  ostr = dynamic_cast<TObjString*>(tokens->At(1));
3506  if (ostr->GetString().IsFloat()) {
3507  dval = ostr->GetString().Atof();
3508  } else {
3509  std::cerr << std::endl << ">> PRunDataHandler::ReadMduAsciiFile **ERROR** line no " << lineNo << ", temperature value is not float/doulbe.";
3510  std::cerr << std::endl << ">> " << line.Data();
3511  std::cerr << std::endl;
3512  success = false;
3513  break;
3514  }
3515  runData.SetTemperature(0, dval, 0.0);
3516 
3517  // clean up tokens
3518  if (tokens) {
3519  delete tokens;
3520  tokens = nullptr;
3521  }
3522  } else if (workStr.BeginsWith("setup:", TString::kIgnoreCase)) {
3523  runData.SetSetup(TString(workStr.Data()+workStr.First(":")+2));
3524  } else if (workStr.BeginsWith("groups:", TString::kIgnoreCase)) {
3525  workStr = TString(workStr.Data()+workStr.First(":")+2);
3526  groups = workStr.Atoi();
3527  if (groups == 0) {
3528  std::cerr << std::endl << ">> PRunDataHandler::ReadMduAsciiFile **ERROR** line no " << lineNo << ", groups is not a number or 0.";
3529  std::cerr << std::endl;
3530  success = false;
3531  break;
3532  }
3533  data.resize(groups);
3534  } else if (workStr.BeginsWith("channels:", TString::kIgnoreCase)) {
3535  workStr = TString(workStr.Data()+workStr.First(":")+2);
3536  channels = workStr.Atoi();
3537  if (channels == 0) {
3538  std::cerr << std::endl << ">> PRunDataHandler::ReadMduAsciiFile **ERROR** line no " << lineNo << ", channels is not a number or 0.";
3539  std::cerr << std::endl;
3540  success = false;
3541  break;
3542  }
3543  } else if (workStr.BeginsWith("resolution:", TString::kIgnoreCase)) {
3544  tokens = workStr.Tokenize(":("); // resolution: val (units)
3545  // check if expected number of tokens present
3546  if (tokens->GetEntries() != 3) {
3547  std::cerr << std::endl << ">> PRunDataHandler::ReadMduAsciiFile **ERROR** line no " << lineNo << ", invalid time resolution entry in header.";
3548  std::cerr << std::endl << line.Data();
3549  std::cerr << std::endl;
3550  success = false;
3551  break;
3552  }
3553  // check if timeresolution value is a number
3554  ostr = dynamic_cast<TObjString*>(tokens->At(1));
3555  if (ostr->GetString().IsFloat()) {
3556  dval = ostr->GetString().Atof();
3557  } else {
3558  std::cerr << std::endl << ">> PRunDataHandler::ReadMduAsciiFile **ERROR** line no " << lineNo << ", time resolution value is not float/doulbe.";
3559  std::cerr << std::endl << ">> " << line.Data();
3560  std::cerr << std::endl;
3561  success = false;
3562  break;
3563  }
3564  // check units, accept (fs), (ps), (ns), (us)
3565  ostr = dynamic_cast<TObjString*>(tokens->At(2));
3566  if (ostr->GetString().Contains("fs"))
3567  unitScaling = 1.0e-6;
3568  else if (ostr->GetString().Contains("ps"))
3569  unitScaling = 1.0e-3;
3570  else if (ostr->GetString().Contains("ns"))
3571  unitScaling = 1.0;
3572  else if (ostr->GetString().Contains("us"))
3573  unitScaling = 1.0e3;
3574  else {
3575  std::cerr << std::endl << ">> PRunDataHandler::ReadMduAsciiFile **ERROR** line no " << lineNo << ", unkown time resolution units.";
3576  std::cerr << std::endl << ">> " << line.Data();
3577  std::cerr << std::endl;
3578  success = false;
3579  break;
3580  }
3581  runData.SetTimeResolution(dval*unitScaling);
3582 
3583  // clean up tokens
3584  if (tokens) {
3585  delete tokens;
3586  tokens = nullptr;
3587  }
3588  } else { // error
3589  std::cerr << std::endl << ">> PRunDataHandler::ReadMduAsciiFile **ERROR** line no " << lineNo << ", illegal header line.";
3590  std::cerr << std::endl;
3591  success = false;
3592  break;
3593  }
3594  } else if (dataTag) {
3595  dataLineCounter++;
3596  tokens = line.Tokenize(" ,\t");
3597  // check if the number of data line entries is correct
3598  if (tokens->GetEntries() != groups+1) {
3599  std::cerr << std::endl << ">> PRunDataHandler::ReadMduAsciiFile **ERROR** found data line with a wrong data format, cannot be handled (line no " << lineNo << ")";
3600  std::cerr << std::endl << ">> line:";
3601  std::cerr << std::endl << ">> " << line.Data();
3602  std::cerr << std::endl;
3603  success = false;
3604  break;
3605  }
3606 
3607  for (Int_t i=1; i<tokens->GetEntries(); i++) {
3608  ostr = dynamic_cast<TObjString*>(tokens->At(i));
3609  data[i-1].push_back(ostr->GetString().Atof());
3610  }
3611 
3612  // clean up tokens
3613  if (tokens) {
3614  delete tokens;
3615  tokens = nullptr;
3616  }
3617  }
3618  }
3619 
3620  f.close();
3621 
3622  // keep data
3623  PRawRunDataSet dataSet;
3624  for (UInt_t i=0; i<data.size(); i++) {
3625  dataSet.Clear();
3626  dataSet.SetHistoNo(i);
3627  dataSet.SetData(data[i]);
3628 
3629  // estimate T0 from maximum of the data
3630  Double_t maxVal = 0.0;
3631  Int_t maxBin = 0;
3632  for (UInt_t j=0; j<data[i].size(); j++) {
3633  if (data[i][j] > maxVal) {
3634  maxVal = data[i][j];
3635  maxBin = j;
3636  }
3637  }
3638  dataSet.SetTimeZeroBinEstimated(maxBin);
3639  runData.SetDataSet(dataSet);
3640  }
3641 
3642  // clean up
3643  for (UInt_t i=0; i<data.size(); i++)
3644  data[i].clear();
3645  data.clear();
3646 
3647  if (dataLineCounter != channels) {
3648  std::cerr << std::endl << ">> PRunDataHandler::ReadMduAsciiFile **WARNING** found " << dataLineCounter << " data bins,";
3649  std::cerr << std::endl << ">> expected " << channels << " according to the header." << std::endl;
3650  }
3651 
3652  fData.push_back(runData);
3653 
3654  return success;
3655 }
3656 
3657 //--------------------------------------------------------------------------
3658 // ReadAsciiFile (private)
3659 //--------------------------------------------------------------------------
3697 {
3698  Bool_t success = true;
3699 
3700  // open file
3701  std::ifstream f;
3702 
3703  // open data-file
3704  f.open(fRunPathName.Data(), std::ifstream::in);
3705  if (!f.is_open()) {
3706  std::cerr << std::endl << ">> PRunDataHandler::ReadAsciiFile **ERROR** Couldn't open data file (" << fRunPathName.Data() << ") for reading, sorry ...";
3707  std::cerr << std::endl;
3708  return false;
3709  }
3710 
3711  PRawRunData runData;
3712 
3713  // init some stuff
3714  runData.fDataNonMusr.SetFromAscii(true);
3715  runData.fDataNonMusr.AppendLabel("??"); // x default label
3716  runData.fDataNonMusr.AppendLabel("??"); // y default label
3717 
3718  runData.SetRunName(fRunName); // keep the run name
3719 
3720  Int_t lineNo = 0;
3721  Char_t instr[512];
3722  TString line, workStr;
3723  Bool_t headerTag = false;
3724  Bool_t dataTag = false;
3725  Double_t x, y, ey;
3726  PDoubleVector xVec, exVec, yVec, eyVec;
3727 
3728  while (!f.eof()) {
3729  f.getline(instr, sizeof(instr));
3730  line = TString(instr);
3731  lineNo++;
3732 
3733  // check if comment line
3734  if (line.BeginsWith("#") || line.BeginsWith("%"))
3735  continue;
3736 
3737  // check if header tag
3738  workStr = line;
3739  workStr.Remove(TString::kLeading, ' '); // remove spaces from the begining
3740  if (workStr.BeginsWith("header", TString::kIgnoreCase)) {
3741  headerTag = true;
3742  dataTag = false;
3743  continue;
3744  }
3745 
3746  // check if data tag
3747  workStr = line;
3748  workStr.Remove(TString::kLeading, ' '); // remove spaces from the beining
3749  if (workStr.BeginsWith("data", TString::kIgnoreCase)) {
3750  headerTag = false;
3751  dataTag = true;
3752  continue;
3753  }
3754 
3755  if (headerTag) {
3756  if (line.IsWhitespace()) // ignore empty lines
3757  continue;
3758  workStr = line;
3759  workStr.Remove(TString::kLeading, ' '); // remove spaces from the beining
3760  if (workStr.BeginsWith("title:", TString::kIgnoreCase)) {
3761  runData.SetRunTitle(TString(workStr.Data()+workStr.First(":")+2));
3762  } else if (workStr.BeginsWith("setup:", TString::kIgnoreCase)) {
3763  runData.SetSetup(TString(workStr.Data()+workStr.First(":")+2));
3764  } else if (workStr.BeginsWith("field:", TString::kIgnoreCase)) {
3765  workStr = TString(workStr.Data()+workStr.First(":")+2);
3766  if (!workStr.IsFloat()) {
3767  std::cerr << std::endl << ">> PRunDataHandler::ReadAsciiFile **ERROR** line no " << lineNo << ", field is not a number.";
3768  std::cerr << std::endl;
3769  success = false;
3770  break;
3771  }
3772  runData.SetField(workStr.Atof());
3773  } else if (workStr.BeginsWith("x-axis-title:", TString::kIgnoreCase)) {
3774  runData.fDataNonMusr.SetLabel(0, TString(workStr.Data()+workStr.First(":")+2));
3775  } else if (workStr.BeginsWith("y-axis-title:", TString::kIgnoreCase)) {
3776  runData.fDataNonMusr.SetLabel(1, TString(workStr.Data()+workStr.First(":")+2));
3777  } else if (workStr.BeginsWith("temp:", TString::kIgnoreCase)) {
3778  workStr = TString(workStr.Data()+workStr.First(":")+2);
3779  if (!workStr.IsFloat()) {
3780  std::cerr << std::endl << ">> PRunDataHandler::ReadAsciiFile **ERROR** line no " << lineNo << ", temperature is not a number.";
3781  std::cerr << std::endl;
3782  success = false;
3783  break;
3784  }
3785  runData.SetTemperature(0, workStr.Atof(), 0.0);
3786  } else if (workStr.BeginsWith("energy:", TString::kIgnoreCase)) {
3787  workStr = TString(workStr.Data()+workStr.First(":")+2);
3788  if (!workStr.IsFloat()) {
3789  std::cerr << std::endl << ">> PRunDataHandler::ReadAsciiFile **ERROR** line no " << lineNo << ", energy is not a number.";
3790  std::cerr << std::endl;
3791  success = false;
3792  break;
3793  }
3794  runData.SetEnergy(workStr.Atof());
3795  runData.SetTransport(PMUSR_UNDEFINED); // just to initialize the variables to some "meaningful" value
3796  } else { // error
3797  std::cerr << std::endl << ">> PRunDataHandler::ReadAsciiFile **ERROR** line no " << lineNo << ", illegal header line.";
3798  std::cerr << std::endl;
3799  success = false;
3800  break;
3801  }
3802  } else if (dataTag) {
3803  if (line.IsWhitespace()) // ignore empty lines
3804  continue;
3805  TObjString *ostr;
3806  TObjArray *tokens;
3807 
3808  // Remove trailing end of line
3809  line.Remove(TString::kTrailing, '\r');
3810 
3811  // check if data have x, y [, error y] structure, and that x, y, and error y are numbers
3812  tokens = line.Tokenize(" ,\t");
3813  // check if the number of data line entries is 2 or 3
3814  if ((tokens->GetEntries() != 2) && (tokens->GetEntries() != 3)) {
3815  std::cerr << std::endl << ">> PRunDataHandler::ReadAsciiFile **ERROR** found data line with a structure different than \"x, y [, error y]\", cannot be handled (line no " << lineNo << ")";
3816  std::cerr << std::endl;
3817  success = false;
3818  break;
3819  }
3820 
3821  // get x
3822  ostr = dynamic_cast<TObjString*>(tokens->At(0));
3823  if (!ostr->GetString().IsFloat()) {
3824  std::cerr << std::endl << ">> PRunDataHandler::ReadAsciiFile **ERROR** line no " << lineNo << ": x = " << ostr->GetString().Data() << " is not a number, sorry.";
3825  std::cerr << std::endl;
3826  success = false;
3827  break;
3828  }
3829  x = ostr->GetString().Atof();
3830 
3831  // get y
3832  ostr = dynamic_cast<TObjString*>(tokens->At(1));
3833  if (!ostr->GetString().IsFloat()) {
3834  std::cerr << std::endl << ">> PRunDataHandler::ReadAsciiFile **ERROR** line no " << lineNo << ": y = " << ostr->GetString().Data() << " is not a number, sorry.";
3835  std::cerr << std::endl;
3836  success = false;
3837  break;
3838  }
3839  y = ostr->GetString().Atof();
3840 
3841  // get error y if present
3842  if (tokens->GetEntries() == 3) {
3843  ostr = dynamic_cast<TObjString*>(tokens->At(2));
3844  if (!ostr->GetString().IsFloat()) {
3845  std::cerr << std::endl << ">> PRunDataHandler::ReadAsciiFile **ERROR** line no " << lineNo << ": error y = " << ostr->GetString().Data() << " is not a number, sorry.";
3846  std::cerr << std::endl;
3847  success = false;
3848  break;
3849  }
3850  ey = ostr->GetString().Atof();
3851  } else {
3852  ey = 1.0;
3853  }
3854 
3855  // clean up tokens
3856  if (tokens) {
3857  delete tokens;
3858  tokens = nullptr;
3859  }
3860 
3861  // keep values
3862  xVec.push_back(x);
3863  exVec.push_back(1.0);
3864  yVec.push_back(y);
3865  eyVec.push_back(ey);
3866 
3867  } else {
3868  std::cerr << std::endl << ">> PRunDataHandler::ReadAsciiFile **ERROR** line no " << lineNo << " neither header nor data line. No idea how to handle it!";
3869  std::cerr << std::endl;
3870  success = false;
3871  break;
3872  }
3873  }
3874 
3875  f.close();
3876 
3877  // keep values
3878  runData.fDataNonMusr.AppendData(xVec);
3879  runData.fDataNonMusr.AppendErrData(exVec);
3880  runData.fDataNonMusr.AppendData(yVec);
3881  runData.fDataNonMusr.AppendErrData(eyVec);
3882 
3883  fData.push_back(runData);
3884 
3885  // clean up
3886  xVec.clear();
3887  exVec.clear();
3888  yVec.clear();
3889  eyVec.clear();
3890 
3891  return success;
3892 }
3893 
3894 //--------------------------------------------------------------------------
3895 // ReadDBFile (private)
3896 //--------------------------------------------------------------------------
4008 {
4009  Bool_t success = true;
4010 
4011  // open file
4012  std::ifstream f;
4013 
4014  // open db-file
4015  f.open(fRunPathName.Data(), std::ifstream::in);
4016  if (!f.is_open()) {
4017  std::cerr << std::endl << ">> PRunDataHandler::ReadDBFile **ERROR** Couldn't open data file (" << fRunPathName.Data() << ") for reading, sorry ...";
4018  std::cerr << std::endl;
4019  return false;
4020  }
4021 
4022  PRawRunData runData;
4023 
4024  runData.fDataNonMusr.SetFromAscii(false);
4025 
4026  Int_t lineNo = 0;
4027  Int_t idx;
4028  Int_t dbTag = -1;
4029  Char_t instr[512];
4030  TString line, workStr;
4031  Double_t val;
4032  Bool_t firstData = true; // needed as a switch to check in which format the data are given.
4033  Bool_t labelledFormat = true; // flag showing if the data are given in row format, or as labelled format (see description above, default is labelled format)
4034  Bool_t dataTagsRead = false; // flag showing if the data tags are alread read
4035 
4036  // variables needed to tokenize strings
4037  TString tstr;
4038  TObjString *ostr;
4039  TObjArray *tokens = nullptr;
4040 
4041  while (!f.eof()) {
4042  // get next line from file
4043  f.getline(instr, sizeof(instr));
4044  line = TString(instr);
4045  lineNo++;
4046 
4047  // check if comment line
4048  if (line.BeginsWith("#") || line.BeginsWith("%"))
4049  continue;
4050 
4051  // ignore empty lines
4052  if (line.IsWhitespace())
4053  continue;
4054 
4055  // check for db specific tags
4056  workStr = line;
4057  workStr.Remove(TString::kLeading, ' '); // remove spaces from the begining
4058  if (workStr.BeginsWith("title", TString::kIgnoreCase)) {
4059  dbTag = 0;
4060  continue;
4061  } else if (workStr.BeginsWith("abstract", TString::kIgnoreCase)) {
4062  dbTag = 1;
4063  continue;
4064  } else if (workStr.BeginsWith("comments", TString::kIgnoreCase)) {
4065  dbTag = 2;
4066  continue;
4067  } else if (workStr.BeginsWith("label", TString::kIgnoreCase)) {
4068  dbTag = 3;
4069  continue;
4070  } else if (workStr.BeginsWith("data", TString::kIgnoreCase) && !dataTagsRead) {
4071  dataTagsRead = true;
4072  dbTag = 4;
4073 
4074  // filter out all data tags
4075  tokens = workStr.Tokenize(" ,\t");
4076  for (Int_t i=1; i<tokens->GetEntries(); i++) {
4077  ostr = dynamic_cast<TObjString*>(tokens->At(i));
4078  runData.fDataNonMusr.AppendDataTag(ostr->GetString());
4079  }
4080 
4081  // clean up tokens
4082  if (tokens) {
4083  delete tokens;
4084  tokens = nullptr;
4085  }
4086  continue;
4087  }
4088 
4089  switch (dbTag) {
4090  case 0: // TITLE
4091  runData.SetRunTitle(workStr);
4092  break;
4093  case 1: // ABSTRACT
4094  // nothing to be done for now
4095  break;
4096  case 2: // COMMENTS
4097  // nothing to be done for now
4098  break;
4099  case 3: // LABEL
4100  runData.fDataNonMusr.AppendLabel(workStr);
4101  break;
4102  case 4: // DATA
4103  // filter out potential start data tag
4104  if (workStr.BeginsWith("\\-e", TString::kIgnoreCase) ||
4105  workStr.BeginsWith("\\e", TString::kIgnoreCase) ||
4106  workStr.BeginsWith("/-e", TString::kIgnoreCase) ||
4107  workStr.BeginsWith("/e", TString::kIgnoreCase)) {
4108  continue;
4109  }
4110 
4111  // check if first real data line
4112  if (firstData) {
4113  // check if data are given just as rows are as labelled columns (see description above)
4114  tokens = workStr.Tokenize(",");
4115  ostr = dynamic_cast<TObjString*>(tokens->At(0));
4116  if (!ostr->GetString().IsFloat()) {
4117  labelledFormat = true;
4118  } else {
4119  labelledFormat = false;
4120  }
4121  // clean up tokens
4122  if (tokens) {
4123  delete tokens;
4124  tokens = nullptr;
4125  }
4126 
4127  // prepare data vector for use
4128  PDoubleVector dummy;
4129  for (UInt_t i=0; i<runData.fDataNonMusr.GetDataTags()->size(); i++) {
4130  runData.fDataNonMusr.AppendData(dummy);
4131  runData.fDataNonMusr.AppendErrData(dummy);
4132  }
4133 
4134  firstData = false;
4135  }
4136 
4137  if (labelledFormat) { // handle labelled formated data
4138  // check if run line
4139  const Char_t *str = workStr.Data();
4140  if (isdigit(str[0])) { // run line
4141  TString run("run");
4142  idx = GetDataTagIndex(run, runData.fDataNonMusr.GetDataTags());
4143  if (idx == -1) {
4144  std::cerr << std::endl << ">> PRunDataHandler::ReadDBFile **ERROR** in line no " << lineNo << ":";
4145  std::cerr << std::endl << ">> " << workStr.Data();
4146  std::cerr << std::endl << ">> found potential run data line without run data tag.";
4147  return false;
4148  }
4149  // split string in tokens
4150  tokens = workStr.Tokenize(","); // line has structure: runNo,,, runTitle
4151  ostr = dynamic_cast<TObjString*>(tokens->At(0));
4152  tstr = ostr->GetString();
4153  if (!tstr.IsFloat()) {
4154  std::cerr << std::endl << ">> PRunDataHandler::ReadDBFile **ERROR** in line no " << lineNo << ":";
4155  std::cerr << std::endl << ">> " << workStr.Data();
4156  std::cerr << std::endl << ">> Expected db-data line with structure: runNo,,, runTitle";
4157  std::cerr << std::endl << ">> runNo = " << tstr.Data() << ", seems to be not a number.";
4158  delete tokens;
4159  return false;
4160  }
4161  val = tstr.Atof();
4162  runData.fDataNonMusr.AppendSubData(idx, val);
4163  runData.fDataNonMusr.AppendSubErrData(idx, 1.0);
4164  } else { // tag = data line
4165  // remove all possible spaces
4166  workStr.ReplaceAll(" ", "");
4167  // split string in tokens
4168  tokens = workStr.Tokenize("=,"); // line has structure: tag = val,err1,err2,
4169  if (tokens->GetEntries() < 3) {
4170  std::cerr << std::endl << ">> PRunDataHandler::ReadDBFile **ERROR** in line no " << lineNo << ":";
4171  std::cerr << std::endl << ">> " << workStr.Data();
4172  std::cerr << std::endl << ">> Expected db-data line with structure: tag = val,err1,err2,\\";
4173  delete tokens;
4174  return false;
4175  }
4176  ostr = dynamic_cast<TObjString*>(tokens->At(0));
4177  tstr = ostr->GetString();
4178  idx = GetDataTagIndex(tstr, runData.fDataNonMusr.GetDataTags());
4179  if (idx == -1) {
4180  std::cerr << std::endl << ">> PRunDataHandler::ReadDBFile **ERROR** in line no " << lineNo << ":";
4181  std::cerr << std::endl << ">> " << workStr.Data();
4182  std::cerr << std::endl << ">> data tag error: " << tstr.Data() << " seems not present in the data tag list";
4183  delete tokens;
4184  return false;
4185  }
4186 
4187  switch (tokens->GetEntries()) {
4188  case 3: // tag = val,,,
4189  ostr = dynamic_cast<TObjString*>(tokens->At(1));
4190  tstr = ostr->GetString();
4191  if (!tstr.IsFloat()) {
4192  std::cerr << std::endl << ">> PRunDataHandler::ReadDBFile **ERROR** in line no " << lineNo << ":";
4193  std::cerr << std::endl << ">> " << workStr.Data();
4194  std::cerr << std::endl << ">> Expected db-data line with structure: tag = val,err1,err2,\\";
4195  std::cerr << std::endl << ">> val = " << tstr.Data() << ", seems to be not a number.";
4196  delete tokens;
4197  return false;
4198  }
4199  val = tstr.Atof();
4200  runData.fDataNonMusr.AppendSubData(idx, val);
4201  runData.fDataNonMusr.AppendSubErrData(idx, 1.0);
4202  break;
4203  case 4: // tag = val,err,,
4204  case 5: // tag = val,err1,err2,
4205  // handle val
4206  ostr = dynamic_cast<TObjString*>(tokens->At(1));
4207  tstr = ostr->GetString();
4208  if (!tstr.IsFloat()) {
4209  std::cerr << std::endl << ">> PRunDataHandler::ReadDBFile **ERROR** in line no " << lineNo << ":";
4210  std::cerr << std::endl << ">> " << workStr.Data();
4211  std::cerr << std::endl << ">> Expected db-data line with structure: tag = val,err1,err2,\\";
4212  std::cerr << std::endl << ">> val = " << tstr.Data() << ", seems to be not a number.";
4213  delete tokens;
4214  return false;
4215  }
4216  val = tstr.Atof();
4217  runData.fDataNonMusr.AppendSubData(idx, val);
4218  // handle err1 (err2 will be ignored for the time being)
4219  ostr = dynamic_cast<TObjString*>(tokens->At(2));
4220  tstr = ostr->GetString();
4221  if (!tstr.IsFloat()) {
4222  std::cerr << std::endl << ">> PRunDataHandler::ReadDBFile **ERROR** in line no " << lineNo << ":";
4223  std::cerr << std::endl << ">> " << workStr.Data();
4224  std::cerr << std::endl << ">> Expected db-data line with structure: tag = val,err1,err2,\\";
4225  std::cerr << std::endl << ">> err1 = " << tstr.Data() << ", seems to be not a number.";
4226  delete tokens;
4227  return false;
4228  }
4229  val = tstr.Atof();
4230  runData.fDataNonMusr.AppendSubErrData(idx, val);
4231  break;
4232  default:
4233  std::cerr << std::endl << ">> PRunDataHandler::ReadDBFile **ERROR** in line no " << lineNo << ":";
4234  std::cerr << std::endl << ">> " << workStr.Data();
4235  std::cerr << std::endl << ">> Expected db-data line with structure: tag = val,err1,err2,\\";
4236  delete tokens;
4237  return false;
4238  }
4239  }
4240 
4241  } else { // handle row formated data
4242  // split string in tokens
4243  tokens = workStr.Tokenize(","); // line has structure: val1, err11, err12, ..., valn, errn1, errn2, runNo, , , , runTitle
4244  if (tokens->GetEntries() != static_cast<Int_t>(3*runData.fDataNonMusr.GetDataTags()->size()+1)) {
4245  std::cerr << std::endl << ">> PRunDataHandler::ReadDBFile **ERROR** in line no " << lineNo << ":";
4246  std::cerr << std::endl << ">> " << workStr.Data();
4247  std::cerr << std::endl << ">> Expected db-data line with structure: val1, err11, err12, ..., valn, errn1, errn2, runNo, , , , runTitle";
4248  std::cerr << std::endl << ">> found = " << tokens->GetEntries() << " tokens, however expected " << 3*runData.fDataNonMusr.GetDataTags()->size()+1;
4249  std::cerr << std::endl << ">> Perhaps there are commas without space inbetween, like 12.3,, 3.2,...";
4250  delete tokens;
4251  return false;
4252  }
4253  // extract data
4254  Int_t j=0;
4255  for (Int_t i=0; i<tokens->GetEntries()-1; i+=3) {
4256  // handle value
4257  ostr = dynamic_cast<TObjString*>(tokens->At(i));
4258  tstr = ostr->GetString();
4259  if (!tstr.IsFloat()) {
4260  std::cerr << std::endl << ">> PRunDataHandler::ReadDBFile **ERROR** in line no " << lineNo << ":";
4261  std::cerr << std::endl << ">> " << workStr.Data();
4262  std::cerr << std::endl << ">> Expected db-data line with structure: val1, err11, err12, ..., valn, errn1, errn2, runNo, , , , runTitle";
4263  std::cerr << std::endl << ">> value=" << tstr.Data() << " seems not to be a number";
4264  delete tokens;
4265  return false;
4266  }
4267  runData.fDataNonMusr.AppendSubData(j, tstr.Atof());
4268 
4269  // handle 1st error if present (2nd will be ignored for now)
4270  ostr = dynamic_cast<TObjString*>(tokens->At(i+1));
4271  tstr = ostr->GetString();
4272  if (tstr.IsWhitespace()) {
4273  runData.fDataNonMusr.AppendSubErrData(j, 1.0);
4274  } else if (tstr.IsFloat()) {
4275  runData.fDataNonMusr.AppendSubErrData(j, tstr.Atof());
4276  } else {
4277  std::cerr << std::endl << ">> PRunDataHandler::ReadDBFile **ERROR** in line no " << lineNo << ":";
4278  std::cerr << std::endl << ">> " << workStr.Data();
4279  std::cerr << std::endl << ">> Expected db-data line with structure: val1, err11, err12, ..., valn, errn1, errn2, runNo, , , , runTitle";
4280  std::cerr << std::endl << ">> error1=" << tstr.Data() << " seems not to be a number";
4281  delete tokens;
4282  return false;
4283  }
4284  j++;
4285  }
4286  }
4287  break;
4288  default:
4289  break;
4290  }
4291  }
4292 
4293  f.close();
4294 
4295  // check that the number of labels == the number of data tags
4296  if (runData.fDataNonMusr.GetLabels()->size() != runData.fDataNonMusr.GetDataTags()->size()) {
4297  std::cerr << std::endl << ">> PRunDataHandler::ReadDBFile **ERROR**";
4298  std::cerr << std::endl << ">> number of LABELS found = " << runData.fDataNonMusr.GetLabels()->size();
4299  std::cerr << std::endl << ">> number of Data tags found = " << runData.fDataNonMusr.GetDataTags()->size();
4300  std::cerr << std::endl << ">> They have to be equal!!";
4301  if (tokens) {
4302  delete tokens;
4303  tokens = nullptr;
4304  }
4305  return false;
4306  }
4307 
4308  // check if all vectors have the same size
4309  for (UInt_t i=1; i<runData.fDataNonMusr.GetData()->size(); i++) {
4310  if (runData.fDataNonMusr.GetData()->at(i).size() != runData.fDataNonMusr.GetData()->at(i-1).size()) {
4311  std::cerr << std::endl << ">> PRunDataHandler::ReadDBFile **ERROR** in line no " << lineNo;
4312  std::cerr << std::endl << ">> label: " << runData.fDataNonMusr.GetDataTags()->at(i-1).Data() << ", number data elements = " << runData.fDataNonMusr.GetData()->at(i-1).size();
4313  std::cerr << std::endl << ">> label: " << runData.fDataNonMusr.GetDataTags()->at(i).Data() << ", number data elements = " << runData.fDataNonMusr.GetData()->at(i).size();
4314  std::cerr << std::endl << ">> They have to be equal!!";
4315  success = false;
4316  break;
4317  }
4318  if (runData.fDataNonMusr.GetErrData()->at(i).size() != runData.fDataNonMusr.GetErrData()->at(i-1).size()) {
4319  std::cerr << std::endl << ">> PRunDataHandler::ReadDBFile **ERROR** in line no " << lineNo;
4320  std::cerr << std::endl << ">> label: " << runData.fDataNonMusr.GetDataTags()->at(i-1).Data() << ", number data elements = " << runData.fDataNonMusr.GetData()->at(i-1).size();
4321  std::cerr << std::endl << ">> label: " << runData.fDataNonMusr.GetDataTags()->at(i).Data() << ", number error data elements = " << runData.fDataNonMusr.GetErrData()->at(i).size();
4322  std::cerr << std::endl << ">> They have to be equal!!";
4323  success = false;
4324  break;
4325  }
4326  }
4327 
4328  // clean up tokens
4329  if (tokens) {
4330  delete tokens;
4331  tokens = nullptr;
4332  }
4333 
4334  // keep run name
4335  runData.SetRunName(fRunName);
4336 
4337  fData.push_back(runData);
4338 
4339  return success;
4340 }
4341 
4342 //--------------------------------------------------------------------------
4343 // ReadDatFile (private)
4344 //--------------------------------------------------------------------------
4354 {
4355  Bool_t success = true;
4356 
4357  // open file
4358  std::ifstream f;
4359 
4360  // open db-file
4361  f.open(fRunPathName.Data(), std::ifstream::in);
4362  if (!f.is_open()) {
4363  std::cerr << std::endl << ">> PRunDataHandler::ReadDatFile **ERROR** Couldn't open data file (" << fRunPathName.Data() << ") for reading, sorry ...";
4364  std::cerr << std::endl;
4365  return false;
4366  }
4367 
4368  PRawRunData runData;
4369  runData.fDataNonMusr.SetFromAscii(false);
4370 
4371  Int_t lineNo = 0;
4372  Char_t instr[4096];
4373  TString line;
4374  Bool_t headerInfo=true;
4375 
4376  // variables needed to tokenize strings
4377  TString tstr;
4378  TObjString *ostr;
4379  TObjArray *tokens = nullptr;
4380 
4381  UInt_t noOfDataSets = 0, noOfEntries = 0;
4382  PBoolVector isData;
4383  Double_t dval;
4384 
4385  while (!f.eof()) {
4386  // get next line from file
4387  f.getline(instr, sizeof(instr));
4388  line = TString(instr);
4389  lineNo++;
4390 
4391  // check if comment line
4392  if (line.BeginsWith("#") || line.BeginsWith("%"))
4393  continue;
4394 
4395  // ignore empty lines
4396  if (line.IsWhitespace())
4397  continue;
4398 
4399  tokens = line.Tokenize(" \t");
4400  if (tokens == nullptr) { // error
4401  std::cerr << std::endl << ">> PRunDataHandler::ReadDatFile **ERROR** couldn't tokenize the line, in lineNo: " << lineNo;
4402  std::cerr << std::endl << ">> line: '" << line << "'.";
4403  std::cerr << std::endl;
4404  return false;
4405  }
4406 
4407  // filter header information
4408  if (headerInfo) {
4409  headerInfo = false;
4410 
4411  // filter out all data tags: this labels are used in the msr-file to select the proper data set
4412  // for the dat-files, label and dataTag are the same
4413  noOfEntries = tokens->GetEntries();
4414  for (Int_t i=0; i<noOfEntries; i++) {
4415  ostr = dynamic_cast<TObjString*>(tokens->At(i));
4416  tstr = ostr->GetString();
4417  if (!tstr.EndsWith("Err", TString::kExact)) {
4418  noOfDataSets++;
4419  isData.push_back(true);
4420  runData.fDataNonMusr.AppendDataTag(tstr);
4421  runData.fDataNonMusr.AppendLabel(tstr);
4422  } else {
4423  isData.push_back(false);
4424  }
4425  }
4426  // set the size for the data
4427  runData.fDataNonMusr.SetSize(noOfDataSets);
4428  } else { // deal with data
4429  if (noOfEntries == 0) {
4430  std::cerr << std::endl << ">> PRunDataHandler::ReadDatFile **ERROR** header information is missing.";
4431  std::cerr << std::endl;
4432  return false;
4433  }
4434  if (tokens->GetEntries() != noOfEntries) { // error
4435  std::cerr << std::endl << ">> PRunDataHandler::ReadDatFile **ERROR** data set with wrong number of entries: " << tokens->GetEntries() << ", should be " << noOfEntries << ".";
4436  std::cerr << std::endl << ">> in line: " << lineNo;
4437  std::cerr << std::endl << ">> line: '" << line << "'.";
4438  std::cerr << std::endl;
4439  return false;
4440  }
4441  // fill data and dataErr sets
4442  UInt_t idx = 0;
4443  for (UInt_t i=0; i<noOfEntries; i++) {
4444  // 1st: check that entry is indeed a number
4445  ostr = dynamic_cast<TObjString*>(tokens->At(i));
4446  tstr = ostr->GetString();
4447  if (!tstr.IsFloat()) { // make sure it is a number
4448  std::cerr << std::endl << ">> PRunDataHandler::ReadDatFile **ERROR** data set entry is not a number: " << tstr.Data();
4449  std::cerr << std::endl << ">> in line: " << lineNo;
4450  std::cerr << std::endl;
4451  return false;
4452  }
4453  dval = tstr.Atof();
4454  if (isData[i]) {
4455  runData.fDataNonMusr.AppendSubData(idx, dval);
4456  idx++;
4457  } else { // error value
4458  if (isData[i-1] == 1) { // Err or PosErr hence keep it
4459  runData.fDataNonMusr.AppendSubErrData(idx-1, dval);
4460  }
4461  }
4462  }
4463  }
4464  // cleanup
4465  if (tokens) {
4466  delete tokens;
4467  tokens = nullptr;
4468  }
4469  }
4470 
4471  f.close();
4472 
4473  // got through all the data sets and if there is NO error vector set it to '0.0'
4474  for (UInt_t i=0; i<noOfDataSets; i++) {
4475  if (runData.fDataNonMusr.GetErrData()->at(i).size() == 0) {
4476  for (UInt_t j=0; j<runData.fDataNonMusr.GetData()->at(i).size(); j++) {
4477  runData.fDataNonMusr.AppendSubErrData(i, 0.0);
4478  }
4479  }
4480  }
4481 
4482  // keep run name
4483  runData.SetRunName(fRunName);
4484 
4485  fData.push_back(runData);
4486 
4487  return success;
4488 }
4489 
4490 //--------------------------------------------------------------------------
4491 // WriteMusrRootFile (private)
4492 //--------------------------------------------------------------------------
4503 {
4504  Bool_t ok = false;
4505  fln = GenerateOutputFileName(fln, ".root", ok);
4506  if (!ok)
4507  return false;
4508 
4510  std::cout << std::endl << ">> PRunDataHandler::WriteMusrRootFile(): writing a root data file (" << fln.Data() << ") ... " << std::endl;
4511 
4512  // generate data file
4513  TFolder *histosFolder;
4514  TFolder *decayAnaModule;
4515  TFolder *runHeader;
4516 
4517  histosFolder = gROOT->GetRootFolder()->AddFolder("histos", "Histograms");
4518  gROOT->GetListOfBrowsables()->Add(histosFolder, "histos");
4519  decayAnaModule = histosFolder->AddFolder("DecayAnaModule", "muSR decay histograms");
4520 
4521  runHeader = gROOT->GetRootFolder()->AddFolder("RunHeader", "MusrRoot Run Header Info");
4522  gROOT->GetListOfBrowsables()->Add(runHeader, "RunHeader");
4523  std::unique_ptr<TMusrRunHeader> header = std::make_unique<TMusrRunHeader>(true);
4524  gROOT->GetListOfBrowsables()->Add(runHeader, "RunHeader");
4525 
4526  // feed header info
4527  TString str, pathName;
4528  Int_t ival;
4529  Double_t dval[2];
4531 
4532  // feed RunInfo
4533  str = fData[0].GetGenericValidatorUrl()->Copy();
4534  header->Set("RunInfo/Generic Validator URL", str);
4535  str = fData[0].GetSpecificValidatorUrl()->Copy();
4536  header->Set("RunInfo/Specific Validator URL", str);
4537  str = fData[0].GetGenerator()->Copy();
4538  header->Set("RunInfo/Generator", str);
4539  str = fData[0].GetFileName()->Copy();
4540  header->Set("RunInfo/File Name", str);
4541  str = fData[0].GetRunTitle()->Copy();
4542  header->Set("RunInfo/Run Title", str);
4543  header->Set("RunInfo/Run Number", fData[0].GetRunNumber());
4544  str = fData[0].GetStartDate()->Copy() + " " + fData[0].GetStartTime()->Copy();
4545  header->Set("RunInfo/Run Start Time", str);
4546  str = fData[0].GetStopDate()->Copy() + " " + fData[0].GetStopTime()->Copy();
4547  header->Set("RunInfo/Run Stop Time", str);
4548  ival = fData[0].GetStopDateTime() - fData[0].GetStartDateTime();
4549  prop.Set("Run Duration", ival, "sec");
4550  header->Set("RunInfo/Run Duration", prop);
4551  str = fData[0].GetLaboratory()->Copy();
4552  header->Set("RunInfo/Laboratory", str);
4553  str = fData[0].GetInstrument()->Copy();
4554  header->Set("RunInfo/Instrument", str);
4555  dval[0] = fData[0].GetMuonBeamMomentum();
4556  prop.Set("Muon Beam Momentum", dval[0], "MeV/c");
4557  header->Set("RunInfo/Muon Beam Momentum", prop);
4558  str = fData[0].GetMuonSpecies()->Copy();
4559  header->Set("RunInfo/Muon Species", str);
4560  str = fData[0].GetMuonSource()->Copy();
4561  header->Set("RunInfo/Muon Source", str);
4562  str = fData[0].GetSetup()->Copy();
4563  header->Set("RunInfo/Setup", str);
4564  str = fData[0].GetComment()->Copy();
4565  header->Set("RunInfo/Comment", str);
4566  str = fData[0].GetSample()->Copy();
4567  header->Set("RunInfo/Sample Name", str);
4568  dval[0] = fData[0].GetTemperature(0);
4569  dval[1] = fData[0].GetTempError(0);
4570  prop.Set("Sample Temperature", MRH_UNDEFINED, dval[0], dval[1], "K");
4571  header->Set("RunInfo/Sample Temperature", prop);
4572  dval[0] = fData[0].GetField();
4573  prop.Set("Sample Magnetic Field", dval[0], "G");
4574  header->Set("RunInfo/Sample Magnetic Field", prop);
4575  header->Set("RunInfo/No of Histos", static_cast<Int_t>(fData[0].GetNoOfHistos()));
4576  dval[0] = fData[0].GetTimeResolution();
4577  prop.Set("Time Resolution", dval[0], "ns");
4578  header->Set("RunInfo/Time Resolution", prop);
4579  header->Set("RunInfo/RedGreen Offsets", fData[0].GetRedGreenOffset());
4580 
4581  // feed DetectorInfo
4582  Int_t histoNo = 0;
4583  PRawRunDataSet *dataSet;
4584  UInt_t size = fData[0].GetNoOfHistos();
4585  for (UInt_t i=0; i<size; i++) {
4586  dataSet = fData[0].GetDataSet(i, false); // i.e. the false means, that i is the index and NOT the histo number
4587  if (dataSet == nullptr) { // something is really wrong
4588  std::cerr << std::endl << ">> PRunDataHandler::WriteMusrRootFile: **ERROR** Couldn't get data set (idx=" << i << ")";
4589  std::cerr << std::endl << ">> something is really wrong!" << std::endl;
4590  return false;
4591  }
4592  histoNo = dataSet->GetHistoNo();
4593  pathName.Form("DetectorInfo/Detector%03d/Name", histoNo);
4594  str = dataSet->GetName();
4595  if (!str.CompareTo("n/a"))
4596  str.Form("Detector%3d", histoNo);
4597  header->Set(pathName, str);
4598  pathName.Form("DetectorInfo/Detector%03d/Histo Number", histoNo);
4599  header->Set(pathName, histoNo);
4600  pathName.Form("DetectorInfo/Detector%03d/Histo Length", histoNo);
4601  header->Set(pathName, static_cast<Int_t>(dataSet->GetData()->size()/fAny2ManyInfo->rebin));
4602  pathName.Form("DetectorInfo/Detector%03d/Time Zero Bin", histoNo);
4603  header->Set(pathName, dataSet->GetTimeZeroBin()/fAny2ManyInfo->rebin);
4604  pathName.Form("DetectorInfo/Detector%03d/First Good Bin", histoNo);
4605  ival = dataSet->GetFirstGoodBin();
4606  header->Set(pathName, static_cast<Int_t>(ival/fAny2ManyInfo->rebin));
4607  pathName.Form("DetectorInfo/Detector%03d/Last Good Bin", histoNo);
4608  ival = dataSet->GetLastGoodBin();
4609  header->Set(pathName, static_cast<Int_t>(ival/fAny2ManyInfo->rebin));
4610  }
4611 
4612  // feed SampleEnvironmentInfo
4613  str = fData[0].GetCryoName()->Copy();
4614  header->Set("SampleEnvironmentInfo/Cryo", str);
4615 
4616  // feed MagneticFieldEnvironmentInfo
4617  str = fData[0].GetMagnetName()->Copy();
4618  header->Set("MagneticFieldEnvironmentInfo/Magnet Name", str);
4619 
4620  // feed BeamlineInfo
4621  str = fData[0].GetBeamline()->Copy();
4622  header->Set("BeamlineInfo/Name", str);
4623 
4624  // feed histos
4625  std::vector<TH1F*> histos;
4626  TH1F *histo = nullptr;
4627  UInt_t length = 0;
4628  if (fAny2ManyInfo->rebin == 1) {
4629  for (UInt_t i=0; i<size; i++) {
4630  dataSet = fData[0].GetDataSet(i, false); // i.e. the false means, that i is the index and NOT the histo number
4631  if (dataSet == nullptr) { // something is really wrong
4632  std::cerr << std::endl << ">> PRunDataHandler::WriteMusrRootFile: **ERROR** Couldn't get data set (idx=" << i << ")";
4633  std::cerr << std::endl << ">> something is really wrong!" << std::endl;
4634  return false;
4635  }
4636  str.Form("hDecay%03d", dataSet->GetHistoNo());
4637  length = dataSet->GetData()->size();
4638  histo = new TH1F(str.Data(), str.Data(), length+1, -0.5, static_cast<Double_t>(length)+0.5);
4639  Int_t sum=0, entries=0;
4640  for (UInt_t j=0; j<length; j++) {
4641  entries = dataSet->GetData()->at(j);
4642  histo->SetBinContent(j+1, entries);
4643  sum += entries;
4644  }
4645  histo->SetEntries(sum);
4646  histos.push_back(histo);
4647  }
4648  } else { // rebin > 1
4649  UInt_t dataRebin = 0;
4650  UInt_t dataCount = 0;
4651  for (UInt_t i=0; i<size; i++) {
4652  dataSet = fData[0].GetDataSet(i, false); // i.e. the false means, that i is the index and NOT the histo number
4653  if (dataSet == nullptr) { // something is really wrong
4654  std::cerr << std::endl << ">> PRunDataHandler::WriteMusrRootFile: **ERROR** Couldn't get data set (idx=" << i << ")";
4655  std::cerr << std::endl << ">> something is really wrong!" << std::endl;
4656  return false;
4657  }
4658  str.Form("hDecay%03d", dataSet->GetHistoNo());
4659  length = dataSet->GetData()->size();
4660  histo = new TH1F(str.Data(), str.Data(), static_cast<Int_t>(length/fAny2ManyInfo->rebin)+1, -0.5, static_cast<Double_t>(static_cast<Int_t>(length/fAny2ManyInfo->rebin))+0.5);
4661  dataCount = 0;
4662  Int_t sum=0, entries=0;
4663  for (UInt_t j=0; j<length; j++) {
4664  if ((j > 0) && (j % fAny2ManyInfo->rebin == 0)) {
4665  dataCount++;
4666  histo->SetBinContent(dataCount, dataRebin);
4667  dataRebin = 0;
4668  }
4669  entries = dataSet->GetData()->at(j);
4670  sum += entries;
4671  dataRebin += static_cast<UInt_t>(entries);
4672  }
4673  histo->SetEntries(sum);
4674  histos.push_back(histo);
4675  }
4676  }
4677 
4678  // add histos to the DecayAnaModule folder
4679  for (UInt_t i=0; i<histos.size(); i++)
4680  decayAnaModule->Add(histos[i]);
4681 
4682  // write file
4683  std::unique_ptr<TFile> fout = std::make_unique<TFile>(fln, "RECREATE", fln);
4684  if (fout == nullptr) {
4685  std::cerr << std::endl << "PRunDataHandler::WriteMusrRootFile(): **ERROR** Couldn't create ROOT file '" << fln << "'" << std::endl;
4686  return false;
4687  }
4688 
4689  fout->cd();
4690  if (header->FillFolder(runHeader))
4691  runHeader->Write();
4692  histosFolder->Write();
4693  fout->Close();
4694 
4695  // cleanup
4696  for (UInt_t i=0; i<histos.size(); i++) {
4697  if (histos[i])
4698  delete histos[i];
4699  }
4700 
4701  // check if root file shall be streamed to stdout
4703  // stream file to stdout
4704  std::ifstream is;
4705  int length=1024;
4706  char *buffer;
4707 
4708  is.open(fln.Data(), std::ios::binary);
4709  if (!is.is_open()) {
4710  std::cerr << std::endl << "PRunDataHandler::WriteMusrRootFile(): **ERROR** Couldn't open the root-file for streaming." << std::endl;
4711  remove(fln.Data());
4712  return false;
4713  }
4714 
4715  // get length of file
4716  is.seekg(0, std::ios::end);
4717  length = is.tellg();
4718  is.seekg(0, std::ios::beg);
4719 
4720  if (length == -1) {
4721  std::cerr << std::endl << "PRunDataHandler::WriteMusrRootFile(): **ERROR** Couldn't determine the root-file size." << std::endl;
4722  remove(fln.Data());
4723  return false;
4724  }
4725 
4726  // allocate memory
4727  buffer = new char [length];
4728 
4729  // read data as a block
4730  while (!is.eof()) {
4731  is.read(buffer, length);
4732  std::cout.write(buffer, length);
4733  }
4734 
4735  is.close();
4736 
4737  delete [] buffer;
4738 
4739  // delete temporary root file
4740  remove(fln.Data());
4741  }
4742 
4743  return true;
4744 }
4745 
4746 //--------------------------------------------------------------------------
4747 // WriteRootFile (private)
4748 //--------------------------------------------------------------------------
4759 {
4760  Bool_t ok = false;
4761  fln = GenerateOutputFileName(fln, ".root", ok);
4762  if (!ok)
4763  return false;
4764 
4765 
4767  std::cout << std::endl << ">> PRunDataHandler::WriteRootFile(): writing a root data file (" << fln.Data() << ") ... " << std::endl;
4768 
4769  // generate data file
4770  TFolder *histosFolder;
4771  TFolder *decayAnaModule;
4772  TFolder *runInfo;
4773 
4774  histosFolder = gROOT->GetRootFolder()->AddFolder("histos", "Histograms");
4775  gROOT->GetListOfBrowsables()->Add(histosFolder, "histos");
4776  decayAnaModule = histosFolder->AddFolder("DecayAnaModule", "muSR decay histograms");
4777 
4778  runInfo = gROOT->GetRootFolder()->AddFolder("RunInfo", "LEM RunInfo");
4779  gROOT->GetListOfBrowsables()->Add(runInfo, "RunInfo");
4780  std::unique_ptr<TLemRunHeader> header = std::make_unique<TLemRunHeader>();
4781  gROOT->GetListOfBrowsables()->Add(runInfo, "RunInfo");
4782 
4783  // feed header info
4784  header->SetRunTitle(fData[0].GetRunTitle()->Data());
4785  header->SetLemSetup(fData[0].GetSetup()->Data());
4786  header->SetRunNumber(fData[0].GetRunNumber());
4787  TString dt = *fData[0].GetStartDate() + "/" + *fData[0].GetStartTime();
4788  header->SetStartTimeString(dt.Data());
4789  dt = *fData[0].GetStopDate() + "/" + *fData[0].GetStopTime();
4790  header->SetStopTimeString(dt.Data());
4791  header->SetStartTime(fData[0].GetStartDateTime());
4792  header->SetStopTime(fData[0].GetStopDateTime());
4793  header->SetModeratorHV(-999.9, 0.0);
4794  header->SetSampleHV(-999.9, 0.0);
4795  header->SetImpEnergy(-999.9);
4796  header->SetSampleTemperature(fData[0].GetTemperature(0), fData[0].GetTempError(0));
4797  header->SetSampleBField(fData[0].GetField(), 0.0);
4798  header->SetTimeResolution(fData[0].GetTimeResolution());
4799  PRawRunDataSet *dataSet = fData[0].GetDataSet(0, false); // i.e. the false means, that i is the index and NOT the histo number
4800  header->SetNChannels(static_cast<UInt_t>(dataSet->GetData()->size()/fAny2ManyInfo->rebin));
4801  header->SetNHist(fData[0].GetNoOfHistos());
4802  header->SetCuts("none");
4803  header->SetModerator("none");
4804 
4805  // feed t0's if possible
4806  UInt_t NoT0s = fData[0].GetNoOfHistos();
4807  if (fData[0].GetNoOfHistos() > NHIST) {
4808  std::cerr << std::endl << ">> PRunDataHandler::WriteRootFile: **WARNING** found more T0's (" << NoT0s << ") than can be handled (" << NHIST << ").";
4809  std::cerr << std::endl << ">> Will only write the first " << NHIST << " T0s!!" << std::endl;
4810  NoT0s = NHIST;
4811  }
4812  Double_t *tt0 = new Double_t[NoT0s];
4813 
4814  for (UInt_t i=0; i<NoT0s; i++) {
4815  dataSet = fData[0].GetDataSet(i, false); // i.e. the false means, that i is the index and NOT the histo number
4816  tt0[i] = dataSet->GetTimeZeroBin()/fAny2ManyInfo->rebin;
4817  }
4818  header->SetTimeZero(tt0);
4819  runInfo->Add(header.get()); // add header to RunInfo folder
4820 
4821  // feed histos
4822  std::vector<TH1F*> histos;
4823  TH1F *histo = nullptr;
4824  Char_t str[32];
4825  UInt_t size = 0;
4826  if (fAny2ManyInfo->rebin == 1) {
4827  for (UInt_t i=0; i<fData[0].GetNoOfHistos(); i++) {
4828  dataSet = fData[0].GetDataSet(i, false); // i.e. the false means, that i is the index and NOT the histo number
4829  if (dataSet == nullptr) { // something is really wrong
4830  std::cerr << std::endl << ">> PRunDataHandler::WriteRootFile: **ERROR** Couldn't get data set (idx=0" << i << ")";
4831  std::cerr << std::endl << ">> something is really wrong!" << std::endl;
4832  return false;
4833  }
4834  size = dataSet->GetData()->size();
4835  snprintf(str, sizeof(str), "hDecay%02d", static_cast<Int_t>(i));
4836  histo = new TH1F(str, str, size+1, -0.5, static_cast<Double_t>(size)+0.5);
4837  for (UInt_t j=0; j<size; j++) {
4838  histo->SetBinContent(j+1, dataSet->GetData()->at(j));
4839  }
4840  histos.push_back(histo);
4841  }
4842  } else { // rebin > 1
4843  UInt_t dataRebin = 0;
4844  UInt_t dataCount = 0;
4845  for (UInt_t i=0; i<fData[0].GetNoOfHistos(); i++) {
4846  dataSet = fData[0].GetDataSet(i, false); // i.e. the false means, that i is the index and NOT the histo number
4847  if (dataSet == nullptr) { // something is really wrong
4848  std::cerr << std::endl << ">> PRunDataHandler::WriteNexusFile: **ERROR** Couldn't get data set (idx=0" << i << ")";
4849  std::cerr << std::endl << ">> something is really wrong!" << std::endl;
4850  return false;
4851  }
4852  size = dataSet->GetData()->size();
4853  snprintf(str, sizeof(str), "hDecay%02d", static_cast<Int_t>(i));
4854  histo = new TH1F(str, str, static_cast<UInt_t>(size/fAny2ManyInfo->rebin)+1, -0.5, static_cast<Double_t>(size)/static_cast<Double_t>(fAny2ManyInfo->rebin)+0.5);
4855  dataCount = 0;
4856  for (UInt_t j=0; j<size; j++) {
4857  if ((j > 0) && (j % fAny2ManyInfo->rebin == 0)) {
4858  dataCount++;
4859  histo->SetBinContent(dataCount, dataRebin);
4860  dataRebin = 0;
4861  }
4862  dataRebin += static_cast<UInt_t>(dataSet->GetData()->at(j));
4863  }
4864  histos.push_back(histo);
4865  }
4866  }
4867 
4868  // add histos to the DecayAnaModule folder
4869  for (UInt_t i=0; i<histos.size(); i++)
4870  decayAnaModule->Add(histos[i]);
4871 
4872  // write file
4873  std::unique_ptr<TFile> fout = std::make_unique<TFile>(fln, "RECREATE", fln);
4874  if (fout == nullptr) {
4875  std::cerr << std::endl << "PRunDataHandler::WriteRootFile(): **ERROR** Couldn't create ROOT file '" << fln << "'" << std::endl;
4876  return false;
4877  }
4878 
4879  fout->cd();
4880  runInfo->Write();
4881  histosFolder->Write();
4882  fout->Close();
4883 
4884  // clean up
4885  for (UInt_t i=0; i<histos.size(); i++) {
4886  delete histos[i];
4887  }
4888  histos.clear();
4889  delete [] tt0;
4890 
4891  // check if root file shall be streamed to stdout
4893  // stream file to stdout
4894  std::ifstream is;
4895  int length=1024;
4896  char *buffer;
4897 
4898  is.open(fln.Data(), std::ios::binary);
4899  if (!is.is_open()) {
4900  std::cerr << std::endl << "PRunDataHandler::WriteRootFile(): **ERROR** Couldn't open the root-file for streaming." << std::endl;
4901  remove(fln.Data());
4902  return false;
4903  }
4904 
4905  // get length of file
4906  is.seekg(0, std::ios::end);
4907  length = is.tellg();
4908  is.seekg(0, std::ios::beg);
4909 
4910  if (length == -1) {
4911  std::cerr << std::endl << "PRunDataHandler::WriteRootFile(): **ERROR** Couldn't determine the root-file size." << std::endl;
4912  remove(fln.Data());
4913  return false;
4914  }
4915 
4916  // allocate memory
4917  buffer = new char [length];
4918 
4919  // read data as a block
4920  while (!is.eof()) {
4921  is.read(buffer, length);
4922  std::cout.write(buffer, length);
4923  }
4924 
4925  is.close();
4926 
4927  delete [] buffer;
4928 
4929  // delete temporary root file
4930  remove(fln.Data());
4931  }
4932 
4933  return true;
4934 }
4935 
4936 //--------------------------------------------------------------------------
4937 // WriteNexusFile (private)
4938 //--------------------------------------------------------------------------
4949 {
4950 #ifdef PNEXUS_ENABLED
4951  Bool_t ok = false;
4952  fln = GenerateOutputFileName(fln, ".nxs", ok);
4953  if (!ok)
4954  return false;
4955 
4957  std::cout << std::endl << ">> PRunDataHandler::WriteNexusFile(): writing a NeXus data file (" << fln.Data() << ") ... " << std::endl;
4958 
4959  // create NeXus object
4960  std::unique_ptr<PNeXus> nxs = std::make_unique<PNeXus>();
4961  if (nxs == nullptr) {
4962  std::cerr << std::endl << ">> PRunDataHandler::WriteNexusFile(): **ERROR** couldn't invoke the NeXus object." << std::endl;
4963  return false;
4964  }
4965 
4966  // set IDF version
4967  nxs->SetIdfVersion(fAny2ManyInfo->idf);
4968 
4969  if (fAny2ManyInfo->idf == 1) {
4970  // fill necessary data structures
4971  nxs->SetFileName(fln.Data());
4972 
4973  // set file creating time
4974  time_t now;
4975  struct tm *tm;
4976  time(&now);
4977  tm = localtime(&now);
4978  std::string str("");
4979  char cstr[128];
4980  strftime(cstr, sizeof(cstr), "%FT%T", tm);
4981  str = std::string(cstr);
4982  nxs->SetFileTime(str);
4983 
4984  nxs->GetEntryIdf1()->SetProgramName("any2many");
4985  nxs->GetEntryIdf1()->SetProgramVersion("$Id$");
4986  nxs->GetEntryIdf1()->SetRunNumber(fData[0].GetRunNumber());
4987  nxs->GetEntryIdf1()->SetTitle(fData[0].GetRunTitle()->Data());
4988  nxs->GetEntryIdf1()->SetNotes("n/a");
4989  nxs->GetEntryIdf1()->SetAnalysis("muonTD");
4990  if (*fData[0].GetLaboratory() != "n/a")
4991  nxs->GetEntryIdf1()->SetLaboratory(fData[0].GetLaboratory()->Data());
4992  if (*fData[0].GetBeamline() != "n/a")
4993  nxs->GetEntryIdf1()->SetBeamline(fData[0].GetBeamline()->Data());
4994  str = std::string(fData[0].GetStartDate()->Data()) + std::string("T") + std::string(fData[0].GetStartTime()->Data());
4995  nxs->GetEntryIdf1()->SetStartTime(str);
4996  str = std::string(fData[0].GetStopDate()->Data()) + std::string("T") + std::string(fData[0].GetStopTime()->Data());
4997  nxs->GetEntryIdf1()->SetStopTime(str);
4998  nxs->GetEntryIdf1()->SetSwitchingState(1);
4999  nxs->GetEntryIdf1()->GetUser()->SetName("n/a");
5000  nxs->GetEntryIdf1()->GetUser()->SetExperimentNumber("n/a");
5001  nxs->GetEntryIdf1()->GetSample()->SetName(fData[0].GetSample()->Data());
5002  nxs->GetEntryIdf1()->GetSample()->SetPhysProp("temperature", fData[0].GetTemperature(0), "Kelvin");
5003  nxs->GetEntryIdf1()->GetSample()->SetPhysProp("magnetic_field", fData[0].GetField(), "Gauss");
5004  nxs->GetEntryIdf1()->GetSample()->SetEnvironment(fData[0].GetSetup()->Data());
5005  nxs->GetEntryIdf1()->GetSample()->SetShape("n/a");
5006  nxs->GetEntryIdf1()->GetSample()->SetMagneticFieldVectorAvailable(0);
5007  if (*fData[0].GetInstrument() != "n/a")
5008  nxs->GetEntryIdf1()->GetInstrument()->SetName(fData[0].GetInstrument()->Data());
5009  nxs->GetEntryIdf1()->GetInstrument()->GetDetector()->SetNumber(fData[0].GetNoOfHistos());
5010  nxs->GetEntryIdf1()->GetInstrument()->GetCollimator()->SetType("n/a");
5011  // calculate the total number of counts
5012  double total_counts = 0;
5013  PRawRunDataSet *dataSet = nullptr;
5014  for (unsigned int i=0; i<fData[0].GetNoOfHistos(); i++) {
5015  dataSet = fData[0].GetDataSet(i, false); // i.e. the false means, that i is the index and NOT the histo number
5016  if (dataSet == nullptr) { // something is really wrong
5017  std::cerr << std::endl << ">> PRunDataHandler::WriteNexusFile: **ERROR** Couldn't get data set (idx=0" << i << ")";
5018  std::cerr << std::endl << ">> something is really wrong!" << std::endl;
5019  return false;
5020  }
5021  for (unsigned int j=0; j<dataSet->GetData()->size(); j++)
5022  total_counts += dataSet->GetData()->at(j);
5023  }
5024  double total_counts_mev = (double) total_counts / 1.0e6;
5025  nxs->GetEntryIdf1()->GetInstrument()->GetBeam()->SetTotalCounts(total_counts_mev);
5026  nxs->GetEntryIdf1()->GetInstrument()->GetBeam()->SetUnits("Mev");
5027 
5028  nxs->GetEntryIdf1()->GetData()->SetTimeResolution(fData[0].GetTimeResolution()*fAny2ManyInfo->rebin, "ns");
5029 
5030  for (unsigned int i=0; i<fData[0].GetNoOfHistos(); i++) {
5031  dataSet = fData[0].GetDataSet(i, false); // i.e. the false means, that i is the index and NOT the histo number
5032  nxs->GetEntryIdf1()->GetData()->SetT0(static_cast<Int_t>(dataSet->GetTimeZeroBin()/fAny2ManyInfo->rebin), i);
5033  nxs->GetEntryIdf1()->GetData()->SetFirstGoodBin(static_cast<Int_t>(dataSet->GetFirstGoodBin()/fAny2ManyInfo->rebin), i);
5034  nxs->GetEntryIdf1()->GetData()->SetLastGoodBin(static_cast<Int_t>(dataSet->GetLastGoodBin()/fAny2ManyInfo->rebin), i);
5035  }
5036 
5037  // feed histos
5038  PUIntVector data;
5039  UInt_t size = 0;
5040  if (fAny2ManyInfo->rebin == 1) {
5041  for (UInt_t i=0; i<fData[0].GetNoOfHistos(); i++) {
5042  dataSet = fData[0].GetDataSet(i, false); // i.e. the false means, that i is the index and NOT the histo number
5043  if (dataSet == nullptr) { // something is really wrong
5044  std::cerr << std::endl << ">> PRunDataHandler::WriteNexusFile: **ERROR** Couldn't get data set (idx=" << i << ")";
5045  std::cerr << std::endl << ">> something is really wrong!" << std::endl;
5046  return false;
5047  }
5048  size = dataSet->GetData()->size();
5049  for (UInt_t j=0; j<size; j++) {
5050  data.push_back((UInt_t)dataSet->GetData()->at(j));
5051  }
5052  nxs->GetEntryIdf1()->GetData()->SetHisto(data, i);
5053  data.clear();
5054  }
5055  } else { // rebin > 1
5056  UInt_t dataRebin = 0;
5057  UInt_t dataCount = 0;
5058  for (UInt_t i=0; i<fData[0].GetNoOfHistos(); i++) {
5059  dataSet = fData[0].GetDataSet(i, false); // i.e. the false means, that i is the index and NOT the histo number
5060  if (dataSet == nullptr) { // something is really wrong
5061  std::cerr << std::endl << ">> PRunDataHandler::WriteNexusFile: **ERROR** Couldn't get data set (idx=" << i << ")";
5062  std::cerr << std::endl << ">> something is really wrong!" << std::endl;
5063  return false;
5064  }
5065  size = dataSet->GetData()->size();
5066  dataCount = 0;
5067  for (UInt_t j=0; j<size; j++) {
5068  if ((j > 0) && (j % fAny2ManyInfo->rebin == 0)) {
5069  dataCount++;
5070  data.push_back(dataRebin);
5071  dataRebin = 0;
5072  }
5073  dataRebin += static_cast<UInt_t>(dataSet->GetData()->at(j));
5074  }
5075  nxs->GetEntryIdf1()->GetData()->SetHisto(data, i);
5076  data.clear();
5077  }
5078  }
5079  } else if (fAny2ManyInfo->idf == 2) {
5080  // fill necessary data structures
5081  nxs->SetFileName(fln.Data());
5082 
5083  // set file creating time
5084  time_t now;
5085  struct tm *tm;
5086  time(&now);
5087  tm = localtime(&now);
5088  std::string str("");
5089  char cstr[128];
5090  strftime(cstr, sizeof(cstr), "%FT%T", tm);
5091  str = std::string(cstr);
5092  nxs->SetFileTime(str);
5093 
5094  // NXroot info
5095  nxs->SetCreator("PSI: any2many");
5096 
5097  // NXentry info
5098  nxs->GetEntryIdf2()->SetDefinition("muonTD");
5099  nxs->GetEntryIdf2()->SetProgramName("any2many");
5100  nxs->GetEntryIdf2()->SetProgramVersion("$Id$");
5101  nxs->GetEntryIdf2()->SetRunNumber(fData[0].GetRunNumber());
5102  nxs->GetEntryIdf2()->SetTitle(fData[0].GetRunTitle()->Data());
5103  str = std::string(fData[0].GetStartDate()->Data()) + std::string("T") + std::string(fData[0].GetStartTime()->Data());
5104  nxs->GetEntryIdf2()->SetStartTime(str);
5105  str = std::string(fData[0].GetStopDate()->Data()) + std::string("T") + std::string(fData[0].GetStopTime()->Data());
5106  nxs->GetEntryIdf2()->SetStopTime(str);
5107 
5108  nxs->GetEntryIdf2()->SetExperimentIdentifier("n/a");
5109 
5110  // NXuser info
5111  nxs->GetEntryIdf2()->GetUser()->SetName("n/a");
5112 
5113  // NXsample info
5114  nxs->GetEntryIdf2()->GetSample()->SetName(fData[0].GetSample()->Data());
5115  nxs->GetEntryIdf2()->GetSample()->SetDescription("n/a");
5116  nxs->GetEntryIdf2()->GetSample()->SetPhysProp("temperature_1", fData[0].GetTemperature(0), "Kelvin");
5117  nxs->GetEntryIdf2()->GetSample()->SetPhysProp("magnetic_field_1", fData[0].GetField(), "Gauss");
5118  nxs->GetEntryIdf2()->GetSample()->SetEnvironmentTemp(fData[0].GetSetup()->Data());
5119  nxs->GetEntryIdf2()->GetSample()->SetEnvironmentField("n/a");
5120 
5121  // here would be the information for NXinstrument. Currently there are not much information to feed this
5122  nxs->GetEntryIdf2()->GetInstrument()->SetName(fData[0].GetInstrument()->Data());
5123 
5124  // NXinstrument/NXsource
5125  nxs->GetEntryIdf2()->GetInstrument()->GetSource()->SetName(fData[0].GetLaboratory()->Data());
5126  nxs->GetEntryIdf2()->GetInstrument()->GetSource()->SetType(fData[0].GetMuonSource()->Data());
5127  nxs->GetEntryIdf2()->GetInstrument()->GetSource()->SetProbe(fData[0].GetMuonSpecies()->Data());
5128 
5129  // NXinstrument/NXbeamline
5130  nxs->GetEntryIdf2()->GetInstrument()->GetBeamline()->SetName(fData[0].GetBeamline()->Data());
5131 
5132  // NXinstrument/NXdetector
5133  nxs->GetEntryIdf2()->GetInstrument()->GetDetector()->SetDescription(fData[0].GetInstrument()->Data()); // assume that this should be the instrument name
5134  nxs->GetEntryIdf2()->GetInstrument()->GetDetector()->SetNoOfPeriods(0); // currently red/green is not distinguished
5135  nxs->GetEntryIdf2()->GetInstrument()->GetDetector()->SetNoOfSpectra(fData[0].GetNoOfHistos());
5136  PRawRunDataSet *dataSet = fData[0].GetDataSet(0, false); // i.e. the false means, that i is the index and NOT the histo number
5137  if (dataSet == nullptr) { // something is really wrong
5138  std::cerr << std::endl << ">> PRunDataHandler::WriteNeXusFile: **ERROR** Couldn't get data set (idx=0)";
5139  std::cerr << std::endl << ">> something is really wrong!" << std::endl;
5140  return false;
5141  }
5142  nxs->GetEntryIdf2()->GetInstrument()->GetDetector()->SetNoOfBins((unsigned int)(dataSet->GetData()->size() / fAny2ManyInfo->rebin));
5143  nxs->GetEntryIdf2()->GetInstrument()->GetDetector()->SetTimeResolution(fData[0].GetTimeResolution()*fAny2ManyInfo->rebin, "ns");
5144  int *histo = nullptr;
5145  int idx = 0;
5146  if (fAny2ManyInfo->rebin == 1) {
5147  histo = new int[fData[0].GetNoOfHistos()*dataSet->GetData()->size()];
5148  idx = 0;
5149  for (int i=0; i<nxs->GetEntryIdf2()->GetInstrument()->GetDetector()->GetNoOfSpectra(); i++) {
5150  dataSet = fData[0].GetDataSet(i, false); // i.e. the false means, that i is the index and NOT the histo number
5151  if (dataSet == nullptr) { // something is really wrong
5152  std::cerr << std::endl << ">> PRunDataHandler::WriteNeXusFile: **ERROR** Couldn't get data set (idx=" << i << ")";
5153  std::cerr << std::endl << ">> something is really wrong!" << std::endl;
5154  return false;
5155  }
5156  for (unsigned int j=0; j<dataSet->GetData()->size(); j++) {
5157  *(histo+idx++) = (int) dataSet->GetData()->at(j);
5158  }
5159  }
5160  nxs->GetEntryIdf2()->GetInstrument()->GetDetector()->SetHistos(histo);
5161  // clean up
5162  if (histo) {
5163  delete [] histo;
5164  histo = nullptr;
5165  }
5166  } else { // rebin > 1
5167  histo = new int[fData[0].GetNoOfHistos()*(int)(dataSet->GetData()->size()/fAny2ManyInfo->rebin)];
5168  int counts = 0;
5169  idx = 0;
5170  for (int i=0; i<nxs->GetEntryIdf2()->GetInstrument()->GetDetector()->GetNoOfSpectra(); i++) {
5171  dataSet = fData[0].GetDataSet(i, false); // i.e. the false means, that i is the index and NOT the histo number
5172  if (dataSet == nullptr) { // something is really wrong
5173  std::cerr << std::endl << ">> PRunDataHandler::WriteNeXusFile: **ERROR** Couldn't get data set (idx=" << i << ")";
5174  std::cerr << std::endl << ">> something is really wrong!" << std::endl;
5175  return false;
5176  }
5177  for (unsigned int j=0; j<dataSet->GetData()->size(); j++) {
5178  if ((j>0) && (j % fAny2ManyInfo->rebin == 0)) {
5179  *(histo+idx++) = counts;
5180  counts = 0;
5181  }
5182  counts += (int) dataSet->GetData()->at(j);
5183  }
5184  }
5185  nxs->GetEntryIdf2()->GetInstrument()->GetDetector()->SetHistos(histo);
5186  // clean up
5187  if (histo) {
5188  delete [] histo;
5189  histo = nullptr;
5190  }
5191  }
5192 
5193  // handle spectrum index
5194  for (int i=0; i<nxs->GetEntryIdf2()->GetInstrument()->GetDetector()->GetNoOfSpectra(); i++)
5195  nxs->GetEntryIdf2()->GetInstrument()->GetDetector()->SetSpectrumIndex(i+1);
5196 
5197  // handle histogram resolution
5198  nxs->GetEntryIdf2()->GetInstrument()->GetDetector()->SetTimeResolution(fData[0].GetTimeResolution()*fAny2ManyInfo->rebin, "ns");
5199 
5200  // handle raw time
5201  std::vector<double> raw_time;
5202  UInt_t size = (unsigned int)(dataSet->GetData()->size() / fAny2ManyInfo->rebin);
5203  for (unsigned int i=0; i<size; i++) {
5204  raw_time.push_back((double)i * fData[0].GetTimeResolution() * fAny2ManyInfo->rebin * 1.0e-3); // since time resolution is given in ns, the factor 1.0e-3 is needed to convert to us
5205  }
5206  nxs->GetEntryIdf2()->GetInstrument()->GetDetector()->SetRawTime(raw_time);
5207  nxs->GetEntryIdf2()->GetInstrument()->GetDetector()->SetRawTimeUnit("micro.second");
5208  nxs->GetEntryIdf2()->GetInstrument()->GetDetector()->SetRawTimeName("time");
5209  raw_time.clear();
5210 
5211  // handle t0
5212  nxs->GetEntryIdf2()->GetInstrument()->GetDetector()->SetT0Tag(2); // i.e. t0[#histo] format
5213  int *t0 = new int[fData[0].GetNoOfHistos()];
5214  int *fgb = new int[fData[0].GetNoOfHistos()];
5215  int *lgb = new int[fData[0].GetNoOfHistos()];
5216  if ((t0==0) || (fgb==0) || (lgb==0)) {
5217  std::cerr << std::endl << ">> PRunDataHandler::WriteNeXusFile: **ERROR** Couldn't allocate memory for t0, fgb, lgb" << std::endl;
5218  return false;
5219  }
5220  for (unsigned int i=0; i<fData[0].GetNoOfHistos(); i++) {
5221  PRawRunDataSet *dataSet = fData[0].GetDataSet(i, false); // i.e. the false means, that i is the index and NOT the histo number
5222  if (dataSet == nullptr) { // something is really wrong
5223  std::cerr << std::endl << ">> PRunDataHandler::WriteNeXusFile: **ERROR** Couldn't get data set (idx=0)";
5224  std::cerr << std::endl << ">> something is really wrong!" << std::endl;
5225  return false;
5226  }
5227  t0[i] = (int)(dataSet->GetTimeZeroBin() / fAny2ManyInfo->rebin);
5228  fgb[i] = (int)(dataSet->GetFirstGoodBin() / fAny2ManyInfo->rebin);
5229  lgb[i] = (int)(dataSet->GetLastGoodBin() / fAny2ManyInfo->rebin);
5230  }
5231  nxs->GetEntryIdf2()->GetInstrument()->GetDetector()->SetT0(t0);
5232  nxs->GetEntryIdf2()->GetInstrument()->GetDetector()->SetFirstGoodBin(fgb);
5233  nxs->GetEntryIdf2()->GetInstrument()->GetDetector()->SetLastGoodBin(lgb);
5234 
5235  // clean up
5236  if (t0) delete [] t0;
5237  if (fgb) delete [] fgb;
5238  if (lgb) delete [] lgb;
5239  } else {
5240  return false;
5241  }
5242 
5243  // filter out the proper file type, i.e. HDF4, HDF5, or XML
5244  char fileType[32];
5245  memset(fileType, '\0', 32);
5246  if (!fAny2ManyInfo->outFormat.CompareTo("nexus1-hdf4", TString::kIgnoreCase) || !fAny2ManyInfo->outFormat.CompareTo("nexus2-hdf4", TString::kIgnoreCase))
5247  strncpy(fileType, "hdf4", sizeof(fileType));
5248  else if (!fAny2ManyInfo->outFormat.CompareTo("nexus1-hdf5", TString::kIgnoreCase) || !fAny2ManyInfo->outFormat.CompareTo("nexus2-hdf5", TString::kIgnoreCase))
5249  strncpy(fileType, "hdf5", sizeof(fileType));
5250  else if (!fAny2ManyInfo->outFormat.CompareTo("nexus1-xml", TString::kIgnoreCase) || !fAny2ManyInfo->outFormat.CompareTo("nexus2-xml", TString::kIgnoreCase))
5251  strncpy(fileType, "xml", sizeof(fileType));
5252  else {
5253  std::cerr << std::endl << ">> PRunDataHandler::WriteNexusFile(): **ERROR** undefined output NeXus format " << fAny2ManyInfo->outFormat.Data() << " found.";
5254  std::cerr << std::endl << ">> Allowed are: hdf4, hdf5, xml" << std::endl;
5255  return false;
5256  }
5257 
5258  // write file
5259  nxs->WriteFile(fln, fileType, fAny2ManyInfo->idf);
5260 #else
5261  std::cout << std::endl << ">> PRunDataHandler::WriteNexusFile(): Sorry, not enabled at configuration level, i.e. --enable-NeXus when executing configure" << std::endl << std::endl;
5262 #endif
5263 
5264  return true;
5265 }
5266 
5267 //--------------------------------------------------------------------------
5268 // WriteWkmFile (private)
5269 //--------------------------------------------------------------------------
5280 {
5281  // check if a LEM nemu file needs to be written
5282  bool lem_wkm_style = false;
5283  if (!fData[0].GetBeamline()->CompareTo("mue4", TString::kIgnoreCase))
5284  lem_wkm_style = true;
5285 
5286  Bool_t ok = false;
5287  if (lem_wkm_style)
5288  fln = GenerateOutputFileName(fln, ".nemu", ok);
5289  else
5290  fln = GenerateOutputFileName(fln, ".wkm", ok);
5291 
5292  if (!ok)
5293  return false;
5294 
5295  TString fileName = fln;
5296 
5298  std::cout << std::endl << ">> PRunDataHandler::WriteWkmFile(): writing a wkm data file (" << fln.Data() << ") ... " << std::endl;
5299 
5300  // write ascii file
5301  std::ofstream fout;
5302  std::streambuf* strm_buffer = nullptr;
5303 
5305  // open data-file
5306  fout.open(fln.Data(), std::ofstream::out);
5307  if (!fout.is_open()) {
5308  std::cerr << std::endl << ">> PRunDataHandler::WriteWkmFile **ERROR** Couldn't open data file (" << fln.Data() << ") for writing, sorry ...";
5309  std::cerr << std::endl;
5310  return false;
5311  }
5312 
5313  // save output buffer of the stream
5314  strm_buffer = std::cout.rdbuf();
5315 
5316  // redirect output into the file
5317  std::cout.rdbuf(fout.rdbuf());
5318  }
5319 
5320  // write header
5321  std::cout << "- WKM data file converted with any2many";
5322  if (lem_wkm_style) {
5323  std::cout << std::endl << "NEMU_Run: " << fData[0].GetRunNumber();
5324  std::cout << std::endl << "nemu_Run: " << fileName.Data();
5325  } else {
5326  std::cout << std::endl << "Run: " << fData[0].GetRunNumber();
5327  }
5328  std::cout << std::endl << "Date: " << fData[0].GetStartTime()->Data() << " " << fData[0].GetStartDate()->Data() << " / " << fData[0].GetStopTime()->Data() << " " << fData[0].GetStopDate()->Data();
5329  if (fData[0].GetRunTitle()->Length() > 0)
5330  std::cout << std::endl << "Title: " << fData[0].GetRunTitle()->Data();
5331  if (fData[0].GetField() != PMUSR_UNDEFINED) {
5332  std::cout << std::endl << "Field: " << fData[0].GetField();
5333  } else {
5334  std::cout << std::endl << "Field: ??";
5335  }
5336  std::cout << std::endl << "Setup: " << fData[0].GetSetup()->Data();
5337  if (fData[0].GetNoOfTemperatures() == 1) {
5338  std::cout << std::endl << "Temp: " << fData[0].GetTemperature(0);
5339  } else if (fData[0].GetNoOfTemperatures() > 1) {
5340  std::cout << std::endl << "Temp(meas1): " << fData[0].GetTemperature(0) << " +- " << fData[0].GetTempError(0);
5341  std::cout << std::endl << "Temp(meas2): " << fData[0].GetTemperature(1) << " +- " << fData[0].GetTempError(1);
5342  } else {
5343  std::cout << std::endl << "Temp: ??";
5344  }
5345  if (lem_wkm_style)
5346  std::cout << std::endl << "TOF(M3S1): nocut";
5347  std::cout << std::endl << "Groups: " << fData[0].GetNoOfHistos();
5348  UInt_t histo0 = 1;
5349  if (fAny2ManyInfo->groupHistoList.size() != 0) { // red/green list found
5350  histo0 = fAny2ManyInfo->groupHistoList[0]+1; // take the first available red/green entry
5351  }
5352  std::cout << std::endl << "Channels: " << static_cast<UInt_t>(fData[0].GetDataBin(histo0)->size()/fAny2ManyInfo->rebin);
5353  std::cout.precision(10);
5354  std::cout << std::endl << "Resolution: " << fData[0].GetTimeResolution()*fAny2ManyInfo->rebin/1.0e3; // ns->us
5355  std::cout.setf(std::ios::fixed,std::ios::floatfield); // floatfield set to fixed
5356 
5357  // write data
5358  UInt_t no_of_bins_per_line = 16;
5359  if (lem_wkm_style)
5360  no_of_bins_per_line = 10;
5361 
5362  PRawRunDataSet *dataSet;
5363 
5364  if (fAny2ManyInfo->rebin == 1) {
5365  for (UInt_t i=0; i<fData[0].GetNoOfHistos(); i++) {
5366  std::cout << std::endl << std::endl;
5367  dataSet = fData[0].GetDataSet(i, false); // i.e. the false means, that i is the index and NOT the histo number
5368  for (UInt_t j=0; j<dataSet->GetData()->size(); j++) {
5369  if ((j > 0) && (j % no_of_bins_per_line == 0))
5370  std::cout << std::endl;
5371  if (lem_wkm_style)
5372  std::cout << std::setw(8) << static_cast<Int_t>(dataSet->GetData()->at(j));
5373  else
5374  std::cout << static_cast<Int_t>(dataSet->GetData()->at(j)) << " ";
5375  }
5376  }
5377  } else { // rebin > 1
5378  Int_t dataRebin = 0;
5379  UInt_t count = 0;
5380  for (UInt_t i=0; i<fData[0].GetNoOfHistos(); i++) {
5381  std::cout << std::endl << std::endl;
5382  count = 0;
5383  dataRebin = 0; // reset rebin
5384  dataSet = fData[0].GetDataSet(i, false); // i.e. the false means, that i is the index and NOT the histo number
5385  for (UInt_t j=0; j<dataSet->GetData()->size(); j++) {
5386  if ((j > 0) && (j % fAny2ManyInfo->rebin == 0)) {
5387  if (lem_wkm_style)
5388  std::cout << std::setw(8) << dataRebin;
5389  else
5390  std::cout << dataRebin << " ";
5391  count++;
5392  dataRebin = 0;
5393  if ((count > 0) && (count % no_of_bins_per_line == 0) && (j != dataSet->GetData()->size()-1))
5394  std::cout << std::endl;
5395  } else {
5396  dataRebin += static_cast<Int_t>(dataSet->GetData()->at(j));
5397  }
5398  }
5399  }
5400  }
5401 
5403  // restore old output buffer
5404  std::cout.rdbuf(strm_buffer);
5405 
5406  fout.close();
5407  }
5408 
5409  return true;
5410 }
5411 
5412 //--------------------------------------------------------------------------
5413 // WritePsiBinFile (private)
5414 //--------------------------------------------------------------------------
5425 {
5426  Bool_t ok = false;
5427  fln = GenerateOutputFileName(fln, ".bin", ok);
5428  if (!ok)
5429  return false;
5430 
5432  std::cout << std::endl << ">> PRunDataHandler::WritePsiBinFile(): writing a psi-bin data file (" << fln.Data() << ") ... " << std::endl;
5433 
5434  MuSR_td_PSI_bin psibin;
5435  int status = 0;
5436 
5437  // fill header information
5438  // run number
5439  psibin.PutRunNumberInt(fData[0].GetRunNumber());
5440  // length of histograms
5441  UInt_t histo0 = 1;
5442  if (fAny2ManyInfo->groupHistoList.size() != 0) { // red/green list found
5443  histo0 = fAny2ManyInfo->groupHistoList[0]+1; // take the first available red/green entry
5444  }
5445  psibin.PutHistoLengthBin(static_cast<int>(fData[0].GetDataBin(histo0)->size()/fAny2ManyInfo->rebin));
5446  // number of histograms
5447  psibin.PutNumberHistoInt(static_cast<int>(fData[0].GetNoOfHistos()));
5448  // run title = sample (10 char) / temp (10 char) / field (10 char) / orientation (10 char)
5449  char cstr[11];
5450  // sample
5451  if (fData[0].GetSample()->Length() > 0)
5452  strncpy(cstr, fData[0].GetSample()->Data(), 10);
5453  else
5454  strcpy(cstr, "??");
5455  cstr[10] = '\0';
5456  psibin.PutSample(cstr);
5457  // temp
5458  if (fData[0].GetNoOfTemperatures() > 0)
5459  snprintf(cstr, 10, "%.1f K", fData[0].GetTemperature(0));
5460  else
5461  strcpy(cstr, "?? K");
5462  cstr[10] = '\0';
5463  psibin.PutTemp(cstr);
5464  // field
5465  if (fData[0].GetField() > 0)
5466  snprintf(cstr, 10, "%.1f G", fData[0].GetField());
5467  else
5468  strcpy(cstr, "?? G");
5469  cstr[10] = '\0';
5470  psibin.PutField(cstr);
5471  // orientation
5472  if (fData[0].GetOrientation()->Length() > 0)
5473  strncpy(cstr, fData[0].GetOrientation()->Data(), 10);
5474  else
5475  strcpy(cstr, "??");
5476  cstr[10] = '\0';
5477  psibin.PutOrient(cstr);
5478  // setup
5479  if (fData[0].GetSetup()->Length() > 0)
5480  strncpy(cstr, fData[0].GetSetup()->Data(), 10);
5481  else
5482  strcpy(cstr, "??");
5483  cstr[10] = '\0';
5484  psibin.PutSetup(cstr);
5485 
5486  // handle PSI-BIN start/stop Time/Date. PSI-BIN requires: Time -> HH:MM:SS, and Date -> DD-MMM-YY
5487  // internally given: Time -> HH:MM:SS, and Date -> YYYY-MM-DD
5488  // run start date
5489  std::vector<std::string> svec;
5490  TString str, date;
5491  TDatime dt;
5492  int year, month, day;
5493 
5494  // 28-Aug-2014, TP: the following line does not work, it generates the wrong date
5495  //dt.Set(fData[0].GetStartDateTime()); //as35
5496  // the following generates the correct date entry
5497  date.Append(*fData[0].GetStartDate());
5498  sscanf((const char*)date.Data(),"%04d-%02d-%02d", &year, &month, &day);
5499  dt.Set(year, month, day, 0, 0, 0);
5500 
5501  date.Form("%02d-", dt.GetDay());
5502  date.Append(GetMonth(dt.GetMonth()));
5503  date.Append("-");
5504  date.Append(GetYear(dt.GetYear()));
5505  strncpy(cstr, date.Data(), 9);
5506  cstr[9] = '\0';
5507  svec.push_back(cstr);
5508  // run start time
5509  strncpy(cstr, fData[0].GetStartTime()->Data(), 8);
5510  cstr[8] = '\0';
5511  svec.push_back(cstr);
5512  psibin.PutTimeStartVector(svec);
5513  svec.clear();
5514 
5515  // run stop date
5516  // 28-Aug-2014, TP: the following line does not work, it generates the wrong date
5517  //dt.Set(fData[0].GetStopDateTime());
5518  // the following generates the correct date entry
5519  date.Clear();
5520  date.Append(*fData[0].GetStopDate());
5521  sscanf((const char*)date.Data(),"%04d-%02d-%02d", &year, &month, &day);
5522  dt.Set(year, month, day, 0, 0, 0);
5523 
5524  date.Form("%02d-", dt.GetDay());
5525  date.Append(GetMonth(dt.GetMonth()));
5526  date.Append("-");
5527  date.Append(GetYear(dt.GetYear()));
5528  strncpy(cstr, date.Data(), 9);
5529  cstr[9] = '\0';
5530  svec.push_back(cstr);
5531  // run stop time
5532  strncpy(cstr, fData[0].GetStopTime()->Data(), 8);
5533  cstr[8] = '\0';
5534  svec.push_back(cstr);
5535  psibin.PutTimeStopVector(svec);
5536  svec.clear();
5537 
5538  // number of measured temperatures
5539  psibin.PutNumberTemperatureInt(fData[0].GetNoOfTemperatures());
5540 
5541  // mean temperatures
5542  std::vector<double> dvec;
5543  for (UInt_t i=0; i<fData[0].GetNoOfTemperatures(); i++)
5544  dvec.push_back(fData[0].GetTemperature(i));
5545  psibin.PutTemperaturesVector(dvec);
5546 
5547  // standard deviation of temperatures
5548  dvec.clear();
5549  for (UInt_t i=0; i<fData[0].GetNoOfTemperatures(); i++)
5550  dvec.push_back(fData[0].GetTempError(i));
5551  psibin.PutDevTemperaturesVector(dvec);
5552 
5553  // write comment
5554  psibin.PutComment(fData[0].GetRunTitle()->Data());
5555 
5556  // write time resolution
5557  psibin.PutBinWidthNanoSec(fData[0].GetTimeResolution()*fAny2ManyInfo->rebin);
5558 
5559  // write scaler dummies
5560  psibin.PutNumberScalerInt(0);
5561 
5562  // feed detector related info like, histogram names, t0, fgb, lgb
5563  Int_t ival = 0;
5564  PRawRunDataSet *dataSet;
5565  UInt_t size = fData[0].GetNoOfHistos();
5566  for (UInt_t i=0; i<size; i++) {
5567  dataSet = fData[0].GetDataSet(i, false); // i.e. the false means, that i is the index and NOT the histo number
5568  if (dataSet == nullptr) { // something is really wrong
5569  std::cerr << std::endl << ">> PRunDataHandler::WritePsiBinFile: **ERROR** Couldn't get data set (idx=" << i << ")";
5570  std::cerr << std::endl << ">> something is really wrong!" << std::endl;
5571  return false;
5572  }
5573 
5574  // detector name
5575  str = dataSet->GetName();
5576  if (!str.CompareTo("n/a"))
5577  str.Form("Detector%3d", i+1);
5578  psibin.PutNameHisto(str.Data(), i);
5579  // time zero bin
5580  ival = static_cast<Int_t>(dataSet->GetTimeZeroBin()/fAny2ManyInfo->rebin);
5581  psibin.PutT0Int(i, ival);
5582  // first good bin
5583  ival = static_cast<Int_t>(dataSet->GetFirstGoodBin()/fAny2ManyInfo->rebin);
5584  psibin.PutFirstGoodInt(i, ival);
5585  // last good bin
5586  ival = static_cast<Int_t>(dataSet->GetLastGoodBin()/fAny2ManyInfo->rebin);
5587  psibin.PutLastGoodInt(i, ival);
5588  }
5589 
5590  // feed histos
5591  std::vector< std::vector<int> > histos;
5592  histos.resize(fData[0].GetNoOfHistos());
5593  UInt_t length = 0;
5594  if (fAny2ManyInfo->rebin == 1) {
5595  for (UInt_t i=0; i<size; i++) {
5596  dataSet = fData[0].GetDataSet(i, false); // i.e. the false means, that i is the index and NOT the histo number
5597  if (dataSet == nullptr) { // something is really wrong
5598  std::cerr << std::endl << ">> PRunDataHandler::WritePsiBinFile: **ERROR** Couldn't get data set (idx=" << i << ")";
5599  std::cerr << std::endl << ">> something is really wrong!" << std::endl;
5600  return false;
5601  }
5602  length = dataSet->GetData()->size();
5603  histos[i].resize(length);
5604  for (UInt_t j=0; j<length; j++) {
5605  histos[i][j] = static_cast<Int_t>(dataSet->GetData()->at(j));
5606  }
5607  }
5608  } else { // rebin > 1
5609  UInt_t dataRebin = 0;
5610  UInt_t dataCount = 0;
5611  for (UInt_t i=0; i<size; i++) {
5612  dataSet = fData[0].GetDataSet(i, false); // i.e. the false means, that i is the index and NOT the histo number
5613  if (dataSet == nullptr) { // something is really wrong
5614  std::cerr << std::endl << ">> PRunDataHandler::WritePsiBinFile: **ERROR** Couldn't get data set (idx=" << i << ")";
5615  std::cerr << std::endl << ">> something is really wrong!" << std::endl;
5616  return false;
5617  }
5618  length = dataSet->GetData()->size();
5619  dataCount = 0;
5620  for (UInt_t j=0; j<length; j++) {
5621  if ((j > 0) && (j % fAny2ManyInfo->rebin == 0)) {
5622  dataCount++;
5623  histos[i].push_back(dataRebin);
5624  dataRebin = 0;
5625  }
5626  dataRebin += static_cast<UInt_t>(dataSet->GetData()->at(j));
5627  }
5628  }
5629  }
5630  status = psibin.PutHistoArrayInt(histos, 2); // tag 2 means: lift histo length restriction on only make sure it is < 32512
5631  if (status != 0) {
5632  std::cerr << std::endl << ">> PRunDataHandler::WritePsiBinFile(): " << psibin.ConsistencyStatus() << std::endl;
5633  return false;
5634  }
5635 
5636  if (!psibin.CheckDataConsistency(2)) {
5637  std::cerr << std::endl << ">> PRunDataHandler::WritePsiBinFile(): " << psibin.ConsistencyStatus() << std::endl;
5638  return false;
5639  }
5640 
5641  // write data to file
5642  status = psibin.Write(fln.Data());
5643 
5644  if (status != 0) {
5645  std::cerr << std::endl << ">> PRunDataHandler::WritePsiBinFile(): " << psibin.WriteStatus() << std::endl;
5646  return false;
5647  }
5648 
5649  return true;
5650 }
5651 
5652 //--------------------------------------------------------------------------
5653 // WriteMudFile (private)
5654 //--------------------------------------------------------------------------
5665 {
5666  Bool_t ok = false;
5667  fln = GenerateOutputFileName(fln, ".msr", ok);
5668  if (!ok)
5669  return false;
5670 
5672  std::cout << std::endl << ">> PRunDataHandler::WriteMudFile(): writing a mud data file (" << fln.Data() << ") ... " << std::endl;
5673 
5674  // generate the mud data file
5675  int fd = MUD_openWrite((char*)fln.Data(), MUD_FMT_TRI_TD_ID);
5676  if (fd == -1) {
5677  std::cerr << std::endl << ">> PRunDataHandler::WriteMudFile(): **ERROR** couldn't open mud data file for write ..." << std::endl;
5678  return false;
5679  }
5680 
5681  // generate header information
5682  char dummy[32], info[128];
5683  strcpy(dummy, "???");
5684  MUD_setRunDesc(fd, MUD_SEC_GEN_RUN_DESC_ID);
5685  MUD_setExptNumber(fd, 0);
5686  MUD_setRunNumber(fd, fData[0].GetRunNumber());
5687  Int_t ival = fData[0].GetStopDateTime()-fData[0].GetStartDateTime();
5688  MUD_setElapsedSec(fd, ival);
5689  MUD_setTimeBegin(fd, fData[0].GetStartDateTime());
5690  MUD_setTimeEnd(fd, fData[0].GetStopDateTime());
5691  MUD_setTitle(fd, (char *)fData[0].GetRunTitle()->Data());
5692  MUD_setLab(fd, (char *)fData[0].GetLaboratory()->Data());
5693  MUD_setArea(fd, (char *)fData[0].GetBeamline()->Data());
5694  MUD_setMethod(fd, (char *)fData[0].GetSetup()->Data());
5695  MUD_setApparatus(fd, (char *)fData[0].GetInstrument()->Data());
5696  MUD_setInsert(fd, dummy);
5697  MUD_setSample(fd, (char *)fData[0].GetSample()->Data());
5698  MUD_setOrient(fd, (char *)fData[0].GetOrientation()->Data());
5699  MUD_setDas(fd, dummy);
5700  MUD_setExperimenter(fd, dummy);
5701  snprintf(info, sizeof(info), "%lf+-%lf (K)", fData[0].GetTemperature(0), fData[0].GetTempError(0));
5702  MUD_setTemperature(fd, info);
5703  snprintf(info, sizeof(info), "%lf", fData[0].GetField());
5704  MUD_setField(fd, info);
5705 
5706  // generate the histograms
5707  MUD_setHists(fd, MUD_GRP_TRI_TD_HIST_ID, fData[0].GetNoOfHistos());
5708 
5709  UInt_t *data, dataSize = fData[0].GetDataSet(0, false)->GetData()->size()/fAny2ManyInfo->rebin + 1;
5710  data = new UInt_t[dataSize];
5711  if (data == nullptr) {
5712  std::cerr << std::endl << ">> PRunDataHandler::WriteMudFile(): **ERROR** couldn't allocate memory for the data ..." << std::endl;
5713  MUD_closeWrite(fd);
5714  return false;
5715  }
5716 
5717  UInt_t noOfEvents = 0, k = 0;
5718  ival = 0;
5719  PRawRunDataSet *dataSet;
5720  for (UInt_t i=0; i<fData[0].GetNoOfHistos(); i++) {
5721  dataSet = fData[0].GetDataSet(i, false); // i.e. the false means, that i is the index and NOT the histo number
5722  if (dataSet == nullptr) { // something is really wrong
5723  std::cerr << std::endl << ">> PRunDataHandler::WriteMudFile: **ERROR** Couldn't get data set (idx=" << i << ")";
5724  std::cerr << std::endl << ">> something is really wrong!" << std::endl;
5725  return false;
5726  }
5727 
5728  // fill data
5729  for (UInt_t j=0; j<dataSize; j++)
5730  data[j] = 0;
5731  noOfEvents = 0;
5732  k = 0;
5733  for (UInt_t j=0; j<dataSet->GetData()->size(); j++) {
5734  if ((j > 0) && (j % fAny2ManyInfo->rebin == 0)) {
5735  data[k] = ival;
5736  noOfEvents += ival;
5737  k++;
5738  ival = 0;
5739  }
5740  ival += static_cast<UInt_t>(dataSet->GetData()->at(j));
5741  }
5742 
5743  // feed data relevant information
5744  // the numbering of the histograms start from '1', hence i+1 needed!!
5745  MUD_setHistType(fd, i+1, MUD_GRP_TRI_TD_HIST_ID);
5746  MUD_setHistNumBytes(fd, i+1, sizeof(data));
5747  MUD_setHistNumBins(fd, i+1, dataSize);
5748  MUD_setHistBytesPerBin(fd, i+1, 0);
5749  MUD_setHistFsPerBin(fd, i+1, static_cast<UINT32>(1.0e6*fAny2ManyInfo->rebin*fData[0].GetTimeResolution())); // time resolution is given in (ns)
5750  MUD_setHistT0_Ps(fd, i+1, static_cast<UINT32>(1.0e3*fData[0].GetTimeResolution()*((dataSet->GetTimeZeroBin()+fAny2ManyInfo->rebin/2)/fAny2ManyInfo->rebin)));
5751  MUD_setHistT0_Bin(fd, i+1, static_cast<UINT32>(dataSet->GetTimeZeroBin()/fAny2ManyInfo->rebin));
5752  MUD_setHistGoodBin1(fd, i+1, static_cast<UINT32>(dataSet->GetFirstGoodBin()/fAny2ManyInfo->rebin));
5753  MUD_setHistGoodBin2(fd, i+1, static_cast<UINT32>(dataSet->GetLastGoodBin()/fAny2ManyInfo->rebin));
5754  MUD_setHistBkgd1(fd, i+1, 0);
5755  MUD_setHistBkgd2(fd, i+1, 0);
5756  MUD_setHistNumEvents(fd, i+1, (UINT32)noOfEvents);
5757  MUD_setHistTitle(fd, i+1, (char *)dataSet->GetName().Data());
5758  REAL64 timeResolution = (fAny2ManyInfo->rebin*fData[0].GetTimeResolution())/1.0e9; // ns -> s
5759  MUD_setHistSecondsPerBin(fd, i+1, timeResolution);
5760 
5761  MUD_setHistData(fd, i+1, data);
5762  }
5763 
5764  MUD_closeWrite(fd);
5765 
5766  delete [] data;
5767 
5768  // check if mud file shall be streamed to stdout
5770  // stream file to stdout
5771  std::ifstream is;
5772  int length=1024;
5773  char *buffer;
5774 
5775  is.open(fln.Data(), std::ios::binary);
5776  if (!is.is_open()) {
5777  std::cerr << std::endl << "PRunDataHandler::WriteMudFile(): **ERROR** Couldn't open the mud-file for streaming." << std::endl;
5778  remove(fln.Data());
5779  return false;
5780  }
5781 
5782  // get length of file
5783  is.seekg(0, std::ios::end);
5784  length = is.tellg();
5785  is.seekg(0, std::ios::beg);
5786 
5787  if (length == -1) {
5788  std::cerr << std::endl << "PRunDataHandler::WriteMudFile(): **ERROR** Couldn't determine the mud-file size." << std::endl;
5789  remove(fln.Data());
5790  return false;
5791  }
5792 
5793  // allocate memory
5794  buffer = new char [length];
5795 
5796  // read data as a block
5797  while (!is.eof()) {
5798  is.read(buffer, length);
5799  std::cout.write(buffer, length);
5800  }
5801 
5802  is.close();
5803 
5804  delete [] buffer;
5805 
5806  // delete temporary root file
5807  remove(fln.Data());
5808  }
5809  return true;
5810 }
5811 
5812 //--------------------------------------------------------------------------
5813 // WriteAsciiFile (private)
5814 //--------------------------------------------------------------------------
5825 {
5826  Bool_t ok = false;
5827  fln = GenerateOutputFileName(fln, ".ascii", ok);
5828  if (!ok)
5829  return false;
5830 
5831  TString fileName = fln;
5832 
5834  std::cout << std::endl << ">> PRunDataHandler::WriteAsciiFile(): writing an ascii data file (" << fln.Data() << ") ... " << std::endl;
5835 
5836  // write ascii file
5837  std::ofstream fout;
5838  std::streambuf* strm_buffer = nullptr;
5839 
5841  // open data-file
5842  fout.open(fln.Data(), std::ofstream::out);
5843  if (!fout.is_open()) {
5844  std::cerr << std::endl << ">> PRunDataHandler::WriteAsciiFile **ERROR** Couldn't open data file (" << fln.Data() << ") for writing, sorry ...";
5845  std::cerr << std::endl;
5846  return false;
5847  }
5848 
5849  // save output buffer of the stream
5850  strm_buffer = std::cout.rdbuf();
5851 
5852  // redirect output into the file
5853  std::cout.rdbuf(fout.rdbuf());
5854  }
5855 
5856  // write header
5857  std::cout << "%*************************************************************************";
5858  std::cout << std::endl << "% file name : " << fileName.Data();
5859  if (fData[0].GetRunTitle()->Length() > 0)
5860  std::cout << std::endl << "% title : " << fData[0].GetRunTitle()->Data();
5861  if (fData[0].GetRunNumber() >= 0)
5862  std::cout << std::endl << "% run number : " << fData[0].GetRunNumber();
5863  if (fData[0].GetSetup()->Length() > 0)
5864  std::cout << std::endl << "% setup : " << fData[0].GetSetup()->Data();
5865  std::cout << std::endl << "% field : " << fData[0].GetField() << " (G)";
5866  if (fData[0].GetStartTime()->Length() > 0)
5867  std::cout << std::endl << "% date : " << fData[0].GetStartTime()->Data() << " " << fData[0].GetStartDate()->Data() << " / " << fData[0].GetStopTime()->Data() << " " << fData[0].GetStopDate()->Data();
5868  if (fData[0].GetNoOfTemperatures() > 0) {
5869  std::cout << std::endl << "% temperature : ";
5870  for (UInt_t i=0; i<fData[0].GetNoOfTemperatures()-1; i++) {
5871  std::cout << fData[0].GetTemperature(i) << "+-" << fData[0].GetTempError(i) << "(K) , ";
5872  }
5873  std::cout << fData[0].GetTemperature(fData[0].GetNoOfTemperatures()-1) << "+-" << fData[0].GetTempError(fData[0].GetNoOfTemperatures()-1) << "(K)";
5874  }
5875  if (fData[0].GetEnergy() != PMUSR_UNDEFINED)
5876  std::cout << std::endl << "% energy : " << fData[0].GetEnergy() << " (keV)";
5877  if (fData[0].GetTransport() != PMUSR_UNDEFINED)
5878  std::cout << std::endl << "% transport : " << fData[0].GetTransport() << " (kV)";
5879  if (fData[0].GetTimeResolution() != PMUSR_UNDEFINED) {
5880  std::cout.precision(10);
5881  std::cout << std::endl << "% time resolution : " << fData[0].GetTimeResolution()*fAny2ManyInfo->rebin << " (ns)";
5882  std::cout.setf(std::ios::fixed,std::ios::floatfield); // floatfield set to fixed
5883  }
5884  PRawRunDataSet *dataSet;
5885  std::cout << std::endl << "% t0 : ";
5886  for (UInt_t i=0; i<fData[0].GetNoOfHistos()-1; i++) {
5887  dataSet = fData[0].GetDataSet(i, false); // i.e. the false means, that i is the index and NOT the histo number
5888  std::cout << static_cast<UInt_t>(dataSet->GetTimeZeroBin()/fAny2ManyInfo->rebin) << ", ";
5889  }
5890  dataSet = fData[0].GetDataSet(fData[0].GetNoOfHistos()-1, false); // i.e. the false means, that i is the index and NOT the histo number
5891  std::cout << static_cast<UInt_t>(dataSet->GetTimeZeroBin()/fAny2ManyInfo->rebin);
5892 
5893  std::cout << std::endl << "% # histos : " << fData[0].GetNoOfHistos();
5894  dataSet = fData[0].GetDataSet(0, false); // i.e. the false means, that i is the index and NOT the histo number
5895  std::cout << std::endl << "% # of bins : " << static_cast<UInt_t>(dataSet->GetData()->size()/fAny2ManyInfo->rebin);
5896  std::cout << std::endl << "%*************************************************************************";
5897 
5898  // write data
5899  UInt_t length = fData[0].GetDataSet(0,false)->GetData()->size();
5900  if (fAny2ManyInfo->rebin == 1) {
5901  for (UInt_t i=0; i<length; i++) {
5902  std::cout << std::endl;
5903  for (UInt_t j=0; j<fData[0].GetNoOfHistos(); j++) {
5904  dataSet = fData[0].GetDataSet(j, false); // i.e. the false means, that i is the index and NOT the histo number
5905  std::cout.width(8);
5906  std::cout << static_cast<Int_t>(dataSet->GetData()->at(i));
5907  }
5908  }
5909  } else {
5910  PUIntVector dataRebin(fData[0].GetNoOfHistos());
5911 
5912  // initialize the dataRebin vector
5913  for (UInt_t i=0; i<dataRebin.size(); i++) {
5914  dataSet = fData[0].GetDataSet(i, false); // i.e. the false means, that i is the index and NOT the histo number
5915  dataRebin[i] = static_cast<UInt_t>(dataSet->GetData()->at(0));
5916  }
5917 
5918  for (UInt_t i=0; i<length; i++) {
5919  if ((i > 0) && ((i % fAny2ManyInfo->rebin) == 0)) {
5920  std::cout << std::endl;
5921  for (UInt_t j=0; j<dataRebin.size(); j++) {
5922  std::cout.width(8);
5923  std::cout << dataRebin[j];
5924  }
5925  // initialize the dataRebin vector for the next pack
5926  for (UInt_t j=0; j<dataRebin.size(); j++) {
5927  dataSet = fData[0].GetDataSet(j, false); // i.e. the false means, that i is the index and NOT the histo number
5928  dataRebin[j] = static_cast<UInt_t>(dataSet->GetData()->at(i));
5929  }
5930  } else {
5931  for (UInt_t j=0; j<fData[0].GetNoOfHistos(); j++) {
5932  dataSet = fData[0].GetDataSet(j, false); // i.e. the false means, that i is the index and NOT the histo number
5933  dataRebin[j] += static_cast<UInt_t>(dataSet->GetData()->at(i));
5934  }
5935  }
5936  }
5937  }
5938 
5939  std::cout << std::endl;
5940 
5942  // restore old output buffer
5943  std::cout.rdbuf(strm_buffer);
5944 
5945  fout.close();
5946  }
5947 
5948  return true;
5949 }
5950 
5951 
5952 //--------------------------------------------------------------------------
5953 // StripWhitespace (private)
5954 //--------------------------------------------------------------------------
5966 {
5967  Char_t *s = nullptr;
5968  Char_t *subs = nullptr;
5969  Int_t i;
5970  Int_t start;
5971  Int_t end;
5972  Int_t size;
5973 
5974  size = static_cast<Int_t>(str.Length());
5975  s = new Char_t[size+1];
5976 
5977  if (!s)
5978  return false;
5979 
5980  for (Int_t i=0; i<size+1; i++)
5981  s[i] = str[i];
5982  s[size] = 0;
5983 
5984  // check for whitespaces at the beginning of the string
5985  i = 0;
5986  while (isblank(s[i]) || iscntrl(s[i])) {
5987  i++;
5988  }
5989  start = i;
5990 
5991  // check for whitespaces at the end of the string
5992  i = strlen(s);
5993  while (isblank(s[i]) || iscntrl(s[i])) {
5994  i--;
5995  }
5996  end = i;
5997 
5998  if (end < start)
5999  return false;
6000 
6001  // make substring
6002  subs = new Char_t[end-start+2];
6003  if (!subs)
6004  return false;
6005 
6006  strncpy(subs, s+start, end-start+1);
6007  subs[end-start+1] = '\0';
6008 
6009  str = TString(subs);
6010 
6011  // clean up
6012  if (subs) {
6013  delete [] subs;
6014  subs = nullptr;
6015  }
6016  if (s) {
6017  delete [] s;
6018  s = nullptr;
6019  }
6020 
6021  return true;
6022 }
6023 
6024 //--------------------------------------------------------------------------
6025 // IsWhitespace (private)
6026 //--------------------------------------------------------------------------
6037 Bool_t PRunDataHandler::IsWhitespace(const Char_t *str)
6038 {
6039  UInt_t i=0;
6040 
6041  while (isblank(str[i]) || iscntrl(str[i])) {
6042  if (str[i] == 0)
6043  break;
6044  i++;
6045  }
6046 
6047  if (i == strlen(str))
6048  return true;
6049  else
6050  return false;
6051 }
6052 
6053 //--------------------------------------------------------------------------
6054 // ToDouble (private)
6055 //--------------------------------------------------------------------------
6066 Double_t PRunDataHandler::ToDouble(TString &str, Bool_t &ok)
6067 {
6068  Char_t *s;
6069  Double_t value;
6070  Int_t size, status;
6071 
6072  ok = true;
6073 
6074  size = static_cast<Int_t>(str.Length());
6075  s = new Char_t[size+1];
6076 
6077  if (!s) {
6078  ok = false;
6079  return 0.0;
6080  }
6081 
6082  // copy string; stupid way but it works
6083  for (Int_t i=0; i<size+1; i++)
6084  s[i] = str[i];
6085  s[size] = 0;
6086 
6087  // extract value
6088  status = sscanf(s, "%lf", &value);
6089  if (status != 1) {
6090  ok = false;
6091  return 0.0;
6092  }
6093 
6094  // clean up
6095  if (s) {
6096  delete [] s;
6097  s = nullptr;
6098  }
6099 
6100  return value;
6101 }
6102 
6103 //--------------------------------------------------------------------------
6104 // ToInt (private)
6105 //--------------------------------------------------------------------------
6116 Int_t PRunDataHandler::ToInt(TString &str, Bool_t &ok)
6117 {
6118  Char_t *s;
6119  Int_t value;
6120  Int_t size, status;
6121 
6122  ok = true;
6123 
6124  size = static_cast<Int_t>(str.Length());
6125  s = new Char_t[size+1];
6126 
6127  if (!s) {
6128  ok = false;
6129  return 0;
6130  }
6131 
6132  // copy string; stupid way but it works
6133  for (Int_t i=0; i<size+1; i++)
6134  s[i] = str[i];
6135  s[size] = 0;
6136 
6137  // extract value
6138  status = sscanf(s, "%d", &value);
6139  if (status != 1) {
6140  ok = false;
6141  return 0;
6142  }
6143 
6144  // clean up
6145  if (s) {
6146  delete [] s;
6147  s = nullptr;
6148  }
6149 
6150  return value;
6151 }
6152 
6153 //--------------------------------------------------------------------------
6154 // GetDataTagIndex (private)
6155 //--------------------------------------------------------------------------
6166 Int_t PRunDataHandler::GetDataTagIndex(TString &str, const PStringVector* dataTags)
6167 {
6168  Int_t result = -1;
6169 
6170  // check all the other possible data tags
6171  for (UInt_t i=0; i<dataTags->size(); i++) {
6172  if (!dataTags->at(i).CompareTo(str, TString::kIgnoreCase)) {
6173  result = i;
6174  break;
6175  }
6176  }
6177 
6178  return result;
6179 }
6180 
6181 
6182 //--------------------------------------------------------------------------
6183 // GenerateOutputFileName (private)
6184 //--------------------------------------------------------------------------
6198 TString PRunDataHandler::GenerateOutputFileName(const TString fileName, const TString extension, Bool_t &ok)
6199 {
6200  TString fln = fileName;
6201  ok = true;
6202 
6203  if (fAny2ManyInfo == nullptr)
6204  return fln;
6205 
6206  // generate output file name if needed
6208  if (fln.Length() == 0) {
6209  ok = false;
6210  fln = GetFileName(extension, ok);
6211  if (!ok) {
6212  fln = "";
6213  return fln;
6214  }
6215  } else {
6216  fln.Prepend(fAny2ManyInfo->outPath);
6217  }
6218 
6219  // make sure that the file name doesn't already exist as a real file
6220  if (!gSystem->AccessPathName(fln)) { // file name already exists!!
6221  Int_t count = 1;
6222  TString newFln;
6223  do {
6224  newFln = fln;
6225  newFln.Insert(newFln.Last('.'), TString::Format(".%d", count++));
6226  } while (!gSystem->AccessPathName(newFln));
6227  std::cerr << std::endl << ">> PRunDataHandler::GenerateOutputFileName **WARNING** needed to modify output filename to " << newFln << ",";
6228  std::cerr << std::endl << ">> due to potential conflict with already existing files." << std::endl << std::endl;
6229  fln = newFln;
6230  }
6231 
6232  // keep the file name if compression is whished
6233  fAny2ManyInfo->outPathFileName.push_back(fln);
6234  } else {
6235  fln = fAny2ManyInfo->outPath + TString("__tmp.") + extension;
6236  }
6237 
6238  return fln;
6239 }
6240 
6241 //--------------------------------------------------------------------------
6242 // GetFileName (private)
6243 //--------------------------------------------------------------------------
6254 TString PRunDataHandler::GetFileName(const TString extension, Bool_t &ok)
6255 {
6256  TString fileName = "";
6257  ok = true;
6258 
6259  Int_t start = fRunPathName.Last('/');
6260  Int_t end = fRunPathName.Last('.');
6261  if (end == -1) {
6262  std::cerr << std::endl << ">> PRunDataHandler::GetFileName(): **ERROR** couldn't generate the output file name ..." << std::endl;
6263  ok = false;
6264  return fileName;
6265  }
6266  // cut out the filename (get rid of the extension, and the path)
6267  Char_t str1[1024], str2[1024];
6268  strncpy(str1, fRunPathName.Data(), sizeof(str1));
6269  for (Int_t i=0; i<end-start-1; i++) {
6270  str2[i] = str1[i+start+1];
6271  }
6272  str2[end-start-1] = 0;
6273 
6274  if (fAny2ManyInfo->inFormat == fAny2ManyInfo->outFormat) { // only rebinning
6275  TString rebinStr;
6276  rebinStr += fAny2ManyInfo->rebin;
6277  fileName = fAny2ManyInfo->outPath + str2 + "_rebin" + rebinStr + extension;
6278  } else { // real conversion
6279  fileName = fAny2ManyInfo->outPath + str2 + extension;
6280  }
6281 
6282  return fileName;
6283 }
6284 
6285 //--------------------------------------------------------------------------
6286 // FileNameFromTemplate (private)
6287 //--------------------------------------------------------------------------
6305 TString PRunDataHandler::FileNameFromTemplate(TString &fileNameTemplate, Int_t run, TString &year, Bool_t &ok)
6306 {
6307  TString result("");
6308 
6309  TObjArray *tok=nullptr;
6310  TObjString *ostr;
6311  TString str;
6312 
6313  // check year string
6314  if ((year.Length() != 2) && (year.Length() != 4)) {
6315  std::cerr << std::endl << ">> PRunDataHandler::FileNameFromTemplate: **ERROR** year needs to be of the format";
6316  std::cerr << std::endl << ">> 'yy' or 'yyyy', found " << year << std::endl;
6317  return result;
6318  }
6319 
6320  // make a short year version 'yy' and a long year version 'yyyy'
6321  TString yearShort="", yearLong="";
6322  if (year.Length() == 2) {
6323  yearShort = year;
6324  yearLong = "20" ;
6325  yearLong += year;
6326  } else {
6327  yearShort = year[2];
6328  yearShort += year[3];
6329  yearLong = year;
6330  }
6331 
6332  // tokenize template string
6333  tok = fileNameTemplate.Tokenize("[]");
6334  if (tok == nullptr) {
6335  std::cerr << std::endl << ">> PRunDataHandler::FileNameFromTemplate: **ERROR** couldn't tokenize template!" << std::endl;
6336  return result;
6337  }
6338  if (tok->GetEntries()==1) {
6339  std::cerr << std::endl << ">> PRunDataHandler::FileNameFromTemplate: **WARNING** template without tags." << std::endl;
6340  }
6341 
6342  // go through the tokens and generate the result string
6343  for (Int_t i=0; i<tok->GetEntries(); i++) {
6344  ostr = dynamic_cast<TObjString*>(tok->At(i));
6345  str = ostr->GetString();
6346 
6347  // check tokens
6348  if (!str.CompareTo("yy", TString::kExact)) { // check for 'yy'
6349  result += yearShort;
6350  } else if (!str.CompareTo("yyyy", TString::kExact)) { // check for 'yyyy'
6351  result += yearLong;
6352  } else if (str.Contains("rr")) { // check for run
6353  // make sure ONLY 'r' are present
6354  Int_t idx;
6355  for (idx=0; idx<str.Length(); idx++) {
6356  if (str[idx] != 'r')
6357  break;
6358  }
6359  if (idx == str.Length()) { // 'r' only
6360  TString runStr("");
6361  char fmt[128];
6362  snprintf(fmt, sizeof(fmt), "%%0%dd", str.Length());
6363  runStr.Form(fmt, run);
6364  result += runStr;
6365  } else { // not only 'r'
6366  result += str;
6367  }
6368  } else { // all parts which are NOT tags
6369  result += str;
6370  }
6371  }
6372 
6373  // clean up
6374  if (tok)
6375  delete tok;
6376 
6377  // everything fine here
6378  ok = true;
6379 
6380  return result;
6381 }
6382 
6383 //--------------------------------------------------------------------------
6384 // DateToISO8601 (private)
6385 //--------------------------------------------------------------------------
6394 bool PRunDataHandler::DateToISO8601(std::string inDate, std::string &iso8601Date)
6395 {
6396  iso8601Date = std::string("");
6397 
6398  struct tm tm;
6399 
6400  // currently only dates like dd-mmm-yy are handled, where dd day of the month, mmm month as abbrivation, e.g. JAN, and yy as the year
6401  if (!strptime(inDate.c_str(), "%d-%b-%y", &tm)) // failed
6402  return false;
6403 
6404  TString str("");
6405  str.Form("%04d-%02d-%02d", 1900+tm.tm_year, tm.tm_mon+1, tm.tm_mday);
6406 
6407  iso8601Date = str.Data();
6408 
6409  return true;
6410 }
6411 
6412 //--------------------------------------------------------------------------
6413 // SplitTimeData (private)
6414 //--------------------------------------------------------------------------
6424 void PRunDataHandler::SplitTimeDate(TString timeData, TString &time, TString &date, Bool_t &ok)
6425 {
6426  struct tm tm;
6427  memset(&tm, 0, sizeof(tm));
6428  strptime(timeData.Data(), "%Y-%m-%d %H:%M:S", &tm);
6429  if (tm.tm_year == 0)
6430  strptime(timeData.Data(), "%Y-%m-%dT%H:%M:S", &tm);
6431 
6432  if (tm.tm_year == 0) {
6433  ok = false;
6434  return;
6435  }
6436 
6437  time = TString::Format("%02d:%02d:%02d", tm.tm_hour, tm.tm_min, tm.tm_sec);
6438  date = TString::Format("%04d-%02d-%02d", tm.tm_year+1900, tm.tm_mon+1, tm.tm_mday);
6439 }
6440 
6441 //--------------------------------------------------------------------------
6442 // GetMonth (private)
6443 //--------------------------------------------------------------------------
6451 TString PRunDataHandler::GetMonth(Int_t month)
6452 {
6453  TString monthString = "???";
6454 
6455  switch (month) {
6456  case 1:
6457  monthString = "JAN";
6458  break;
6459  case 2:
6460  monthString = "FEB";
6461  break;
6462  case 3:
6463  monthString = "MAR";
6464  break;
6465  case 4:
6466  monthString = "APR";
6467  break;
6468  case 5:
6469  monthString = "MAY";
6470  break;
6471  case 6:
6472  monthString = "JUN";
6473  break;
6474  case 7:
6475  monthString = "JUL";
6476  break;
6477  case 8:
6478  monthString = "AUG";
6479  break;
6480  case 9:
6481  monthString = "SEP";
6482  break;
6483  case 10:
6484  monthString = "OCT";
6485  break;
6486  case 11:
6487  monthString = "NOV";
6488  break;
6489  case 12:
6490  monthString = "DEC";
6491  break;
6492  default:
6493  break;
6494  }
6495 
6496  return monthString;
6497 }
6498 
6499 //--------------------------------------------------------------------------
6500 // GetYear (private)
6501 //--------------------------------------------------------------------------
6509 TString PRunDataHandler::GetYear(Int_t year)
6510 {
6511  TString yearString = "??";
6512  Int_t yy=0;
6513 
6514  if (year < 2000) {
6515  yy = year - 1900;
6516  if (yy > 0)
6517  yearString.Form("%02d", yy);
6518  } else {
6519  yy = year - 2000;
6520  if ((yy >= 0) && (yy <= 99))
6521  yearString.Form("%02d", yy);
6522  }
6523 
6524  return yearString;
6525 }
#define A2M_ASCII
virtual void SetInstrument(const TString &str)
Definition: PMusr.h:449
virtual TString * GetBeamline(UInt_t idx=0)
Definition: PMusr.cpp:1274
virtual Bool_t StripWhitespace(TString &str)
virtual void SetTimeZeroBinEstimated(Double_t tzb)
Definition: PMusr.h:337
virtual ~PRunDataHandler()
virtual Double_t ToDouble(TString &str, Bool_t &ok)
virtual void SetRunTitle(const TString str)
Definition: PMusr.h:456
#define A2M_MUSR_ROOT
virtual void SetSpecificValidatorUrl(const TString &str)
Definition: PMusr.h:443
virtual Bool_t ReadWkmFile()
TString year
holds the information about the year to be used
Definition: PMusr.h:831
virtual void SetStopDate(const TString str)
Definition: PMusr.h:462
PRawRunDataList fData
keeping all the raw data
virtual void Init(const Int_t tag=0)
virtual Bool_t ReadRootFile()
virtual void TestFileName(TString &runName, const TString &ext)
virtual void SetMagnetName(const TString str)
Definition: PMusr.h:464
virtual void SetMuonSource(const TString &str)
Definition: PMusr.h:450
virtual const Double_t GetTimeResolution()
Definition: PMusr.h:429
virtual TString GetYear(Int_t month)
std::vector< PMsrRunBlock > PMsrRunList
Definition: PMusr.h:754
virtual Bool_t WriteMudFile(TString fln="")
virtual Double_t GetValue() const
virtual TString GenerateOutputFileName(const TString fileName, const TString extension, Bool_t &ok)
virtual TString GetFileName(const TString extension, Bool_t &ok)
virtual TString * GetFileFormat(UInt_t idx=0)
Definition: PMusr.cpp:1358
virtual void SetBeamline(TString &str, Int_t idx=-1)
Definition: PMusr.cpp:1291
#define A2M_MUD
virtual void SetEnergy(const Double_t dval)
Definition: PMusr.h:472
virtual void SetStopTime(const TString str)
Definition: PMusr.h:461
virtual void SetGenericValidatorUrl(const TString &str)
Definition: PMusr.h:442
virtual Bool_t WriteMusrRootFile(TString fln="")
virtual void SetRingAnode(const UInt_t idx, const Double_t dval)
Definition: PMusr.cpp:681
TString fRunPathName
current path file name
virtual const std::vector< PDoubleVector > * GetData()
Definition: PMusr.h:293
virtual Bool_t ReadFilesMsr()
TString outFileName
holds the output file name
Definition: PMusr.h:835
virtual void SetTempError(const UInt_t idx, const Double_t errTemp)
Definition: PMusr.cpp:716
virtual void SplitTimeDate(TString timeDate, TString &time, TString &date, Bool_t &ok)
virtual void SetRunName(const TString &str)
Definition: PMusr.h:454
virtual void Clear()
Definition: PMusr.cpp:261
#define PRH_PPC_OFFSET
TString fRunName
current run name
virtual PRawRunData * GetRunData(const TString &runName)
virtual Bool_t WritePsiBinFile(TString fln="")
virtual const PStringVector * GetLabels()
Definition: PMusr.h:291
UInt_t compressionTag
0=no compression, 1=gzip compression, 2=bzip2 compression
Definition: PMusr.h:839
virtual Bool_t FileExistsCheck(PMsrRunBlock &runInfo, const UInt_t idx)
virtual TString GetMonth(Int_t month)
TString inFormat
holds the information about the input data file format
Definition: PMusr.h:827
virtual Bool_t WriteNexusFile(TString fln="")
virtual void SetStartTime(const TString str)
Definition: PMusr.h:458
Bool_t fAllDataAvailable
flag indicating if all data sets could be read
#define PRH_MUSR_ROOT
PAny2ManyInfo * fAny2ManyInfo
pointer to the any2many data structure
virtual Bool_t ReadDatFile()
virtual void SetRunNumber(const Int_t &val)
Definition: PMusr.h:455
virtual void SetFileFormat(TString &str, Int_t idx=-1)
Definition: PMusr.cpp:1375
virtual void SetBeamline(const TString &str)
Definition: PMusr.h:448
std::vector< Int_t > PIntVector
Definition: PMusr.h:178
virtual void AppendSubErrData(const UInt_t idx, const Double_t dval)
Definition: PMusr.cpp:229
PStringVector fDataPath
vector containing all the search paths where to look for data files
virtual void SetMuonSpinAngle(const Double_t dval)
Definition: PMusr.h:453
std::vector< UInt_t > PUIntVector
Definition: PMusr.h:172
virtual void SetDataSet(PRawRunDataSet &dataSet, UInt_t idx=-1)
Definition: PMusr.h:477
#define A2M_ROOT
virtual void AppendLabel(const TString str)
Definition: PMusr.h:298
std::vector< Double_t > PDoubleVector
Definition: PMusr.h:196
PStringVector outPathFileName
holds the out path/file name
Definition: PMusr.h:836
virtual PMsrRunList * GetMsrRunList()
Definition: PMsrHandler.h:63
virtual Bool_t SetRunData(PRawRunData *data, UInt_t idx=0)
#define MRH_UNDEFINED
virtual void SetStartDate(const TString str)
Definition: PMusr.h:459
virtual const PStringVector * GetDataTags()
Definition: PMusr.h:292
virtual void SetComment(const TString &str)
Definition: PMusr.h:445
virtual void SetGenerator(const TString &str)
Definition: PMusr.h:444
virtual Int_t GetLastGoodBin()
Definition: PMusr.h:328
virtual TString GetName()
Definition: PMusr.h:323
PIntVector groupHistoList
holds the histo group list offset (used to define for MusrRoot files, what to be exported) ...
Definition: PMusr.h:833
virtual bool DateToISO8601(std::string inDate, std::string &iso8601Date)
virtual Bool_t WriteAsciiFile(TString fln="")
virtual void SetRedGreenOffset(PIntVector &ivec)
Definition: PMusr.h:476
virtual TString * GetMsrFileDirectoryPath()
Definition: PMsrHandler.h:69
virtual Bool_t IsWhitespace(const Char_t *str)
virtual Int_t GetHistoNo()
Definition: PMusr.h:324
virtual TString GetUnit() const
TString outTemplate
holds the output file template
Definition: PMusr.h:830
virtual void SetSample(const TString str)
Definition: PMusr.h:467
virtual Bool_t ReadMduAsciiFile()
virtual Bool_t WriteData(TString fileName="")
virtual const std::vector< PDoubleVector > * GetErrData()
Definition: PMusr.h:294
virtual void SetVersion(const TString &str)
Definition: PMusr.h:441
virtual Int_t ToInt(TString &str, Bool_t &ok)
#define PHR_INIT_ALL
virtual Bool_t ReadAsciiFile()
TString fFileFormat
keeps the file format if explicitly given
virtual void SetTimeZeroBin(Double_t tzb)
Definition: PMusr.h:336
virtual void ClearTemperature()
Definition: PMusr.h:469
virtual void SetMuonBeamMomentum(const Double_t dval)
Definition: PMusr.h:452
virtual void Set(TString label, Double_t demand, Double_t value, Double_t error, TString unit, TString description=TString("n/a"))
virtual void SetName(TString str)
Definition: PMusr.h:334
#define PHR_INIT_MSR
virtual Double_t GetError() const
UInt_t rebin
holds the number of bins to be packed
Definition: PMusr.h:838
virtual void SetInstitute(TString &str, Int_t idx=-1)
Definition: PMusr.cpp:1333
virtual void AppendDataTag(const TString str)
Definition: PMusr.h:300
virtual void SetHistoNo(Int_t no)
Definition: PMusr.h:335
virtual Bool_t ReadNexusFile()
#define PMUSR_UNDEFINED
Definition: PMusr.h:89
virtual void SetSize(const UInt_t size)
Definition: PMusr.cpp:162
virtual void SetFirstGoodBin(Int_t fgb)
Definition: PMusr.h:338
virtual void SetStartDateTime(const time_t val)
Definition: PMusr.h:460
virtual void SetOrientation(const TString str)
Definition: PMusr.h:468
virtual Bool_t WriteRootFile(TString fln="")
#define A2M_WKM
virtual void SetFileName(const TString &str)
Definition: PMusr.h:446
virtual Bool_t ReadMudFile()
#define A2M_NEXUS
TString compressFileName
holds the name of the outputfile name in case of compression is used
Definition: PMusr.h:840
TString outPath
holds the output path
Definition: PMusr.h:837
virtual void AppendSubData(const UInt_t idx, const Double_t dval)
Definition: PMusr.cpp:209
virtual void SetData(PDoubleVector data)
Definition: PMusr.h:342
Bool_t useStandardOutput
flag showing if the converted shall be sent to the standard output
Definition: PMusr.h:826
virtual void SetStopDateTime(const time_t val)
Definition: PMusr.h:463
virtual PDoubleVector * GetData()
Definition: PMusr.h:331
#define PRH_LEM_ROOT
virtual void AppendErrData(const PDoubleVector &data)
Definition: PMusr.h:302
virtual Bool_t ReadPsiBinFile()
virtual TString * GetInstitute(UInt_t idx=0)
Definition: PMusr.cpp:1316
virtual Bool_t ReadDBFile()
virtual void SetLabel(const UInt_t idx, const TString str)
Definition: PMusr.cpp:189
virtual void SetCryoName(const TString str)
Definition: PMusr.h:466
virtual void SetSetup(const TString str)
Definition: PMusr.h:457
#define A2M_PSIBIN
virtual Bool_t WriteWkmFile(TString fln="")
PUIntVector runList
holds the run number list to be converted
Definition: PMusr.h:832
virtual TString FileNameFromTemplate(TString &fileNameTemplate, Int_t run, TString &year, Bool_t &ok)
virtual Int_t GetFirstGoodBin()
Definition: PMusr.h:327
virtual void SetLaboratory(const TString &str)
Definition: PMusr.h:447
virtual void SetField(const Double_t dval)
Definition: PMusr.h:465
#define A2M_PSIMDU
virtual void SetLastBkgBin(Int_t lbb)
Definition: PMusr.h:341
virtual const PDoubleVector * GetDataBin(const UInt_t histoNo)
Definition: PMusr.h:438
virtual Int_t GetDataTagIndex(TString &str, const PStringVector *fLabels)
virtual const UInt_t GetNoOfHistos()
Definition: PMusr.h:436
std::vector< TString > PStringVector
Definition: PMusr.h:214
virtual void SetLastGoodBin(Int_t lgb)
Definition: PMusr.h:339
virtual void AppendData(const PDoubleVector &data)
Definition: PMusr.h:301
virtual void SetFromAscii(const Bool_t bval)
Definition: PMusr.h:296
return status
virtual Double_t GetTimeZeroBin()
Definition: PMusr.h:325
PNonMusrRawRunData fDataNonMusr
keeps all ascii- or db-file info in case of nonMusr fit
Definition: PMusr.h:479
#define A2M_UNDEFINED
virtual Bool_t ReadWriteFilesList()
PMsrHandler * fMsrInfo
pointer to the msr-file handler
virtual void ConvertData()
#define PHR_INIT_ANY2MANY
virtual void SetMuonSpecies(const TString &str)
Definition: PMusr.h:451
virtual void SetTimeResolution(const Double_t dval)
Definition: PMusr.h:475
UInt_t idf
IDF version for NeXus files.
Definition: PMusr.h:841
virtual void SetTemperature(const UInt_t idx, const Double_t temp, const Double_t errTemp)
Definition: PMusr.cpp:698
TString outFormat
holds the information about the output data file format
Definition: PMusr.h:828
virtual void SetTransport(const Double_t dval)
Definition: PMusr.h:473
virtual TString * GetRunName(UInt_t idx=0)
Definition: PMusr.cpp:1232
virtual void ReadData()
PStringVector inFileName
holds the file name of the input data file
Definition: PMusr.h:834
virtual Bool_t FileAlreadyRead(TString runName)
virtual void SetFirstBkgBin(Int_t fbb)
Definition: PMusr.h:340
std::vector< Bool_t > PBoolVector
Definition: PMusr.h:166
#define POST_PILEUP_HISTO_OFFSET
Definition: PMusr.h:86
TString inTemplate
holds the input file template
Definition: PMusr.h:829