TRestAxionOpticsMirror.cxx

/******************** REST disclaimer ***********************************
 * This file is part of the REST software framework.                     *
 *                                                                       *
 * Copyright (C) 2016 GIFNA/TREX (University of Zaragoza)                *
 * For more information see http://gifna.unizar.es/trex                  *
 *                                                                       *
 * REST is free software: you can redistribute it and/or modify          *
 * it under the terms of the GNU General Public License as published by  *
 * the Free Software Foundation, either version 3 of the License, or     *
 * (at your option) any later version.                                   *
 *                                                                       *
 * REST is distributed in the hope that it will be useful,               *
 * but WITHOUT ANY WARRANTY; without even the implied warranty of        *
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the          *
 * GNU General Public License for more details.                          *
 *                                                                       *
 * You should have a copy of the GNU General Public License along with   *
 * REST in $REST_PATH/LICENSE.                                           *
 * If not, see http://www.gnu.org/licenses/.                             *
 * For the list of contributors see $REST_PATH/CREDITS.                  *
 *************************************************************************/
 
//
 
#include "TRestAxionOpticsMirror.h"
 
#include <TAxis.h>
#include <TGraph.h>
#include <TH1F.h>
#include <TLatex.h>
#include <TLegend.h>
 
using namespace std;
 
ClassImp(TRestAxionOpticsMirror);
 
TRestAxionOpticsMirror::TRestAxionOpticsMirror() : TRestMetadata() { Initialize(); }
 
TRestAxionOpticsMirror::TRestAxionOpticsMirror(const char* cfgFileName, string name)
    : TRestMetadata(cfgFileName) {
    RESTDebug << "Entering TRestAxionOpticsMirror constructor( cfgFileName, name )" << RESTendl;
 
    Initialize();
 
    LoadConfigFromFile(fConfigFileName, name);
 
    if (GetVerboseLevel() >= TRestStringOutput::REST_Verbose_Level::REST_Info) PrintMetadata();
}
 
TRestAxionOpticsMirror::~TRestAxionOpticsMirror() {}
 
void TRestAxionOpticsMirror::Initialize() {
    SetSectionName(this->ClassName());
    SetLibraryVersion(LIBRARY_VERSION);
 
    fHenkeKeys.clear();
 
    if (fMirrorType == "Single") {
        fHenkeKeys["Ldensity"] = "-1";
        fHenkeKeys["Layer"] = fLayerTop;
        fHenkeKeys["Thick"] = fLayerThicknessTop;
        fHenkeKeys["Sigma1"] = fSigmaTop;
        fHenkeKeys["Sigma2"] = "0";
    }
 
    if (fMirrorType == "Bilayer") {
        fHenkeKeys["Tdensity"] = "-1";
        fHenkeKeys["Tlayer"] = fLayerTop;
        fHenkeKeys["Thickt"] = fLayerThicknessTop;
        fHenkeKeys["Sigmat"] = fSigmaTop;
        fHenkeKeys["Bdensity"] = "-1";
        fHenkeKeys["Blayer"] = fLayerBottom;
        fHenkeKeys["Thickb"] = fLayerThicknessBottom;
        fHenkeKeys["Sigmab"] = fSigmaBottom;
        fHenkeKeys["Sigmas"] = "0";
    }
 
    fHenkeKeys["Substrate"] = fSubstrate;
    fHenkeKeys["Sdensity"] = "-1";
    fHenkeKeys["Pol"] = "0";
    fHenkeKeys["Scan"] = "Angle";
    fHenkeKeys["Min"] = "0";
    fHenkeKeys["Max"] = "4.5";
    fHenkeKeys["Npts"] = "450";
    fHenkeKeys["temp"] = "Energy+%28eV%29";
    fHenkeKeys["Fixed"] = "200";
    fHenkeKeys["Plot"] = "LinLog";
    fHenkeKeys["Output"] = "Plot";
}
 
