Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Remove Parabolic and Hyperbolic mirror restrictions and add x_sep #81

Merged
merged 41 commits into from
Feb 23, 2024
Merged

Remove Parabolic and Hyperbolic mirror restrictions and add x_sep #81

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
41 commits
Select commit Hold shift + click to select a range
9815117
Adding x_sep in TRestAxionTrueWolterOptics.h
jovoy Aug 22, 2023
dbcc47c
Adding x_sep to TRestAxionWolterOptics.h
jovoy Aug 22, 2023
92fda9d
Adding x_sep to TRestAxionTrueWolterOptics.cxx
jovoy Aug 22, 2023
350dd66
[pre-commit.ci] auto fixes from pre-commit.com hooks
pre-commit-ci[bot] Aug 22, 2023
2aefce7
Adding x_sep to TRestAxionWolterOptics.cxx
jovoy Aug 22, 2023
941945d
Correct variable spelling
jovoy Aug 22, 2023
7deb7d9
Correcting variable spelling
jovoy Aug 22, 2023
e6ed715
Adding x_sep for hyperbolic and parabolic mirror interaction
jovoy Aug 22, 2023
5077ff5
[pre-commit.ci] auto fixes from pre-commit.com hooks
pre-commit-ci[bot] Aug 22, 2023
dc5541f
Adding x_sep to entrance and exit position in TRestAxionTrueWolterOpt…
jovoy Aug 22, 2023
ad5b10d
Fixed x_sep into a vector in TRestAxionTrueWolterOptics.cxx
jovoy Aug 22, 2023
a4dc161
[pre-commit.ci] auto fixes from pre-commit.com hooks
pre-commit-ci[bot] Aug 22, 2023
c69bc30
Adding x_sep to entrance and exit position in TRestAxionWolterOptics.h
jovoy Aug 22, 2023
b967ccb
Correcting x_sep into a vector in TRestAxionWolterOptics.cxx
jovoy Aug 22, 2023
f44569d
[pre-commit.ci] auto fixes from pre-commit.com hooks
pre-commit-ci[bot] Aug 22, 2023
4540f6c
Added x_sep in first and second mirror interaction in TRestAxionWolte…
jovoy Aug 22, 2023
40ca3be
[pre-commit.ci] auto fixes from pre-commit.com hooks
pre-commit-ci[bot] Aug 22, 2023
c9b4362
Removed variables from TRestPhysics
jovoy Aug 28, 2023
923a11f
Corrected a sign switch
jovoy Aug 28, 2023
04f2cee
Correction of vectors in definition of x_sep
jovoy Aug 28, 2023
6adc342
[pre-commit.ci] auto fixes from pre-commit.com hooks
pre-commit-ci[bot] Aug 28, 2023
30aff7f
X-sep vector correction
jovoy Aug 28, 2023
8b00dea
[pre-commit.ci] auto fixes from pre-commit.com hooks
pre-commit-ci[bot] Aug 28, 2023
d45f68f
README.md dummy change to test pipelines
jgalan Aug 29, 2023
6c11abb
Merge branch 'master' into jovoy-mirror_fix
jovoy Dec 15, 2023
ae3d5a9
Merge branch 'master' into jovoy-mirror_fix
jgalan Dec 15, 2023
b2eaf97
Merge branch 'master' into jovoy-mirror_fix
jovoy Feb 20, 2024
320b038
Comment update
jovoy Feb 20, 2024
0015cc3
Add R2 and R4
jovoy Feb 21, 2024
a5836d0
R2 and R4 in cone shaped optics
jovoy Feb 21, 2024
3d2ff54
Have x_sep be 0 for R2 equals R3
jovoy Feb 21, 2024
cc7e1e4
[pre-commit.ci] auto fixes from pre-commit.com hooks
pre-commit-ci[bot] Feb 21, 2024
a402d7f
Add x_sep equals 0 in cone shaped
jovoy Feb 21, 2024
bfb5535
[pre-commit.ci] auto fixes from pre-commit.com hooks
pre-commit-ci[bot] Feb 21, 2024
f971d3b
Update x_sep description true Wolter
jovoy Feb 22, 2024
2799fde
[pre-commit.ci] auto fixes from pre-commit.com hooks
pre-commit-ci[bot] Feb 22, 2024
d46604e
Update Wolter description of x_sep
jovoy Feb 22, 2024
701c676
[pre-commit.ci] auto fixes from pre-commit.com hooks
pre-commit-ci[bot] Feb 22, 2024
d688e78
Change alpha to beta in hyperbolic fuctions
jovoy Feb 22, 2024
60f5503
Update commentation
jovoy Feb 22, 2024
ccc31a6
[pre-commit.ci] auto fixes from pre-commit.com hooks
pre-commit-ci[bot] Feb 22, 2024
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
1 change: 1 addition & 0 deletions README.md
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -65,6 +65,7 @@ It is necessary to enable the `REST_MPFR` cmake option in order to benefit from
### Publications

