This example macro produces the correlation between normalized transition strength B(E2)/A expressed in Weisskopf unit and the ratio between the energy of the first 4 \(^+\) and 2 \(^+\) states for all even-even nuclei. The color scale and the shape of the points reflect the quadrupole deformation \(\beta_2\) of the ground state while the size of each point in proportional to A \(^{2/3}\).
#include "TH2F.h"
#include "TStyle.h"
#include "TGraph2D.h"
#include "TGraph.h"
#include "TLatex.h"
#include "TLegend.h"
#include "TView.h"
#include "TEllipse.h"
#include "TMath.h"
#include "TArrow.h"
#include "TROOT.h"
#include "tkmanager.h"
#include "tknuclear_chart.h"
#define COLORS 20
TH2F* setup_histo(double ts = .05);
void draw_legend(int palette[COLORS], double betamax, double ts = .05);
void draw_deformation()
{
gROOT->SetStyle("tkn-histo");
gStyle->SetOptStat(0);
gStyle->SetOptTitle(0);
glog.set_warnings(false);
int cs = 1300;
TCanvas *c = new TCanvas("c","c",cs,cs*.7);
c->SetLeftMargin(0.15);
c->SetBottomMargin(0.15);
c->SetRightMargin(0.05);
c->SetTopMargin(0.05);
c->SetTickx();
c->SetTicky();
c->SetBorderMode(0);
c->SetBorderSize(0);
TH2F *h = setup_histo();
h->Draw();
Double_t r[] = {0.0, 0.0, 1.0, 1.0};
Double_t g[] = {0.0, 1.0, 1.0, 0.0};
Double_t b[] = {1.0, 1.0, 0.0, 0.0};
Double_t stop[] = {0.0, 0.3, 0.7, 1.0};
int colu=TColor::CreateGradientColorTable(4, stop, r, g, b, COLORS);
int palette[COLORS];
for(int k=0;k<COLORS;k++) palette[k]=colu+k;
double betamax = .6;
int col=0;
for(
const auto &nuc : dtm.
get_nuclei([](
auto nuc) {return (nuc->get_z()%2==0 && nuc->get_n()%2==0);})) {
auto lev_scheme = nuc->get_level_scheme();
auto lvl2 = lev_scheme->get_level("2+1",true);
auto lvl4 = lev_scheme->get_level("4+1",true);
if(!lvl2||!lvl4) continue;
if(!lvl2->has_property("energy")||!lvl4->has_property("energy")) continue;
double e2 = lvl2->get_energy();
double e4 = lvl4->get_energy();
double r42 = e4/e2;
auto gamma = lev_scheme->get_decay<
tkgammadecay>(
"2+1->0+1",
false);
if(!gamma) continue;
if(BE2W<=0||std::isnan(BE2W)) continue;
auto quad_def = nuc->get("quadrupoleDeformation");
if(!quad_def) continue;
double beta = quad_def->get_value();
if(beta>betamax) continue;
col = palette[TMath::Nint((COLORS-1)*beta/betamax)];
double radius = .008*sqrt((1+beta)*pow(nuc->get_a(),2./3.));
TEllipse* el = new TEllipse(r42,BE2W/nuc->get_a(),radius,radius*1.2*(1+beta));
el->SetLineWidth(0);
el->SetFillColorAlpha(col,.7);
el->SetLineColorAlpha(col,0);
el->Draw("same");
}
draw_legend(palette, betamax);
}
TH2F *setup_histo(double ts)
{
TH2F *h = new TH2F("h","h",100,1,3.5,100,-0.1,2);
h->GetXaxis()->SetNdivisions(209);
h->GetYaxis()->SetNdivisions(209);
h->GetXaxis()->SetLabelSize(ts);
h->GetYaxis()->SetLabelSize(ts);
h->GetXaxis()->SetTitleSize(ts+.01);
h->GetYaxis()->SetTitleSize(ts+.01);
h->GetXaxis()->SetTickLength(.02);
h->GetYaxis()->SetTickLength(.02*.7);
h->GetXaxis()->SetTitleOffset(1.1);
h->GetYaxis()->SetTitle("B(E2)/A (in Weisskopf unit)");
h->GetXaxis()->SetTitle("E(4^{+})/E(2^{+})");
return h;
}
void draw_legend(int palette[], double betamax, double ts)
{
TEllipse* el = 0;
TLatex* txt = 0;
double yy = 1.7;
double dyy = .02;
int col = 0;
double val =betamax;
el = new TEllipse(1.2,yy,.05,.05*1.2*(1+val));
col = palette[TMath::Nint((COLORS-1)*val/betamax)];
el->SetFillColorAlpha(col,.7);
el->SetLineColorAlpha(kBlack,0);
el->Draw("same");
txt = new TLatex(1.3,yy, "0.6");
txt->SetTextAlign(12);
txt->SetTextFont(132);
txt->SetTextSize(ts);
txt->Draw("same");
yy-= .05*1.2*(1+val)+dyy;
val = 0.4;
yy-= .05*1.2*(1+val);
el = new TEllipse(1.2,yy,.05,.05*1.2*(1+val));
col = palette[TMath::Nint((COLORS-1)*val/betamax)];
el->SetFillColorAlpha(col,.7);
el->SetLineColorAlpha(kBlack,0);
el->Draw("same");
txt = new TLatex(1.3,yy, "0.4");
txt->SetTextAlign(12);
txt->SetTextFont(132);
txt->SetTextSize(ts);
txt->Draw("same");
yy-= .05*1.2*(1+val)+dyy;
val = 0.2;
yy-= .05*1.2*(1+val);
el = new TEllipse(1.2,yy,.05,.05*1.2*(1+val));
col = palette[TMath::Nint((COLORS-1)*val/betamax)];
el->SetFillColorAlpha(col,.7);
el->SetLineColorAlpha(kBlack,0);
el->Draw("same");
txt = new TLatex(1.3,yy, "0.2");
txt->SetTextAlign(12);
txt->SetTextFont(132);
txt->SetTextSize(ts);
txt->Draw("same");
yy-= .05*1.2*(1+val)+dyy;
val = 0.0;
yy-= .05*1.2*(1+val);
el = new TEllipse(1.2,yy,.05,.05*1.2*(1+val));
col = palette[TMath::Nint((COLORS-1)*val/betamax)];
el->SetFillColorAlpha(col,.7);
el->SetLineColorAlpha(kBlack,0);
el->Draw("same");
txt = new TLatex(1.3,yy, "0.0");
txt->SetTextAlign(12);
txt->SetTextFont(132);
txt->SetTextSize(ts);
txt->Draw("same");
}
Stores information on a gamma-ray decay.
double get_trans_prob(bool _elec, int _L, bool _WU)
returns the gamma gamma transition probability value (ex: BE2 in Weiksopf units, get_trans_prob(1,...
Manages the database loading and provides access to the physics properties.
vector< shared_ptr< tknucleus > > get_nuclei(std::function< bool(shared_ptr< tknucleus >)>const &_selection)
return a vector containing all the known nuclei filtered by the lambda function