void TRestAxionOpticsMirror::LoadTables() {
    Initialize();
 
    string mirrorFile = SearchFile(GetReflectivityFilename());
    string windowFile = SearchFile(GetTransmissionFilename());
 
    if (mirrorFile != "" && windowFile != "") {
        TRestTools::ReadBinaryTable(mirrorFile, fReflectivityTable, 901);
        TRestTools::ReadBinaryTable(windowFile, fTransmissionTable, 901);
        return;
    }
 
    cout << "--------------------------------- INFO ------------------------------" << endl;
    cout << "The optics material properties were not found in the database" << endl;
    cout << "Downloading from Henke database. This process will take a long time." << endl;
    cout << "After producing these files, please, contribute them to the " << endl;
    cout << "rest-for-physics/axionlib-data inside the optics directory ..." << endl;
    cout << "---------------------------------- INFO ------------------------------" << endl;
 
    fReflectivityTable.clear();
    fTransmissionTable.clear();
    map<Int_t, std::vector<Float_t>> reflectivity;
    map<Int_t, std::vector<Float_t>> transmission;
    for (double n = 0; n < 9; n = n + 4.5) {
        fHenkeKeys["Min"] = DoubleToString(n);
        fHenkeKeys["Max"] = DoubleToString(n + 4.5);
 
        cout.flush();
        vector<Float_t> reflect;
        vector<Float_t> transm;
        reflect.clear();
        transm.clear();
 
        cout << "Scanning angles between " << n << " and " << n + 4.5 << ".";
        for (int e = 30; e <= 15000; e += 30) {
            fHenkeKeys["Fixed"] = IntegerToString(e);
            cout << ".";
            cout.flush();
            string fname = DownloadHenkeFile();
 
            std::vector<std::vector<Double_t>> data;
            if (!TRestTools::ReadASCIITable(fname, data, 2)) {
                RESTError << "TRestAxionOpticsMirror. Error reading HenkeFile table" << RESTendl;
            }
 
            // we skip the last point if we are not at the latest angles file
            Int_t N = data.size() - 1;
            if (n == 4.5) N = N + 1;
 
            for (int m = 0; m < N; m++) reflectivity[e].push_back((Float_t)data[m][1]);
            for (int m = 0; m < N; m++) transmission[e].push_back((Float_t)data[m][2]);
        }
        cout << endl;
    }
    cout << "Number of angles stored : " << reflectivity[30].size() << endl;
 
    for (const auto& x : reflectivity) {
        fReflectivityTable.push_back(x.second);
    }
    for (const auto& x : transmission) {
        fTransmissionTable.push_back(x.second);
    }
 
    ExportTables();
}
 
std::string TRestAxionOpticsMirror::GetReflectivityFilename() {
    string fnameR = "Reflectivity_" + fMirrorType + "_" + fLayerTop + "_" + fLayerThicknessTop + "_" +
                    fSubstrate + "_" + fSigmaTop + ".N901f";
    if (fMirrorType == "Bilayer")
        fnameR = "Reflectivity_" + fMirrorType + "_" + fLayerTop + "_" + fLayerThicknessTop + "_" +
                 fSigmaTop + "_" + fLayerBottom + "_" + fLayerThicknessBottom + "_" + fSigmaBottom + "_" +
                 fSubstrate + ".N901f";
    return fnameR;
}
 
std::string TRestAxionOpticsMirror::GetTransmissionFilename() {
    string fnameT = "Transmission_" + fMirrorType + "_" + fLayerTop + "_" + fLayerThicknessTop + "_" +
                    fSubstrate + "_" + fSigmaTop + ".N901f";
    if (fMirrorType == "Bilayer")
        fnameT = "Transmission_" + fMirrorType + "_" + fLayerTop + "_" + fLayerThicknessTop + "_" +
                 fSigmaTop + "_" + fLayerBottom + "_" + fLayerThicknessBottom + "_" + fSigmaBottom + "_" +
                 fSubstrate + ".N901f";
    return fnameT;
}
 
Int_t TRestAxionOpticsMirror::ExportTables() {
    if (fReflectivityTable.size() == 0) {
        RESTError << "Nothing to export!" << RESTendl;
        return 1;
    }
 
    string path = REST_USER_PATH + "/export/";
 
    if (!TRestTools::fileExists(path)) {
        std::cout << "Creating path: " << path << std::endl;
        int z = system(("mkdir -p " + path).c_str());
        if (z != 0) RESTError << "Problem creating directory: " << path << RESTendl;
    }
 
    string fnameR = GetReflectivityFilename();
    TRestTools::ExportBinaryTable(path + fnameR, fReflectivityTable);
    RESTInfo << "Reflectivity table generated at: " << path + fnameR << RESTendl;
 
    string fnameT = GetTransmissionFilename();
    TRestTools::ExportBinaryTable(path + fnameT, fTransmissionTable);
    RESTInfo << "Transmission table generated at: " << path + fnameT << RESTendl;
 
    return 0;
}
 
