File: plotff.py

package info (click to toggle)
libformfactor 0.3.1-1
  • links: PTS, VCS
  • area: main
  • in suites: forky, sid, trixie
  • size: 1,288 kB
  • sloc: cpp: 17,289; python: 382; makefile: 15
file content (138 lines) | stat: -rw-r--r-- 3,460 bytes parent folder | download | duplicates (3)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
# script to plot the form factor
# after the ROOT window with plots appear,
# just save the picture in that format which you need

import ROOT
from libBornAgainCore import *

# define a form factor, I recommend a 10--20 nm diameter
#ff = FormFactorFullSphere(10.0*nanometer)
#ff =  FormFactorSphere(10.0*nanometer, 13.0*nanometer)
ff =  FormFactorPyramid(13*nanometer, 10.0*nanometer, 60*degree)


# volume
# I suggest, that the form factor evaluated for q=0 gives the volume
# please correct if not
zero = cvector_t(0,0,0)
V = abs(ff.evaluate_for_q(zero))


# number of bins for qx, qy, qz
# I recommend to put at least 400 for a better picture quality
nqy = 400
nqz = 400
nqx = 400

# minimum and maximum values for qx, qy qz
qymin = -2.0
qymax = 2.0
qzmin = -2.0
qzmax = 2.0
qxmin = -2.0
qxmax = 2.0

# first and last bin for the projection slice
# if number of bins=400, then fbin should be set to 200 and lbin to 201
#fbin=200
#lbin=201

# step for qx, qy, qz
stepqy = (qymax - qymin)/(nqy-1)
stepqz = (qzmax - qzmin)/(nqz-1)
stepqx = (qxmax - qxmin)/(nqx-1)

# define style of the diagram
ROOT.gROOT.SetStyle("Plain")
ROOT.gStyle.SetOptStat(0);
ROOT.gStyle.SetOptTitle(0);
ROOT.gStyle.SetLabelSize(0.06, "xyz");
ROOT.gStyle.SetTitleSize(0.06, "xyz");
#ROOT.gStyle.SetPadRightMargin(0.1)
ROOT.gStyle.SetPadLeftMargin(0.18)
ROOT.gStyle.SetPadBottomMargin(0.18)

t = ROOT.TText()
# create ROOT histograms
hist = ROOT.TH2D("hist","Sphere:H=R",nqy,qymin,qymax, nqz, qzmin, qzmax)
hist2 = ROOT.TH2D("hist2","Sphere:H=R",nqx,qxmin,qxmax, nqy, qymin, qymax)

# and fill them with the values
for i in range(nqy):
    qy = qymin + i*stepqy
    for j in range(nqz):
        qz = qzmin + j*stepqz
        k = cvector_t(0,qy,qz)
        hist.Fill(qy,qz,(abs(ff.evaluate_for_q(k)))**2+1)

for i in range(nqy):
    qy = qymin + i*stepqy
    for j in range(nqx):
        qx = qxmin + j*stepqx
        k = cvector_t(qx,qy,0)
        hist2.Fill(qx,qy,(abs(ff.evaluate_for_q(k)))**2+1)


# create a ROOT canvas and put all plots on it
c = ROOT.TCanvas("FormfactorSphere","Formfactor Sphere", 1000,500)

c.Divide(2,1)
hist.SetMinimum(1)
hist.SetMaximum(V**2)
hist2.SetMinimum(1)
hist2.SetMaximum(V**2)
hist.SetContour(50)
hist2.SetContour(50)
#hist.GetZaxis().SetTitle(" |F|^{2} ")
#hist2.GetZaxis().SetTitle(" |F|^{2} ")


c.cd(1)
ROOT.gPad.SetLogz()
hist.GetXaxis().SetTitle(" q_{y} [nm^{-1}] ")
hist.GetXaxis().CenterTitle()
hist.GetXaxis().SetTitleOffset(1.1)
hist.GetYaxis().SetTitle(" q_{z} [nm^{-1}] ")
hist.GetYaxis().CenterTitle()
hist.GetYaxis().SetTitleOffset(1.4)
hist.GetZaxis().SetTitleOffset(1.3)
hist.Draw("col")
#t.SetNDC(1)
#t.SetTextSize(8.0e-2)
#t.DrawText(0.1, 0.01, "a")

c.cd(2)
ROOT.gPad.SetLogz()
hist2.SetMinimum(1)
hist2.GetXaxis().SetTitle(" q_{x} [nm^{-1}] ")
hist2.GetXaxis().CenterTitle()
hist2.GetXaxis().SetTitleOffset(1.1)
hist2.GetYaxis().SetTitle(" q_{y} [nm^{-1}] ")
hist2.GetYaxis().CenterTitle()
hist2.GetYaxis().SetTitleOffset(1.4)
hist2.GetZaxis().SetTitleOffset(1.3)
hist2.Draw("col")
#t.DrawText(0.1, 0.004, "b")

#c.cd(3)
#ROOT.gPad.SetLogy()
#py = hist.ProjectionX("py", fbin, lbin, 'o')
#py.GetYaxis().SetTitle(" |F|^{2} ")
#py.GetYaxis().SetTitleOffset(1.3)
#py.Draw()

#t.DrawText(0.1, 0.01, "c")

#c.cd(4)
#ROOT.gPad.SetLogy()
#px = hist2.ProjectionX("px", fbin, lbin, 'o')
#px.GetYaxis().SetTitle(" |F|^{2} ")
#px.GetYaxis().SetTitleOffset(1.3)
#px.Draw()
#t.DrawText(0.1, 0.004, "d")

c.Update()
ROOT.gApplication.Run()