REST-for-Physics  v2.3
Rare Event Searches ToolKit for Physics
TRestAxionSolarFlux.cxx
1 /******************** REST disclaimer ***********************************
2  * This file is part of the REST software framework. *
3  * *
4  * Copyright (C) 2016 GIFNA/TREX (University of Zaragoza) *
5  * For more information see http://gifna.unizar.es/trex *
6  * *
7  * REST is free software: you can redistribute it and/or modify *
8  * it under the terms of the GNU General Public License as published by *
9  * the Free Software Foundation, either version 3 of the License, or *
10  * (at your option) any later version. *
11  * *
12  * REST is distributed in the hope that it will be useful, *
13  * but WITHOUT ANY WARRANTY; without even the implied warranty of *
14  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
15  * GNU General Public License for more details. *
16  * *
17  * You should have a copy of the GNU General Public License along with *
18  * REST in $REST_PATH/LICENSE. *
19  * If not, see http://www.gnu.org/licenses/. *
20  * For the list of contributors see $REST_PATH/CREDITS. *
21  *************************************************************************/
22 
84 
85 #include "TRestAxionSolarFlux.h"
86 using namespace std;
87 
88 ClassImp(TRestAxionSolarFlux);
89 
94 
109 TRestAxionSolarFlux::TRestAxionSolarFlux(const char* cfgFileName, string name) : TRestMetadata(cfgFileName) {
110  RESTDebug << "Entering TRestAxionSolarFlux constructor( cfgFileName, name )" << RESTendl;
111  RESTDebug << "File: " << cfgFileName << " Name: " << name << RESTendl;
112 }
113 
118 
123  SetLibraryVersion(LIBRARY_VERSION);
124 
125  fTablesLoaded = false;
126  if (LoadTables()) fTablesLoaded = true;
127 
128  if (!fRandom) {
129  delete fRandom;
130  fRandom = nullptr;
131  }
132 
133  if (fRandom != nullptr) {
134  delete fRandom;
135  fRandom = nullptr;
136  }
137 
138  fRandom = new TRandom3(fSeed);
139  fSeed = fRandom->TRandom::GetSeed();
140 }
141 
146 TH1F* TRestAxionSolarFlux::GetFluxHistogram(string fname, Double_t binSize) {
147  string fullPathName = SearchFile(fname);
148 
149  std::vector<std::vector<Double_t>> fluxData;
150  TRestTools::ReadASCIITable(fullPathName, fluxData, 3);
151 
152  TH2F* originalHist =
153  new TH2F(Form("FluxTable_%s", GetName()), "", 100, 0., 1., (Int_t)(20. / binSize), 0., 20.);
154 
155  for (const auto& data : fluxData) {
156  Double_t r = 0.005 + data[0];
157  Double_t en = data[1] - 0.005;
158  Double_t flux = data[2] * binSize; // flux in cm-2 s-1 bin-1
159 
160  originalHist->Fill(r, en, flux);
161  }
162 
163  return (TH1F*)originalHist->ProjectionY();
164 }
165 
170 TCanvas* TRestAxionSolarFlux::DrawFluxFile(string fname, Double_t binSize) {
171  if (fCanvas != nullptr) {
172  delete fCanvas;
173  fCanvas = nullptr;
174  }
175  fCanvas = new TCanvas("canv", "This is the canvas title", 1400, 1200);
176  fCanvas->Draw();
177 
178  TPad* pad1 = new TPad("pad1", "This is pad1", 0.01, 0.02, 0.99, 0.97);
179  pad1->Draw();
180 
181  fCanvas->cd();
182  pad1->cd();
183 
184  GetFluxHistogram(fname, binSize)->Draw("hist");
185 
186  return fCanvas;
187 }
188 
193  if (fCanvas != nullptr) {
194  delete fCanvas;
195  fCanvas = nullptr;
196  }
197  fCanvas = new TCanvas("canv", "This is the canvas title", 1200, 500);
198  fCanvas->Draw();
199 
200  TPad* pad1 = new TPad("pad1", "This is pad1", 0.01, 0.02, 0.99, 0.97);
201  pad1->Draw();
202 
203  pad1->cd();
204  pad1->SetLogy();
205  pad1->SetRightMargin(0.09);
206  pad1->SetLeftMargin(0.15);
207  pad1->SetBottomMargin(0.15);
208 
209  TH1F* ht = GetEnergySpectrum();
210  ht->SetLineColor(kBlack);
211  ht->SetFillStyle(4050);
212  ht->SetFillColor(kBlue - 10);
213 
214  ht->Draw("hist");
215 
216  return fCanvas;
217 }
218 
224 
225  RESTMetadata << " - Coupling type : " << fCouplingType << RESTendl;
226  RESTMetadata << " - Coupling strength : " << fCouplingStrength << RESTendl;
227  RESTMetadata << "--------" << RESTendl;
228  RESTMetadata << " - Random seed : " << fSeed << RESTendl;
229  RESTMetadata << "--------" << RESTendl;
230 }
A metadata class to load tabulated solar axion fluxes.
void Initialize()
It is required in order to load solar flux tables into memory.
TRestAxionSolarFlux()
Default constructor.
TCanvas * fCanvas
A canvas pointer for drawing.
~TRestAxionSolarFlux()
Default destructor.
virtual Bool_t LoadTables()=0
It defines how to read the solar tables at the inhereted class.
Double_t fCouplingStrength
Axion coupling strength.
TCanvas * DrawFluxFile(std::string fname, Double_t binSize=0.01)
It draws the contents of a .flux file. This method just receives the name of the ....
TH1F * GetFluxHistogram(std::string fname, Double_t binSize=0.01)
It builds a histogram using the contents of the .flux file given in the argument.
std::string fCouplingType
Axion coupling. Defines coupling type and strength.
Int_t fSeed
Seed used in random generator.
TRandom3 * fRandom
Random number generator.
virtual TH1F * GetEnergySpectrum(Double_t m=0)=0
It returns an energy integrated spectrum in cm-2 s-1 keV-1.
void PrintMetadata()
Prints on screen the information about the metadata members of TRestAxionSolarFlux.
Bool_t fTablesLoaded
A metadata member to control if this class has been initialized.
virtual TCanvas * DrawSolarFlux()
It draws the contents of a .flux file. This method just receives the.
A base class for any REST metadata class.
Definition: TRestMetadata.h:74
virtual void PrintMetadata()
Implemented it in the derived metadata class to print out specific metadata information.
endl_t RESTendl
Termination flag object for TRestStringOutput.
void SetLibraryVersion(TString version)
Set the library version of this metadata class.
std::string SearchFile(std::string filename)
Search files in current directory and directories specified in "searchPath" section.
static int ReadASCIITable(std::string fName, std::vector< std::vector< Double_t >> &data, Int_t skipLines=0, std::string separator="\t")
Reads an ASCII file containing a table with values.
Definition: TRestTools.cxx:577