std::string TRestAxionOpticsMirror::DownloadHenkeFile() {
    string perlName = "laymir.pl";
    if (fMirrorType == "Bilayer") perlName = "bimir.pl";
    string url = "https://henke.lbl.gov/cgi-bin/" + perlName;
    string result = TRestTools::POSTRequest(url, fHenkeKeys);
 
    size_t start = result.find("HREF=\"") + 6;
    size_t length = result.find(".dat\">") + 4 - start;
    result = result.substr(start, length);
 
    return TRestTools::DownloadRemoteFile("https://henke.lbl.gov/" + result);
}
 
Double_t TRestAxionOpticsMirror::GetReflectivity(const Double_t angle, const Double_t energy) {
    if (fReflectivityTable.size() == 0) LoadTables();
 
    Double_t en = energy;
    if (en < 0.030) {
        RESTWarning << "Energy is below 30eV! It should be between 30eV and 15keV" << RESTendl;
        RESTWarning << "Setting energy to 30eV" << RESTendl;
        en = 0.030;
    }
 
    if (en > 15) {
        RESTWarning << "Energy is above 15keV! It should be between 30eV and 15keV" << RESTendl;
        RESTWarning << "Setting energy to 15keV" << RESTendl;
        en = 15;
    }
 
    Int_t lowEnBin = (Int_t)((en - 0.03) / 0.03);
    Double_t deltaE = (en - (Double_t)(lowEnBin + 1) * 0.03) / 0.03;  // between 0 and 1
 
    Double_t ang = angle;
    if (ang < 0.0) {
        RESTWarning << "Angle is below 0 degrees! It should be between 0 and 9 degrees" << RESTendl;
        RESTWarning << "Setting angle to 0 degrees" << RESTendl;
        ang = 0.0;
    }
 
    if (ang > 9) {
        RESTWarning << "Angle is above 9 degrees! It should be between 0 and 9 degrees" << RESTendl;
        RESTWarning << "Setting angle to 9 degrees" << RESTendl;
        ang = 9;
    }
 
    Int_t lowAngBin = (Int_t)((ang) / 0.01);
    Double_t deltaAng = (ang - (Double_t)(lowAngBin)*0.01) / 0.01;  // between 0 and 1
 
    Double_t REnLowAngLow = fReflectivityTable[lowEnBin][lowAngBin];
    Double_t REnLowAngHi = fReflectivityTable[lowEnBin][lowAngBin + 1];
    Double_t REnHiAngLow = fReflectivityTable[lowEnBin + 1][lowAngBin];
    Double_t REnHiAngHi = fReflectivityTable[lowEnBin + 1][lowAngBin + 1];
 
    // We have renormalized the grid to unity and now we apply the equation z = f(x,y)
    // where x is associated with the energy, y is associated with the angle
    // z = f(x,y) = (1−x)(1−y) * v00+x(1−y) * v10+(1−x)y * v01+xy * v11
    // So that, for example, when x=1 and v=1 we get v11=REnHiAngHi
    return (1 - deltaE) * (1 - deltaAng) * REnLowAngLow + deltaE * (1 - deltaAng) * REnHiAngLow +
           (1 - deltaE) * deltaAng * REnLowAngHi + deltaE * deltaAng * REnHiAngHi;
}
 
