CMS 3D CMS Logo

bigModule.py
Go to the documentation of this file.
1 ##########################################################################
2 # Creates histograms of the modules of one structure. and returns them as
3 # a list of PlotData objects.
4 
5 
6 import logging
7 
8 import ROOT
9 ROOT.PyConfig.IgnoreCommandLineOptions = True
10 ROOT.gROOT.SetBatch()
11 
12 import Alignment.MillePedeAlignmentAlgorithm.mpsvalidate.subModule as mpsv_subModule
13 import Alignment.MillePedeAlignmentAlgorithm.mpsvalidate.classes as mpsv_classes
14 import Alignment.MillePedeAlignmentAlgorithm.mpsvalidate.style as mpsv_style
15 
16 
17 def plot(MillePedeUser, alignables, config):
18  logger = logging.getLogger("mpsvalidate")
19 
20  alignables.create_list(MillePedeUser)
21 
22  # number of bins to start
23  numberOfBins = 10000
24 
25  ######################################################################
26  # initialize data hierarchy
27  # plots[mode][struct]
28  #
29 
30  plots = []
31  # loop over mode
32  for modeNumber, mode in enumerate(["xyz", "rot", "dist"]):
33  plots.append([])
34  # loop over structures
35  for structNumber, struct in enumerate(alignables.structures):
36  plots[modeNumber].append(mpsv_classes.PlotData(mode))
37 
38  # initialize histograms
39  for modeNumber, mode in enumerate(["xyz", "rot", "dist"]):
40  for structNumber, struct in enumerate(alignables.structures):
41  # use a copy for shorter name
42  plot = plots[modeNumber][structNumber]
43 
44  for i in range(3):
45  if (mode == "xyz"):
46  plot.histo.append(ROOT.TH1F("{0} {1} {2}".format(struct.get_name(), plot.xyz[
47  i], mode), "", numberOfBins, -1000, 1000))
48  else:
49  plot.histo.append(ROOT.TH1F("{0} {1} {2}".format(struct.get_name(), plot.xyz[
50  i], mode), "", numberOfBins, -0.1, 0.1))
51 
52  if (plot.unit!=""):
53  plot.histo[i].SetXTitle("#Delta"+plot.xyz[i]+" ["+plot.unit+"]")
54  else:
55  plot.histo[i].SetXTitle("#Delta"+plot.xyz[i])
56  plot.histo[i].SetYTitle("number of alignables")
57  plot.histo[i].GetXaxis().SetTitleOffset(0.8)
58  plot.histoAxis.append(plot.histo[i].GetXaxis())
59 
60  # add labels
61  plot.title = ROOT.TPaveLabel(
62  0.1, 0.8, 0.9, 0.9, "Module: {0} {1}".format(struct.get_name(), mode))
63  plot.text = ROOT.TPaveText(0.05, 0.1, 0.95, 0.75)
64  plot.text.SetTextAlign(12)
65  plot.text.SetTextSizePixels(20)
66 
67  # save copy
68  plots[modeNumber][structNumber] = plot
69 
70  ######################################################################
71  # fill histogram
72  #
73 
74  for line in MillePedeUser:
75  # is module ?
76  if (line.ObjId == 1):
77  for modeNumber, mode in enumerate(["xyz", "rot", "dist"]):
78  for structNumber, struct in enumerate(alignables.structures):
79  # use a copy for shorter name
80  plot = plots[modeNumber][structNumber]
81 
82  # module in struct ?
83  if (struct.contains_detid(line.Id)):
84  for i in range(3):
85  if (abs(line.Par[plot.data[i]]) != 999999):
86  # transform xyz data from cm to #mu m
87  if (mode == "xyz"):
88  plot.histo[i].Fill(
89  10000 * line.Par[plot.data[i]])
90  else:
91  plot.histo[i].Fill(line.Par[plot.data[i]])
92 
93  # save copy
94  plots[modeNumber][structNumber] = plot
95 
96  ######################################################################
97  # find the best range
98  #
99 
100  for modeNumber, mode in enumerate(["xyz", "rot", "dist"]):
101  for structNumber, struct in enumerate(alignables.structures):
102  # use a copy for shorter name
103  plot = plots[modeNumber][structNumber]
104 
105  for i in range(3):
106  # get first and last bin with content and chose the one which
107  # has a greater distance to the center
108  if (abs(numberOfBins / 2 - plot.histo[i].FindFirstBinAbove()) > abs(plot.histo[i].FindLastBinAbove() - numberOfBins / 2)):
109  plot.maxBinShift[i] = abs(
110  numberOfBins / 2 - plot.histo[i].FindFirstBinAbove())
111  # set the maxShift value
112  plot.maxShift[i] = plot.histo[i].GetBinCenter(
113  plot.histo[i].FindFirstBinAbove())
114  else:
115  plot.maxBinShift[i] = abs(
116  plot.histo[i].FindLastBinAbove() - numberOfBins / 2)
117  # set the maxShift value
118  plot.maxShift[i] = plot.histo[i].GetBinCenter(
119  plot.histo[i].FindLastBinAbove())
120  # skip empty histogram
121  if (abs(plot.maxBinShift[i]) == numberOfBins / 2 + 1):
122  plot.maxBinShift[i] = 0
123 
124  # three types of ranges
125 
126  # 1. multiple of standard dev
127  if (config.rangemode == "stddev"):
128  for i in range(3):
129  if (plot.histo[i].GetEntries() != 0 and plot.histo[i].GetStdDev() != 0):
130  # if the plotrange is much bigger than the standard
131  # deviation use config.widthstdev * StdDev als Range
132  if (max(plot.maxShift) / plot.histo[i].GetStdDev() > config.defpeak):
133  # corresponding bin config.widthstdev*StdDev
134  binShift = int(plot.histo[i].FindBin(
135  config.widthstddev * plot.histo[i].GetStdDev()) - numberOfBins / 2)
136  else:
137  binShift = max(plot.maxBinShift)
138 
139  # save used binShift
140  plot.binShift[i] = binShift
141 
142  # 2. show all
143  if (config.rangemode == "all"):
144  for i in range(3):
145  plot.binShift[i] = plot.maxBinShift[i]
146 
147  # 3. use given ranges
148  if (config.rangemode == "given"):
149  for i in range(3):
150  if (mode == "xyz"):
151  valuelist = config.rangexyzM
152  if (mode == "rot"):
153  valuelist = config.rangerotM
154  if (mode == "dist"):
155  valuelist = config.rangedistM
156 
157  for value in valuelist:
158  # maximum smaller than given value
159  if (abs(plot.maxShift[i]) < value):
160  binShift = value
161  break
162  # if not possible, force highest
163  if (abs(plot.maxShift[i]) > valuelist[-1]):
164  binShift = valuelist[-1]
165  # calculate binShift
166  plot.binShift[i] = int(
167  binShift / plot.histo[i].GetBinWidth(1))
168 
169  # all plot the same range
170  if (config.samerange == 1):
171  for i in range(3):
172  plot.binShift[i] = max(plot.binShift)
173 
174  # save used range
175  for i in range(3):
176  plot.usedRange[i] = plot.binShift[i]
177 
178  # count entries which are not shown anymore
179  for i in range(3):
180  # bin 1 to begin of histogram
181  for j in range(1, numberOfBins / 2 - plot.binShift[i]):
182  plot.hiddenEntries[i] += plot.histo[i].GetBinContent(j)
183  # from the end of shown bins to the end of histogram
184  for j in range(numberOfBins / 2 + plot.binShift[i], plot.histo[i].GetNbinsX()):
185  plot.hiddenEntries[i] += plot.histo[i].GetBinContent(j)
186 
187  # apply new range
188  for i in range(3):
189  if (plot.histo[i].GetEntries() != 0):
190  # merge bins, ca. 100 should be visible in the resulting
191  # plot
192  mergeNumberBins = plot.binShift[i]
193  # skip empty histogram
194  if (mergeNumberBins != 0):
195  # the 2*maxBinShift bins should shrink to 100 bins
196  mergeNumberBins = int(
197  2. * mergeNumberBins / config.numberofbins)
198  # the total number of bins should be dividable by the
199  # bins shrinked together
200  if (mergeNumberBins == 0):
201  mergeNumberBins = 1
202  while (numberOfBins % mergeNumberBins != 0 and mergeNumberBins != 1):
203  mergeNumberBins -= 1
204 
205  # Rebin and save new created histogram and axis
206  plot.histo[i] = plot.histo[i].Rebin(mergeNumberBins)
207  plot.histoAxis[i] = plot.histo[i].GetXaxis()
208 
209  # set view range. it is important to note that the number of bins have changed with the rebinning
210  # the total number and the number of shift must be
211  # corrected with / mergeNumberBins
212  plot.histoAxis[i].SetRange(int(numberOfBins / (2 * mergeNumberBins) - plot.binShift[
213  i] / mergeNumberBins), int(numberOfBins / (2 * mergeNumberBins) + plot.binShift[i] / mergeNumberBins))
214 
215  # error if shift is bigger than limit
216  limit = config.limit[mode]
217  for i in range(3):
218  # skip empty
219  if (plot.histo[i].GetEntries() > 0):
220  if (plot.unit!=""):
221  plot.text.AddText("max. shift {0}: {1:.2} {2}".format(plot.xyz[i], plot.maxShift[i], plot.unit))
222  if (abs(plot.maxShift[i]) > limit):
223  plot.text.AddText("! {0} shift bigger than {1} {2}".format(plot.xyz[i], limit, plot.unit))
224  else:
225  plot.text.AddText("max. shift {0}: {1:.2}".format(plot.xyz[i], plot.maxShift[i]))
226  if (abs(plot.maxShift[i]) > limit):
227  plot.text.AddText("! {0} shift bigger than {1}".format(plot.xyz[i], limit))
228 
229  if (plot.hiddenEntries[i] != 0):
230  plot.text.AddText("! {0}: {1} outlier !".format(
231  plot.xyz[i], int(plot.hiddenEntries[i])))
232 
233  # save copy
234  plots[modeNumber][structNumber] = plot
235 
236  ######################################################################
237  # make the plots
238  #
239 
240  # show the skewness in the legend
241  ROOT.gStyle.SetOptStat("emrs")
242 
243  for modeNumber, mode in enumerate(["xyz", "rot", "dist"]):
244  for structNumber, struct in enumerate(alignables.structures):
245  # use a copy for shorter name
246  plot = plots[modeNumber][structNumber]
247 
248  canvas = ROOT.TCanvas("canvasModules{0}_{1}".format(
249  struct.get_name(), mode), "Parameter", 300, 0, 800, 600)
250  canvas.Divide(2, 2)
251 
252  canvas.cd(1)
253  plot.title.Draw()
254  plot.text.Draw()
255 
256  # draw identification
257  ident = mpsv_style.identification(config)
258  ident.Draw()
259 
260  # is there any plot?
261  plotNumber = 0
262 
263  # loop over coordinates
264  for i in range(3):
265  if(plot.histo[i].GetEntries() > 0):
266  plotNumber += 1
267  canvas.cd(i + 2)
268  mpsv_style.setstatsize(canvas, plot.histo[i], config)
269  plot.histo[i].DrawCopy()
270 
271  if (plotNumber == 0):
272  break
273 
274  canvas.Update()
275 
276  # save as pdf
277  canvas.Print(
278  "{0}/plots/pdf/modules_{1}_{2}.pdf".format(config.outputPath, mode, struct.get_name()))
279 
280  # export as png
281  image = ROOT.TImage.Create()
282  image.FromPad(canvas)
283  image.WriteImage(
284  "{0}/plots/png/modules_{1}_{2}.png".format(config.outputPath, mode, struct.get_name()))
285 
286  # add to output list
287  output = mpsv_classes.OutputData(plottype="mod", name=struct.get_name(),
288  parameter=mode, filename="modules_{0}_{1}".format(mode, struct.get_name()))
289  config.outputList.append(output)
290 
291  ######################################################################
292  # make plots with substructure
293  #
294 
295  if config.showsubmodule:
296  alignables.create_children_list()
297  for modeNumber, mode in enumerate(["xyz", "rot", "dist"]):
298  for structNumber, struct in enumerate(alignables.structures):
299  # use a copy for shorter name
300  plot = plots[modeNumber][structNumber]
301 
302  mpsv_subModule.plot(MillePedeUser, alignables,
303  mode, struct, plot, config)
void Fill(HcalDetId &id, double val, std::vector< TH2F > &depth)
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
def plot(MillePedeUser, alignables, config)
Definition: bigModule.py:17