Skip to content

Commit

Permalink
[pre-commit.ci] auto fixes from pre-commit.com hooks
Browse files Browse the repository at this point in the history
for more information, see https://pre-commit.ci
  • Loading branch information
pre-commit-ci[bot] committed Apr 11, 2024
1 parent fd357f2 commit 8032860
Show file tree
Hide file tree
Showing 6 changed files with 191 additions and 209 deletions.
14 changes: 7 additions & 7 deletions inc/TRestAxionMagnetModel.h
Original file line number Diff line number Diff line change
Expand Up @@ -23,21 +23,21 @@
#ifndef REST_TRestAxionMagnetModel
#define REST_TRestAxionMagnetModel

#include <TRestMetadata.h>
#include <TRestAxionMagneticFit.h>
#include <TRestMetadata.h>

/// A class to fit and store the fit parameters from the magnetic field map
class TRestAxionMagnetModel: public TObject {
class TRestAxionMagnetModel : public TObject {
private:
std::vector <TRestAxionMagneticFit> fLinearFit;
std::vector<TRestAxionMagneticFit> fLinearFit;

public:
virtual std::vector<Double_t> GetComponentAlongPath(Int_t axis, TVector3 from, TVector3 to, Double_t dl = 1.,
Int_t Nmax = 0) = 0;
virtual std::vector<Double_t> GetComponentAlongPath(Int_t axis, TVector3 from, TVector3 to,
Double_t dl = 1., Int_t Nmax = 0) = 0;

void Test(Double_t X = 340, Double_t Y = 0, Double_t dl = 50);
void Test(Double_t X = 340, Double_t Y = 0, Double_t dl = 50);

const TRestAxionMagneticFit &GetFit( size_t n ) const{ return fLinearFit[n]; }
const TRestAxionMagneticFit& GetFit(size_t n) const { return fLinearFit[n]; }

ClassDef(TRestAxionMagnetModel, 1);
};
Expand Down
7 changes: 2 additions & 5 deletions inc/TRestAxionMagneticField.h
Original file line number Diff line number Diff line change
Expand Up @@ -23,7 +23,8 @@
#ifndef _TRestAxionMagneticField
#define _TRestAxionMagneticField


#include <TRestAxionMagnetModel.h>
#include <TRestMetadata.h>

#include <iostream>

Expand All @@ -33,10 +34,6 @@
#include "TVector3.h"
#include "TVectorD.h"

#include <TRestMetadata.h>

#include <TRestAxionMagnetModel.h>

/// A structure to define the properties and store the field data of a single magnetic volume inside
/// TRestAxionMagneticField
struct MagneticFieldVolume {
Expand Down
4 changes: 2 additions & 2 deletions inc/TRestAxionMagneticFit.h
Original file line number Diff line number Diff line change
Expand Up @@ -124,8 +124,8 @@ class TRestAxionMagneticFit {
void LoadData(const std::vector<Double_t>& z, const std::vector<Double_t>& Bx,
const std::vector<Double_t>& By, const std::vector<Double_t>& Bz);

Double_t GetChi2_X();
Double_t GetChi2_X2() { return fChiSquareBx; }
Double_t GetChi2_X();
Double_t GetChi2_X2() { return fChiSquareBx; }

TRestAxionMagneticFit();
~TRestAxionMagneticFit();
Expand Down
252 changes: 121 additions & 131 deletions macros/REST_Axion_PaperPlots.C
100755 → 100644
Original file line number Diff line number Diff line change
@@ -1,6 +1,6 @@
#include <TRestTools.h>
//*******************************************************************************************************
//*** Description:
//*** Description:
//***
//*** Arguments by default are (in order):
//***
Expand All @@ -20,147 +20,137 @@
const Int_t samples = 1000000;
std::string fluxesString = "LennertHoofPrimakoff:LennertHoofABC";

int SolarPlots()
{
int SolarPlots() {
std::vector<std::string> fluxes = TRestTools::GetOptions(fluxesString);

TCanvas *c1 = new TCanvas("c1", "My canvas", 800, 350);
c1->GetFrame()->SetBorderSize(6);
c1->GetFrame()->SetBorderMode(-1);

TPad *pad1 = new TPad("pad1", "This is pad1", 0.01, 0.02, 0.99, 0.97);
pad1->Divide(2, 1);
pad1->Draw();

TRestAxionSolarQCDFlux *solarFlux1 = new TRestAxionSolarQCDFlux("fluxes.rml", fluxes[0]);
solarFlux1->Initialize();
solarFlux1->PrintMetadata();

TRestAxionSolarQCDFlux *solarFlux2 = new TRestAxionSolarQCDFlux("fluxes.rml", fluxes[1]);
solarFlux2->Initialize();
solarFlux2->PrintMetadata();

Int_t enBins = 2000;
Double_t enMax = 20;

//// Producing histograms ///
TH2D *EvsRHisto1 = new TH2D("FluxEvsR1", "", enBins, 0, enMax, 100, 0, 1);
for (int n = 0; n < samples; n++ )
{
std::pair<Double_t,Double_t> x = solarFlux1->GetRandomEnergyAndRadius();
EvsRHisto1->Fill(x.first, x.second);
}

TH2D *EvsRHisto2 = new TH2D("FluxEvsR2", "", enBins, 0, enMax, 100, 0, 1);
for (int n = 0; n < samples; n++ )
{
std::pair<Double_t,Double_t> x = solarFlux2->GetRandomEnergyAndRadius();
EvsRHisto2->Fill(x.first, x.second);
}

TRandom3 *rnd = new TRandom3(0);
TH2D *solarDisk = new TH2D("solar_disk", "SolarDisk", 120, -1.2, 1.2, 120, -1.2, 1.2);
for (int n = 0; n < samples; n++ )
{
Double_t angle = rnd->Rndm() * 2 * TMath::Pi();
std::pair<Double_t,Double_t> x = solarFlux1->GetRandomEnergyAndRadius();

solarDisk->Fill(x.second * TMath::Cos(angle), x.second * TMath::Sin(angle));
}

TH1D *enSpt1 = EvsRHisto1->ProjectionX();
TH1D *rSpt1 = EvsRHisto1->ProjectionY();

TH1D *enSpt2 = EvsRHisto2->ProjectionX();
TH1D *rSpt2 = EvsRHisto2->ProjectionY();

pad1->cd(1);
pad1->cd(1)->SetRightMargin(0.12);
pad1->cd(1)->SetLeftMargin(0.15);
pad1->cd(1)->SetBottomMargin(0.15);

EvsRHisto1->SetStats(0);
EvsRHisto1->GetXaxis()->SetTitle("Energy [keV]");
EvsRHisto1->GetXaxis()->SetTitleSize(0.05);
EvsRHisto1->GetXaxis()->SetLabelSize(0.05);
EvsRHisto1->GetYaxis()->SetTitle("Solar radius");
EvsRHisto1->GetYaxis()->SetTitleSize(0.05);
EvsRHisto1->GetYaxis()->SetLabelSize(0.05);
EvsRHisto1->Draw("colz");

pad1->cd(2);
//pad1->cd(2)->SetLogy();
pad1->cd(2)->SetRightMargin(0.09);
pad1->cd(2)->SetLeftMargin(0.15);
pad1->cd(2)->SetBottomMargin(0.15);

enSpt1->GetYaxis()->SetTitleSize(0.05);
enSpt1->SetStats(0);
enSpt1->GetXaxis()->SetTitle("Energy [keV]");
enSpt1->GetXaxis()->SetTitleSize(0.05);
enSpt1->GetXaxis()->SetLabelSize(0.05);
enSpt1->GetYaxis()->SetTitle("Counts [keV^{-1}]");
enSpt1->GetYaxis()->SetTitleSize(0.05);
enSpt1->GetYaxis()->SetLabelSize(0.05);
//enSpt1->SetFillStyle(4050);
//enSpt1->SetFillColor(kBlue - 9);
enSpt1->SetFillColorAlpha(kBlue,0.15);
enSpt1->SetLineColor(kBlack);
enSpt1->Scale(enBins/enMax);
TCanvas* c1 = new TCanvas("c1", "My canvas", 800, 350);
c1->GetFrame()->SetBorderSize(6);
c1->GetFrame()->SetBorderMode(-1);

TPad* pad1 = new TPad("pad1", "This is pad1", 0.01, 0.02, 0.99, 0.97);
pad1->Divide(2, 1);
pad1->Draw();

TRestAxionSolarQCDFlux* solarFlux1 = new TRestAxionSolarQCDFlux("fluxes.rml", fluxes[0]);
solarFlux1->Initialize();
solarFlux1->PrintMetadata();

TRestAxionSolarQCDFlux* solarFlux2 = new TRestAxionSolarQCDFlux("fluxes.rml", fluxes[1]);
solarFlux2->Initialize();
solarFlux2->PrintMetadata();

Int_t enBins = 2000;
Double_t enMax = 20;

//// Producing histograms ///
TH2D* EvsRHisto1 = new TH2D("FluxEvsR1", "", enBins, 0, enMax, 100, 0, 1);
for (int n = 0; n < samples; n++) {
std::pair<Double_t, Double_t> x = solarFlux1->GetRandomEnergyAndRadius();
EvsRHisto1->Fill(x.first, x.second);
}

TH2D* EvsRHisto2 = new TH2D("FluxEvsR2", "", enBins, 0, enMax, 100, 0, 1);
for (int n = 0; n < samples; n++) {
std::pair<Double_t, Double_t> x = solarFlux2->GetRandomEnergyAndRadius();
EvsRHisto2->Fill(x.first, x.second);
}

TRandom3* rnd = new TRandom3(0);
TH2D* solarDisk = new TH2D("solar_disk", "SolarDisk", 120, -1.2, 1.2, 120, -1.2, 1.2);
for (int n = 0; n < samples; n++) {
Double_t angle = rnd->Rndm() * 2 * TMath::Pi();
std::pair<Double_t, Double_t> x = solarFlux1->GetRandomEnergyAndRadius();

solarDisk->Fill(x.second * TMath::Cos(angle), x.second * TMath::Sin(angle));
}

TH1D* enSpt1 = EvsRHisto1->ProjectionX();
TH1D* rSpt1 = EvsRHisto1->ProjectionY();

TH1D* enSpt2 = EvsRHisto2->ProjectionX();
TH1D* rSpt2 = EvsRHisto2->ProjectionY();

pad1->cd(1);
pad1->cd(1)->SetRightMargin(0.12);
pad1->cd(1)->SetLeftMargin(0.15);
pad1->cd(1)->SetBottomMargin(0.15);

EvsRHisto1->SetStats(0);
EvsRHisto1->GetXaxis()->SetTitle("Energy [keV]");
EvsRHisto1->GetXaxis()->SetTitleSize(0.05);
EvsRHisto1->GetXaxis()->SetLabelSize(0.05);
EvsRHisto1->GetYaxis()->SetTitle("Solar radius");
EvsRHisto1->GetYaxis()->SetTitleSize(0.05);
EvsRHisto1->GetYaxis()->SetLabelSize(0.05);
EvsRHisto1->Draw("colz");

pad1->cd(2);
// pad1->cd(2)->SetLogy();
pad1->cd(2)->SetRightMargin(0.09);
pad1->cd(2)->SetLeftMargin(0.15);
pad1->cd(2)->SetBottomMargin(0.15);

enSpt1->GetYaxis()->SetTitleSize(0.05);
enSpt1->SetStats(0);
enSpt1->GetXaxis()->SetTitle("Energy [keV]");
enSpt1->GetXaxis()->SetTitleSize(0.05);
enSpt1->GetXaxis()->SetLabelSize(0.05);
enSpt1->GetYaxis()->SetTitle("Counts [keV^{-1}]");
enSpt1->GetYaxis()->SetTitleSize(0.05);
enSpt1->GetYaxis()->SetLabelSize(0.05);
// enSpt1->SetFillStyle(4050);
// enSpt1->SetFillColor(kBlue - 9);
enSpt1->SetFillColorAlpha(kBlue, 0.15);
enSpt1->SetLineColor(kBlack);
enSpt1->Scale(enBins / enMax);
enSpt1->GetYaxis()->SetRangeUser(0., 500000);
enSpt1->Draw("HIST FL");
//enSpt2->SetFillColor(kRed - 9);
enSpt2->SetFillColorAlpha(kRed,0.15);
enSpt2->SetLineColor(kBlack);
enSpt2->Scale(enBins/enMax);
enSpt2->Draw("HIST SAME");

TLegend* legend = new TLegend(0.87, 0.75, 0.55, 0.9);
legend->AddEntry(enSpt1, "Primakoff", "lf");
legend->AddEntry(enSpt2, "ABC", "lf");
legend->Draw("same");
enSpt1->Draw("HIST FL");
// enSpt2->SetFillColor(kRed - 9);
enSpt2->SetFillColorAlpha(kRed, 0.15);
enSpt2->SetLineColor(kBlack);
enSpt2->Scale(enBins / enMax);
enSpt2->Draw("HIST SAME");

TLegend* legend = new TLegend(0.87, 0.75, 0.55, 0.9);
legend->AddEntry(enSpt1, "Primakoff", "lf");
legend->AddEntry(enSpt2, "ABC", "lf");
legend->Draw("same");

c1->Print("/tmp/solarFlux.pdf");
return 0;
return 0;
}

int MagnetPlots()
{
TRestAxionMagneticField *field = new TRestAxionMagneticField( "fields.rml", "babyIAXO_2024_cutoff" );
TCanvas *c = field->DrawTracks( TVector3(0,0,8000), 100 );
c->Print("/tmp/tracks.pdf");
return 0;
int MagnetPlots() {
TRestAxionMagneticField* field = new TRestAxionMagneticField("fields.rml", "babyIAXO_2024_cutoff");
TCanvas* c = field->DrawTracks(TVector3(0, 0, 8000), 100);
c->Print("/tmp/tracks.pdf");
return 0;
}

int OpticsPlots()
{
TCanvas *c1 = new TCanvas("c1", "My canvas", 800, 350);
TRestAxionTrueWolterOptics *optics = new TRestAxionTrueWolterOptics("xmmTrueWolter.rml", "xmm" );
TPad *pad = optics->DrawDensityMaps(7000, 4.2, 0.0, samples, 7500, 9, true);
c1->Print("/tmp/optics.pdf");

TCanvas *c2 = new TCanvas("c2", "My canvas", 600, 800);
TPad *pad2 = optics->DrawParticleTracks(0.001, 40, true);
c2->Print("/tmp/optics2.pdf");

return 0;
int OpticsPlots() {
TCanvas* c1 = new TCanvas("c1", "My canvas", 800, 350);
TRestAxionTrueWolterOptics* optics = new TRestAxionTrueWolterOptics("xmmTrueWolter.rml", "xmm");
TPad* pad = optics->DrawDensityMaps(7000, 4.2, 0.0, samples, 7500, 9, true);
c1->Print("/tmp/optics.pdf");

TCanvas* c2 = new TCanvas("c2", "My canvas", 600, 800);
TPad* pad2 = optics->DrawParticleTracks(0.001, 40, true);
c2->Print("/tmp/optics2.pdf");

return 0;
}

int MirrorPlots()
{
TRestAxionOpticsMirror *mirror = new TRestAxionOpticsMirror("mirrors.rml", "default");
TCanvas *c = mirror->DrawOpticsProperties("[4,7,10](0,1){0.73,0.79,0.9,0.96}:[0.25,0.5,0.75,1](0,10){0.73,0.79,0.9,0.96}", 1.e-5, 1.e-3, true);
c->Print("/tmp/mirrors.pdf");
return 0;
int MirrorPlots() {
TRestAxionOpticsMirror* mirror = new TRestAxionOpticsMirror("mirrors.rml", "default");
TCanvas* c = mirror->DrawOpticsProperties(
"[4,7,10](0,1){0.73,0.79,0.9,0.96}:[0.25,0.5,0.75,1](0,10){0.73,0.79,0.9,0.96}", 1.e-5, 1.e-3, true);
c->Print("/tmp/mirrors.pdf");
return 0;
}

int REST_Axion_PaperPlots(std::string plotsString = "magnetPlots" )
{
if( plotsString == "magnetPlots" )
MagnetPlots();
if( plotsString == "opticsPlots" )
OpticsPlots();
if( plotsString == "mirrorPlots" )
MirrorPlots();
return 0;
int REST_Axion_PaperPlots(std::string plotsString = "magnetPlots") {
if (plotsString == "magnetPlots") MagnetPlots();
if (plotsString == "opticsPlots") OpticsPlots();
if (plotsString == "mirrorPlots") MirrorPlots();
return 0;
}
Loading

0 comments on commit 8032860

Please sign in to comment.