Double_t TRestAxionOpticsMirror::GetTransmission(const Double_t angle, const Double_t energy) {
    if (fTransmissionTable.size() == 0) LoadTables();
 
    Double_t en = energy;
    if (en < 0.030) {
        RESTWarning << "Energy is below 30eV! It should be between 30eV and 15keV" << RESTendl;
        RESTWarning << "Setting energy to 30eV" << RESTendl;
        en = 0.030;
    }
 
    if (en > 15) {
        RESTWarning << "Energy is above 15keV! It should be between 30eV and 15keV" << RESTendl;
        RESTWarning << "Setting energy to 15keV" << RESTendl;
        en = 15;
    }
 
    Int_t lowEnBin = (Int_t)((en - 0.03) / 0.03);
    Double_t deltaE = (en - (Double_t)(lowEnBin + 1) * 0.03) / 0.03;  // between 0 and 1
 
    Double_t ang = angle;
    if (ang < 0.0) {
        RESTWarning << "Angle is below 0 degrees! It should be between 0 and 9 degrees" << RESTendl;
        RESTWarning << "Setting angle to 0 degrees" << RESTendl;
        ang = 0.0;
    }
 
    if (ang > 9) {
        RESTWarning << "Angle is above 9 degrees! It should be between 0 and 9 degrees" << RESTendl;
        RESTWarning << "Setting angle to 9 degrees" << RESTendl;
        ang = 9;
    }
 
    Int_t lowAngBin = (Int_t)((ang) / 0.01);
    Double_t deltaAng = (ang - (Double_t)(lowAngBin)*0.01) / 0.01;  // between 0 and 1
 
    Double_t REnLowAngLow = fTransmissionTable[lowEnBin][lowAngBin];
    Double_t REnLowAngHi = fTransmissionTable[lowEnBin][lowAngBin + 1];
    Double_t REnHiAngLow = fTransmissionTable[lowEnBin + 1][lowAngBin];
    Double_t REnHiAngHi = fTransmissionTable[lowEnBin + 1][lowAngBin + 1];
 
    // We have renormalized the grid to unity and now we apply the equation z = f(x,y)
    // where x is associated with the energy, y is associated with the angle
    // z = f(x,y) = (1−x)(1−y) * v00+𝑥(1−y) * v10+(1−𝑥)y * v01+xy * v11
    // So that, for example, when x=1 and v=1 we get v11=REnHiAngHi
 
    return (1 - deltaE) * (1 - deltaAng) * REnLowAngLow + deltaE * (1 - deltaAng) * REnHiAngLow +
           (1 - deltaE) * deltaAng * REnLowAngHi + deltaE * deltaAng * REnHiAngHi;
}
 
void TRestAxionOpticsMirror::PrintMetadata() {
    TRestMetadata::PrintMetadata();
 
    RESTMetadata << "Mirror type: " << fMirrorType << RESTendl;
    if (fMirrorType == "Single") {
        RESTMetadata << "Layer material: " << fLayerTop << RESTendl;
        RESTMetadata << "Layer thickness: " << fLayerThicknessTop << " nm" << RESTendl;
        RESTMetadata << "Layer roughness: " << fSigmaTop << "nm" << RESTendl;
    }
 
    if (fMirrorType == "Bilayer") {
        RESTMetadata << "Top layer material: " << fLayerTop << RESTendl;
        RESTMetadata << "Top layer thickness: " << fLayerThicknessTop << " nm" << RESTendl;
        RESTMetadata << "Top layer roughness: " << fSigmaTop << "nm" << RESTendl;
        RESTMetadata << "Bottom layer material: " << fLayerBottom << RESTendl;
        RESTMetadata << "Bottom layer thickness: " << fLayerThicknessBottom << " nm" << RESTendl;
        RESTMetadata << "Bottom layer roughness: " << fSigmaBottom << "nm" << RESTendl;
    }
    RESTMetadata << "Substrate material: " << fSubstrate << RESTendl;
    RESTMetadata << "+++++++++++++++++++++++++++++++++++++++++++++++++" << RESTendl;
}
 
