TkN 2.1
Toolkit for Nuclei
deformation_island.C
1#include "TMultiGraph.h"
2#include "TGraph.h"
3#include "tkmanager.h"
4#include "TROOT.h"
5#include "TCanvas.h"
6#include "TAxis.h"
7#include "TColor.h"
8#include "TLegend.h"
9#include "TLegendEntry.h"
10#include "TLine.h"
11
12using namespace tkn;
13
14const char* ccs[7] = {"#1ddbdc","#1aabbf","#1a83a5","#1d6190","#1e427a","#182562","#120952"};
15int cols[7];
16int marks[7] = {20,33,21,47,22,34,23};
17double masize[7] = {1.3,1.8,1,1.4,1.4,1.4,1.4};
18
19TGraph* get_graph(int zz);
20TLegend* get_legend();
21TLine* get_line(TMultiGraph* mg);
22
23void draw_deformation_island(){
24 gROOT->SetStyle("tkn-histo");
25
26 for(int ii=0; ii<7; ii++) cols[ii] = TColor::GetColor(ccs[ii]);
27 tkmanager dtm;
28
29 TMultiGraph* mgr = new TMultiGraph;
30 TMultiGraph* mgs = new TMultiGraph;
31 TMultiGraph* mgb2 = new TMultiGraph;
32 TMultiGraph* mgr42 = new TMultiGraph;
33
34 for(int zz=32; zz<=45; zz+=2)
35 {
36 TGraph* gr = get_graph(zz);
37 TGraph* gs = get_graph(zz);
38 TGraph* gb2 = get_graph(zz);
39 TGraph* gr42 = get_graph(zz);
40
41 for(int nn=48; nn<=68; nn+=2)
42 {
43 auto nuc = dtm.get_nucleus(zz,nn+zz);
44 if(!nuc) continue;
45
46 if(nuc->has_property("radius"))
47 gr->SetPoint(gr->GetN(),nuc->get_n(), nuc->get("radius")->get_value());
48
49 if(nuc->has_property("twoNeutronSeparationEnergy")&&nuc->get("twoNeutronSeparationEnergy")->get_info()==tkn::tkproperty::kKnown)
50 gs->SetPoint(gs->GetN(),nuc->get_n(), nuc->get("twoNeutronSeparationEnergy")->get_value()*0.001);
51
52 auto lev_scheme = nuc->get_level_scheme();
53 auto gamma = lev_scheme->get_decay<tkgammadecay>("2+1->0+1",false);
54
55 if(!gamma){}
56 else{
57 auto BE2W = gamma->get_trans_prob(true,2,true);
58 if(BE2W<=0||std::isnan(BE2W)) {}
59 else gb2->SetPoint(gb2->GetN(),nuc->get_n(),BE2W/(nuc->get_a()));
60 }
61
62 auto lvl2 = lev_scheme->get_level("2+1",false);
63 if(!lvl2||!lvl2->has_property("energy")) {}
64 else{
65 double e2 = lvl2->get_energy();
66 gr42->SetPoint(gr42->GetN(),nuc->get_n(),e2*0.001);
67 }
68 }
69
70 if(gr->GetN()>1) mgr->Add(gr,"PL");
71 else delete gr;
72
73 if(gs->GetN()>1) mgs->Add(gs,"PL");
74 else delete gs;
75
76 if(gb2->GetN()>1) mgb2->Add(gb2,"PL");
77 else delete gb2;
78
79 if(gr42->GetN()>1) mgr42->Add(gr42,"PL");
80 else delete gr42;
81 }
82
83 TCanvas* cc = new TCanvas("Deformation","Deformation island",1200,1000);
84 cc->Divide(2,2);
85
86 cc->cd(1);
87 mgs->Draw("PAL");
88 mgs->GetYaxis()->SetTitle("S_{2n} (MeV)");
89 mgs->GetXaxis()->SetTitle("N");
90 get_line(mgs)->Draw();
91
92 cc->cd(2);
93 mgr->Draw("PAL");
94 mgr->GetYaxis()->SetTitle("r_{ch} (fm)");
95 mgr->GetXaxis()->SetTitle("N");
96 get_line(mgr)->Draw();
97
98 cc->cd(3);
99 mgb2->Draw("PAL");
100 mgb2->GetYaxis()->SetTitle("B(E2)/A (W.u./A)");
101 mgb2->GetXaxis()->SetTitle("N");
102 get_legend()->Draw();
103 get_line(mgb2)->Draw();
104
105 cc->cd(4);
106 mgr42->Draw("PAL");
107 mgr42->GetYaxis()->SetTitle("E(2^{+}) (MeV)");
108 mgr42->GetXaxis()->SetTitle("N");
109 get_line(mgr42)->Draw();
110
111 return;
112}
113
114
115TGraph* get_graph(int zz)
116{
117 TGraph* gg = new TGraph;
118 gg->SetMarkerStyle(marks[zz/2-16]);
119 gg->SetMarkerSize(masize[zz/2-16]);
120 gg->SetMarkerColor(cols[zz/2-16]);
121 gg->SetLineColor(cols[zz/2-16]);
122 gg->SetName(Form("Z=%d",zz));
123
124 return gg;
125}
126
127TLine* get_line(TMultiGraph* mg)
128{
129 TLine *ll = new TLine(60,mg->GetYaxis()->GetXmin(),60,mg->GetYaxis()->GetXmax());
130 ll->SetLineStyle(9);
131 return ll;
132}
133
134
135TLegend* get_legend()
136{
137 TLegend *leg = new TLegend(0.203499,0.35485,0.542606,0.937585,NULL,"brNDC");
138 leg->SetBorderSize(0);
139 leg->SetTextFont(132);
140 leg->SetLineColor(1);
141 leg->SetLineStyle(1);
142 leg->SetLineWidth(1);
143 leg->SetFillColor(0);
144 leg->SetFillStyle(1001);
145
146 TLegendEntry *entry=0;
147 Int_t ci;
148
149 entry=leg->AddEntry("Z=44","Z=44","lpf");
150 entry->SetFillStyle(1000);
151 ci = TColor::GetColor("#120952");
152 entry->SetLineColor(ci);
153 entry->SetLineStyle(1);
154 entry->SetLineWidth(1);
155 entry->SetMarkerColor(ci);
156 entry->SetMarkerStyle(20);
157 entry->SetMarkerSize(1.3);
158 entry->SetTextFont(132);
159
160 entry=leg->AddEntry("Z=42","Z=42","lpf");
161 entry->SetFillStyle(1000);
162 ci = TColor::GetColor("#182562");
163 entry->SetLineColor(ci);
164 entry->SetLineStyle(1);
165 entry->SetLineWidth(1);
166 entry->SetMarkerColor(ci);
167 entry->SetMarkerStyle(20);
168 entry->SetMarkerSize(1.3);
169 entry->SetTextFont(132);
170
171 entry=leg->AddEntry("Z=40","Z=40","lpf");
172 entry->SetFillStyle(1000);
173 ci = TColor::GetColor("#1e427a");
174 entry->SetLineColor(ci);
175 entry->SetLineStyle(1);
176 entry->SetLineWidth(1);
177 entry->SetMarkerColor(ci);
178 entry->SetMarkerStyle(20);
179 entry->SetMarkerSize(1.3);
180 entry->SetTextFont(132);
181
182 entry=leg->AddEntry("Z=38","Z=38","lpf");
183 entry->SetFillStyle(1000);
184 ci = TColor::GetColor("#1d6190");
185 entry->SetLineColor(ci);
186 entry->SetLineStyle(1);
187 entry->SetLineWidth(1);
188 entry->SetMarkerColor(ci);
189 entry->SetMarkerStyle(20);
190 entry->SetMarkerSize(1.3);
191 entry->SetTextFont(132);
192
193 entry=leg->AddEntry("Z=36","Z=36","lpf");
194 entry->SetFillStyle(1000);
195 ci = TColor::GetColor("#1a83a5");
196 entry->SetLineColor(ci);
197 entry->SetLineStyle(1);
198 entry->SetLineWidth(1);
199 entry->SetMarkerColor(ci);
200 entry->SetMarkerStyle(20);
201 entry->SetMarkerSize(1.3);
202 entry->SetTextFont(132);
203
204 entry=leg->AddEntry("Z=34","Z=34","lpf");
205 entry->SetFillStyle(1000);
206 ci = TColor::GetColor("#1aabbf");
207 entry->SetLineColor(ci);
208 entry->SetLineStyle(1);
209 entry->SetLineWidth(1);
210 entry->SetMarkerColor(ci);
211 entry->SetMarkerStyle(20);
212 entry->SetMarkerSize(1.3);
213 entry->SetTextFont(132);
214
215 entry=leg->AddEntry("Z=32","Z=32","lpf");
216 entry->SetFillStyle(1000);
217 ci = TColor::GetColor("#1ddbdc");
218 entry->SetLineColor(ci);
219 entry->SetLineStyle(1);
220 entry->SetLineWidth(1);
221 entry->SetMarkerColor(ci);
222 entry->SetMarkerStyle(20);
223 entry->SetMarkerSize(1.3);
224 entry->SetTextFont(132);
225
226 return leg;
227}
Stores information on a gamma-ray decay.
Definition: tkdecay.h:131
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,...
Definition: tkdecay.cpp:262
Manages the database loading and provides access to the physics properties.
Definition: tkmanager.h:55
shared_ptr< tknucleus > get_nucleus(const tkstring &_nuc)
return a shared pointer to a nucleus from its name
Definition: tkmanager.h:126
Definition: tklog.cpp:39