This repository makes use of the following published codes:

- K. Altenmuller et al, REST-for-Physics, a ROOT-based framework for event oriented data analysis and combined Monte Carlo response, [Computer Physics Communications 273, April 2022, 108281](https://doi.org/10.1016/j.cpc.2021.108281).
- S.Hoof, J.Jaeckel, T.J.Lennert, Quantifying uncertainties in the solar axion flux and their impact on determining axion model parameters, [JCAP09(2021)006](https://doi.org/10.1088/1475-7516/2021/09/006).
- T.Kittelmann, E.Klinkby, E.B.Knudsen, P.Willendrup, X.X.Cai, K.Kanaki, Monte Carlo Particle Lists: MCPL, Computer Physics Communications 218 (2017) 17–42](http://dx.doi.org/10.17632/cby92vsv5g.1).
27 changes: 25 additions & 2 deletions inc/TRestAxionTrueWolterOptics.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -35,9 +35,15 @@ class TRestAxionTrueWolterOptics : public TRestAxionOptics {
/// Entrance radius R1 in mm. See schematic figure.
std::vector<Double_t> fR1; //!

/// Radius R2 in mm. See schematic figure.
std::vector<Double_t> fR2; //!

/// Radius R3 in mm. See schematic figure.
std::vector<Double_t> fR3; //!

/// Radius R4 in mm. See schematic figure.
std::vector<Double_t> fR4; //!

/// Radius R5 in mm. See schematic figure.
std::vector<Double_t> fR5; //!

Expand All @@ -47,6 +53,11 @@ class TRestAxionTrueWolterOptics : public TRestAxionOptics {
/// Mirror thickness in mm. See schematic figure.
std::vector<Double_t> fThickness; //!

/// Distance between mirror stacks in mm. See schematic figure.
/// This is here calculated using the functions from
/// https://backend.orbit.dtu.dk/ws/portalfiles/portal/122353510/phdthesis_for_DTU_orbit.pdf.
std::vector<Double_t> fXSep; //!

/// The spider structure to be used as an optical opaque mask (common to all planes)
TRestSpiderMask* fSpiderMask = nullptr; //<

Expand Down Expand Up @@ -83,12 +94,24 @@ class TRestAxionTrueWolterOptics : public TRestAxionOptics {
return r;
}

/// It returns a vector with the values of R2
std::vector<Double_t> GetR2() {
std::vector<Double_t> r = TRestTools::GetColumnFromTable(fOpticsData, 1);
return r;
}

/// It returns a vector with the values of R3
std::vector<Double_t> GetR3() {
std::vector<Double_t> r = TRestTools::GetColumnFromTable(fOpticsData, 2);
return r;
}

/// It returns a vector with the values of R4
std::vector<Double_t> GetR4() {
std::vector<Double_t> r = TRestTools::GetColumnFromTable(fOpticsData, 3);
return r;
}

/// It returns a vector with the values of R5
std::vector<Double_t> GetR5() {
std::vector<Double_t> r = TRestTools::GetColumnFromTable(fOpticsData, 4);
Expand All @@ -110,13 +133,13 @@ class TRestAxionTrueWolterOptics : public TRestAxionOptics {

/// It returns the entrance Z-position defined by the optical axis.
Double_t GetEntrancePositionZ() override {
if (fCosAlpha.size() > 0) return -fMirrorLength * fCosAlpha[0];
if (fCosAlpha.size() > 0) return -fMirrorLength * fCosAlpha[0] - 0.5 * fXSep[0];
return 0;
}

/// It returns the exit Z-position defined by the optical axis
Double_t GetExitPositionZ() override {
if (fCosAlpha_3.size() > 0) return fMirrorLength * fCosAlpha_3[0];
if (fCosAlpha_3.size() > 0) return fMirrorLength * fCosAlpha_3[0] + 0.5 * fXSep[0];
return 0;
}

Expand Down
27 changes: 25 additions & 2 deletions inc/TRestAxionWolterOptics.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -35,9 +35,15 @@ class TRestAxionWolterOptics : public TRestAxionOptics {
/// Entrance radius R1 in mm. See schematic figure.
std::vector<Double_t> fR1; //!

/// Radius R2 in mm. See schematic figure.
std::vector<Double_t> fR2; //!

/// Radius R3 in mm. See schematic figure.
std::vector<Double_t> fR3; //!

/// Radius R4 in mm. See schematic figure.
std::vector<Double_t> fR4; //!

/// Radius R5 in mm. See schematic figure.
std::vector<Double_t> fR5; //!

Expand All @@ -47,6 +53,11 @@ class TRestAxionWolterOptics : public TRestAxionOptics {
/// Mirror thickness in mm. See schematic figure.
std::vector<Double_t> fThickness; //!

/// Distance between mirror stacks in mm. See schematic figure.
/// This is here calculated using the functions from
/// https://backend.orbit.dtu.dk/ws/portalfiles/portal/122353510/phdthesis_for_DTU_orbit.pdf.
std::vector<Double_t> fXSep; //!

/// The spider structure to be used as an optical opaque mask (common to all planes)
TRestSpiderMask* fSpiderMask = nullptr; //<

Expand Down Expand Up @@ -83,12 +94,24 @@ class TRestAxionWolterOptics : public TRestAxionOptics {
return r;
}

/// It returns a vector with the values of R2
std::vector<Double_t> GetR2() {
std::vector<Double_t> r = TRestTools::GetColumnFromTable(fOpticsData, 1);
return r;
}

/// It returns a vector with the values of R3
std::vector<Double_t> GetR3() {
std::vector<Double_t> r = TRestTools::GetColumnFromTable(fOpticsData, 2);
return r;
}

/// It returns a vector with the values of R4
std::vector<Double_t> GetR4() {
std::vector<Double_t> r = TRestTools::GetColumnFromTable(fOpticsData, 3);
return r;
}

/// It returns a vector with the values of R5
std::vector<Double_t> GetR5() {
std::vector<Double_t> r = TRestTools::GetColumnFromTable(fOpticsData, 4);
Expand All @@ -110,13 +133,13 @@ class TRestAxionWolterOptics : public TRestAxionOptics {

/// It returns the entrance Z-position defined by the optical axis.
Double_t GetEntrancePositionZ() override {
if (fCosAlpha.size() > 0) return -fMirrorLength * fCosAlpha[0];
if (fCosAlpha.size() > 0) return -fMirrorLength * fCosAlpha[0] - 0.5 * fXSep[0];
return 0;
}

/// It returns the exit Z-position defined by the optical axis
Double_t GetExitPositionZ() override {
if (fCosAlpha_3.size() > 0) return fMirrorLength * fCosAlpha_3[0];
if (fCosAlpha_3.size() > 0) return fMirrorLength * fCosAlpha_3[0] + 0.5 * fXSep[0];
return 0;
}

Expand Down
29 changes: 21 additions & 8 deletions src/TRestAxionTrueWolterOptics.cxx
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -169,11 +169,22 @@ void TRestAxionTrueWolterOptics::Initialize() {
SetLibraryVersion(LIBRARY_VERSION);

fR1 = GetR1();
fR2 = GetR2();
fR3 = GetR3();
fR4 = GetR4();
fR5 = GetR5();
fAlpha = GetAlpha();
fThickness = GetThickness();

fXSep.clear();
for (unsigned int n = 0; n < fAlpha.size(); n++)
if (fR2[n] == fR3[n]) {
fXSep.push_back(0);
} else {
fXSep.push_back(2 * (fR1[n] - fR3[n] - fMirrorLength * TMath::Sin(fAlpha[n])) /
TMath::Tan(fAlpha[n]));
}

if (fAlpha.size() == 0) return;

fCosAlpha.clear();
Expand Down Expand Up @@ -275,10 +286,11 @@ Int_t TRestAxionTrueWolterOptics::FirstMirrorReflection(const TVector3& pos, con
TVector3 vertex(0, 0, fFrontVertex[mirror]);

//// Reflection on first mirror
fFirstInteractionPosition = REST_Physics::GetParabolicVectorIntersection(
pos, dir, fAlpha[mirror], fR3[mirror], fMirrorLength); // should add this: TVector3(0, 0, -1), vertex
fFirstInteractionPosition =
REST_Physics::GetParabolicVectorIntersection(pos, dir, fAlpha[mirror], fR3[mirror]);

if (fFirstInteractionPosition.Z() < GetEntrancePositionZ() || fFirstInteractionPosition.Z() > 0) {
if (fFirstInteractionPosition.Z() < GetEntrancePositionZ() ||
fFirstInteractionPosition.Z() > -(0.5 * fXSep[mirror])) {
RESTDebug << "TRestAxionTrueWolterOptics::FirstMirrorReflection. No interaction!" << RESTendl;
fFirstInteractionPosition = REST_Physics::MoveByDistance(pos, dir, fMirrorLength / 2.);
fMiddleDirection = fEntranceDirection;
Expand Down Expand Up @@ -320,13 +332,14 @@ Int_t TRestAxionTrueWolterOptics::SecondMirrorReflection(const TVector3& pos, co

TVector3 vertex(0, 0, fBackVertex[mirror]);
Double_t focal = fR3[mirror] / TMath::Tan(4 * fAlpha[mirror]);
Copy link
Member

@Vindaar Vindaar Feb 22, 2024 •
edited
Loading

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Is fR3 the middle point along the extensions of mirrors 1 and 2 /or/ the virtual reflection point? If it is the former it's not the right equation, if focal is supposed to be the focal length of the telescope.

(this is how you end up with a "1500mm" LLNL telescope, which actually has a 1530mm focal length)

Or is this equation for some reason correct as is for a true Wolter optic?

edit:
Given that the docstrings for fR3 reference a schematic, which I assume is this one from here:

https://rest-for-physics.github.io/framework/classTRestAxionTrueWolterOptics.html

it seems to be the former still.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

So in the document for the true Wolter optics they call this the z-Koordinate of the detector/ focal point. I thought this would be roughly the same for each layer. What do you mean by "virtual reflection point"?

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

So for these true Wolter functions it is not the focal length of the whole optics but of the the single hyperbolic mirror as far as I understand.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

llnl_layers_explanation.pdf

The virtual point is indicated in this schematic I made for my thesis.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

llnl_layers_explanation

here as an SVG so that it's visible embedded.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

spiegel_parameter.pdf

Here is the true Wolter calculation I got from Vadim's colleague and it does state that what I call "focal" is calculated with R3. Interesting point though to note about the virtual reflection point!

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Ok, good to know. Then I suppose this is indeed correct for a true Wolter optic.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

But it's not the whole optics focal spot, that is true and confusing^^

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

But the PDF you sent explicitly says 'Fokallänge' for x_0 - x_f?

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Yeah, but the calculation is right there and they clearly mean R3 with r0

Double_t beta = 3 * fAlpha[mirror];

//// Reflection on first mirror
//// Reflection on second mirror
fSecondInteractionPosition =
REST_Physics::GetHyperbolicVectorIntersection(pos, dir, fAlpha[mirror], fR3[mirror], fMirrorLength,
focal); // should add this: TVector3(0, 0, -1), vertex,
REST_Physics::GetHyperbolicVectorIntersection(pos, dir, beta, fR3[mirror], focal);

if (fSecondInteractionPosition.Z() > GetExitPositionZ() || fSecondInteractionPosition.Z() < 0) {
if (fSecondInteractionPosition.Z() > GetExitPositionZ() ||
fSecondInteractionPosition.Z() < (0.5 * fXSep[mirror])) {
RESTDebug << "TRestAxionTrueWolterOptics::SecondMirrorReflection. No interaction!" << RESTendl;
fSecondInteractionPosition = REST_Physics::MoveByDistance(pos, dir, fMirrorLength / 2.);
fExitDirection = fMiddleDirection;
Expand All @@ -335,7 +348,7 @@ Int_t TRestAxionTrueWolterOptics::SecondMirrorReflection(const TVector3& pos, co
}

TVector3 hyperNormal =
REST_Physics::GetHyperbolicNormal(fSecondInteractionPosition, fAlpha[mirror], fR3[mirror], focal);
REST_Physics::GetHyperbolicNormal(fSecondInteractionPosition, beta, fR3[mirror], focal);

fExitDirection = GetVectorReflection(fMiddleDirection, hyperNormal);

Expand Down
19 changes: 16 additions & 3 deletions src/TRestAxionWolterOptics.cxx
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -169,11 +169,22 @@ void TRestAxionWolterOptics::Initialize() {
SetLibraryVersion(LIBRARY_VERSION);

fR1 = GetR1();
fR2 = GetR2();
fR3 = GetR3();
fR4 = GetR4();
fR5 = GetR5();
fAlpha = GetAlpha();
fThickness = GetThickness();

fXSep.clear();
for (unsigned int n = 0; n < fAlpha.size(); n++)
if (fR2[n] == fR3[n]) {
fXSep.push_back(0);
} else {
fXSep.push_back(2 * (fR1[n] - fR3[n] - fMirrorLength * TMath::Sin(fAlpha[n])) /
TMath::Tan(fAlpha[n]));
}

if (fAlpha.size() == 0) return;

fCosAlpha.clear();
Expand Down Expand Up @@ -278,7 +289,8 @@ Int_t TRestAxionWolterOptics::FirstMirrorReflection(const TVector3& pos, const T
fFirstInteractionPosition =
pos + dir * REST_Physics::GetConeVectorIntersection(pos, dir, TVector3(0, 0, -1), vertex, cosA);

if (fFirstInteractionPosition.Z() < GetEntrancePositionZ() || fFirstInteractionPosition.Z() > 0) {
if (fFirstInteractionPosition.Z() < GetEntrancePositionZ() ||
fFirstInteractionPosition.Z() > -(0.5 * fXSep[mirror])) {
RESTDebug << "TRestAxionWolterOptics::FirstMirrorReflection. No interaction!" << RESTendl;
fFirstInteractionPosition = REST_Physics::MoveByDistance(pos, dir, fMirrorLength / 2.);
fMiddleDirection = fEntranceDirection;
Expand Down Expand Up @@ -319,11 +331,12 @@ Int_t TRestAxionWolterOptics::SecondMirrorReflection(const TVector3& pos, const
TVector3 vertex(0, 0, fBackVertex[mirror]);
Double_t cosA = fCosAlpha_3[mirror];

//// Reflection on first mirror
//// Reflection on second mirror
fSecondInteractionPosition =
pos + dir * REST_Physics::GetConeVectorIntersection(pos, dir, TVector3(0, 0, -1), vertex, cosA);

if (fSecondInteractionPosition.Z() > GetExitPositionZ() || fSecondInteractionPosition.Z() < 0) {
if (fSecondInteractionPosition.Z() > GetExitPositionZ() ||
fSecondInteractionPosition.Z() < (0.5 * fXSep[mirror])) {
RESTDebug << "TRestAxionWolterOptics::SecondMirrorReflection. No interaction!" << RESTendl;
fSecondInteractionPosition = REST_Physics::MoveByDistance(pos, dir, fMirrorLength / 2.);
fExitDirection = fMiddleDirection;
Expand Down
Loading