TCanvas* TRestAxionOpticsMirror::DrawOpticsProperties(std::string options, Double_t lowRange,
                                                      Double_t lowRange2) {
    if (fReflectivityTable.size() == 0) LoadTables();
 
    std::vector<string> optList = TRestTools::GetOptions(options);
 
    if (optList.size() == 0)
        optList = TRestTools::GetOptions(
            "[1,4,7,10](0,9){0.6,0.68,0.9,0.88}:[0.25,0.5,0.75,1](0,10){0.2,0.2,0.45,0.45}");
 
    if (optList.size() != 2) {
        RESTError << "TRestAxionOpticsMirror::DrawOpticsProperties. Wrong arguments!" << RESTendl;
        return fCanvas;
    }
 
    std::vector<double> energies = StringToElements(optList[0], "[", ",", "]");
    std::vector<double> aRange = StringToElements(optList[0], "(", ",", ")");
    std::vector<double> eLegendCoords = StringToElements(optList[0], "{", ",", "}");
 
    std::vector<double> angles = StringToElements(optList[1], "[", ",", "]");
    std::vector<double> eRange = StringToElements(optList[1], "(", ",", ")");
    std::vector<double> aLegendCoords = StringToElements(optList[1], "{", ",", "}");
 
    if (eRange[0] < 0.03) eRange[0] = 0.03;
    if (eRange[1] > 15) eRange[1] = 15;
    if (aRange[0] < 0.0) aRange[0] = 0.0;
    if (aRange[1] > 9) aRange[1] = 9;
 
    //   Double_t lowReflec = TRestTools::GetMinValueFromTable(fReflectivityTable);
    //   Double_t highReflec = TRestTools::GetMaxValueFromTable(fReflectivityTable);
 
    if (fCanvas != NULL) {
        delete fCanvas;
        fCanvas = NULL;
    }
    fCanvas = new TCanvas("canv", "This is the canvas title", 1400, 1200);
    fCanvas->Draw();
 
    TPad* pad1 = new TPad("pad1", "This is pad1", 0.01, 0.02, 0.99, 0.97);
    pad1->Divide(2, 2);
    pad1->Draw();
 
    pad1->cd(1);
    pad1->cd(1)->SetLogy();
    pad1->cd(1)->SetRightMargin(0.09);
    pad1->cd(1)->SetLeftMargin(0.15);
    pad1->cd(1)->SetBottomMargin(0.15);
 
    std::vector<TGraph*> ref_vs_ang_graph;
 
    for (unsigned int n = 0; n < energies.size(); n++) {
        string grname = "gr" + IntegerToString(n);
        TGraph* gr = new TGraph();
        gr->SetName(grname.c_str());
        for (double a = aRange[0]; a <= aRange[1]; a += (aRange[1] - aRange[0]) / 10000.) {
            gr->SetPoint(gr->GetN(), a, GetReflectivity(a, energies[n]));
        }
        gr->SetLineColor(49 - n * 3);
        gr->SetLineWidth(2);
        ref_vs_ang_graph.push_back(gr);
    }
 
    ref_vs_ang_graph[0]->GetXaxis()->SetLimits(aRange[0], aRange[1]);
    ref_vs_ang_graph[0]->GetHistogram()->SetMaximum(1);
    ref_vs_ang_graph[0]->GetHistogram()->SetMinimum(lowRange);
 
    ref_vs_ang_graph[0]->GetXaxis()->SetTitle("Angle [degrees]");
    ref_vs_ang_graph[0]->GetXaxis()->SetTitleSize(0.05);
    ref_vs_ang_graph[0]->GetXaxis()->SetLabelSize(0.05);
    ref_vs_ang_graph[0]->GetYaxis()->SetTitle("Reflectivity");
    ref_vs_ang_graph[0]->GetYaxis()->SetTitleOffset(1.5);
    ref_vs_ang_graph[0]->GetYaxis()->SetTitleSize(0.05);
    ref_vs_ang_graph[0]->GetYaxis()->SetLabelSize(0.05);
    pad1->cd(1)->SetLogy();
    ref_vs_ang_graph[0]->Draw("AL");
    for (unsigned int n = 1; n < energies.size(); n++) ref_vs_ang_graph[n]->Draw("L");
 
    Double_t lx1 = 0.6, ly1 = 0.75, lx2 = 0.9, ly2 = 0.95;
    if (eLegendCoords.size() > 0) {
        lx1 = eLegendCoords[0];
        ly1 = eLegendCoords[1];
        lx2 = eLegendCoords[2];
        ly2 = eLegendCoords[3];
    }
    TLegend* legend = new TLegend(lx1, ly1, lx2, ly2);
 
    legend->SetTextSize(0.03);
    legend->SetHeader("Energies", "C");  // option "C" allows to center the header
    for (unsigned int n = 0; n < energies.size(); n++) {
        std::string lname = "gr" + IntegerToString(n);
        std::string ltitle = DoubleToString(energies[n]) + " keV";
 
        legend->AddEntry(lname.c_str(), ltitle.c_str(), "l");
    }
    legend->Draw();
 
    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);
 
    std::vector<TGraph*> ref_vs_en_graph;
 
    for (unsigned int n = 0; n < angles.size(); n++) {
        string grname = "agr" + IntegerToString(n);
        TGraph* gr = new TGraph();
        gr->SetName(grname.c_str());
        for (double e = eRange[0]; e <= eRange[1]; e += (eRange[1] - eRange[0]) / 10000.) {
            gr->SetPoint(gr->GetN(), e, GetReflectivity(angles[n], e));
        }
        gr->SetLineColor(49 - n * 3);
        gr->SetLineWidth(2);
        ref_vs_en_graph.push_back(gr);
    }
 
    ref_vs_en_graph[0]->GetXaxis()->SetLimits(eRange[0], eRange[1]);
    ref_vs_en_graph[0]->GetHistogram()->SetMaximum(1);
    ref_vs_en_graph[0]->GetHistogram()->SetMinimum(lowRange2);
 
    ref_vs_en_graph[0]->GetXaxis()->SetTitle("Energy [keV]");
    ref_vs_en_graph[0]->GetXaxis()->SetTitleSize(0.05);
    ref_vs_en_graph[0]->GetXaxis()->SetLabelSize(0.05);
    ref_vs_en_graph[0]->GetYaxis()->SetTitle("Reflectivity");
    ref_vs_en_graph[0]->GetYaxis()->SetTitleOffset(1.5);
    ref_vs_en_graph[0]->GetYaxis()->SetTitleSize(0.05);
    ref_vs_en_graph[0]->GetYaxis()->SetLabelSize(0.05);
    pad1->cd(2)->SetLogy();
    ref_vs_en_graph[0]->Draw("AL");
    for (unsigned int n = 1; n < angles.size(); n++) ref_vs_en_graph[n]->Draw("L");
 
    if (aLegendCoords.size() > 0) {
        lx1 = aLegendCoords[0];
        ly1 = aLegendCoords[1];
        lx2 = aLegendCoords[2];
        ly2 = aLegendCoords[3];
    }
    TLegend* legendA = new TLegend(lx1, ly1, lx2, ly2);
    legendA->SetTextSize(0.03);
    legendA->SetHeader("Angles", "C");  // option "C" allows to center the header
    for (unsigned int n = 0; n < angles.size(); n++) {
        std::string lname = "agr" + IntegerToString(n);
        std::string ltitle = DoubleToString(angles[n]) + " degrees";
 
        legendA->AddEntry(lname.c_str(), ltitle.c_str(), "l");
    }
    legendA->Draw();
 
 
    pad1->cd(3);
    pad1->cd(3)->SetRightMargin(0.09);
    pad1->cd(3)->SetLeftMargin(0.15);
    pad1->cd(3)->SetBottomMargin(0.15);
    // ref_vs_ang_graph[0]->GetHistogram()->SetMinimum(0);
    ref_vs_ang_graph[0]->Draw("AL");
    for (unsigned int n = 1; n < energies.size(); n++) ref_vs_ang_graph[n]->Draw("L");
 
    pad1->cd(4);
    pad1->cd(4)->SetRightMargin(0.09);
    pad1->cd(4)->SetLeftMargin(0.15);
    pad1->cd(4)->SetBottomMargin(0.15);
 
    // ref_vs_en_graph[0]->GetHistogram()->SetMinimum(0);
    ref_vs_en_graph[0]->Draw("AL");
    for (unsigned int n = 1; n < angles.size(); n++) ref_vs_en_graph[n]->Draw("L");
 
    fCanvas->cd();  // c1 is the TCanvas
    TPad* pad5 = new TPad("all", "all", 0, 0, 1, 1);
    pad5->SetFillStyle(4000);  // transparent
    pad5->Draw();
    pad5->cd();
    TLatex* lat = new TLatex();
 
    std::string title;
    if (fMirrorType == "Single") {
        title = "Mirror type:" + fMirrorType + " Substrate: " + fSubstrate + ". Layer: " + fLayerTop +
                " Thickness: " + fLayerThicknessTop + "nm Roughness: " + fSigmaTop + "nm.";
        lat->SetTextSize(0.02);
    }
    if (fMirrorType == "Bilayer") {
        title = "Mirror type:" + fMirrorType + " Substrate: " + fSubstrate + ". TOP Layer: " + fLayerTop +
                " Thickness: " + fLayerThicknessTop + "nm Roughness: " + fSigmaTop +
                "nm. \nBOTTOM Layer: " + fLayerBottom + " Thickness: " + fLayerThicknessBottom +
                "nm Roughness: " + fSigmaBottom + " nm.";
        lat->SetTextSize(0.015);
    }
    lat->SetTextAlign(22);
    lat->DrawLatexNDC(.5, .95, title.c_str());
    return fCanvas;
}