Coverage for /usr/local/lib/python3.11/dist-packages/grond/problems/plot.py: 14%

1import math 

2import logging 

3import random 

4 

5import numpy as num 

6from matplotlib import cm, patches, colors as mcolors 

7 

8from pyrocko.guts import Tuple, Float, Int, String, List, Bool, StringChoice 

9 

10from pyrocko.plot import mpl_margins, mpl_color, mpl_init, mpl_graph_color 

11from pyrocko.plot import beachball, hudson 

12 

13from grond.plot.config import PlotConfig 

14from grond.plot.section import SectionPlotConfig, SectionPlot 

15from grond.plot.collection import PlotItem 

16from grond import meta, core, stats 

17from matplotlib import pyplot as plt 

18 

19logger = logging.getLogger('grond.problem.plot') 

20 

21guts_prefix = 'grond' 

22 

23 

24def cluster_label(icluster, perc): 

25 if icluster == -1: 

26 return 'Unclust. (%.0f%%)' % perc 

27 else: 

28 return 'Cluster %i (%.0f%%)' % (icluster, perc) 

29 

30 

31def cluster_color(icluster): 

32 if icluster == -1: 

33 return mpl_color('aluminium3') 

34 else: 

35 return mpl_graph_color(icluster) 

36 

37 

38def fixlim(lo, hi): 

39 if lo == hi: 

40 return lo - 1.0, hi + 1.0 

41 else: 

42 return lo, hi 

43 

44 

45def eigh_sorted(mat): 

46 evals, evecs = num.linalg.eigh(mat) 

47 iorder = num.argsort(evals) 

48 return evals[iorder], evecs[:, iorder] 

49 

50 

51class JointparPlot(PlotConfig): 

52 ''' 

53 Source problem parameter's tradeoff plots. 

54 ''' 

55 

56 name = 'jointpar' 

57 size_cm = Tuple.T(2, Float.T(), default=(20., 20.)) 

58 misfit_cutoff = Float.T(optional=True) 

59 ibootstrap = Int.T(optional=True) 

60 color_parameter = String.T(default='misfit') 

61 exclude = List.T(String.T()) 

62 include = List.T(String.T()) 

63 show_ellipses = Bool.T(default=False) 

64 nsubplots = Int.T(default=6) 

65 show_ticks = Bool.T(default=False) 

66 show_reference = Bool.T(default=True) 

67 show_axes = Bool.T(default=False) 

68 

69 def make(self, environ): 

70 cm = environ.get_plot_collection_manager() 

71 history = environ.get_history(subset='harvest') 

72 optimiser = environ.get_optimiser() 

73 sref = 'Dark gray boxes mark reference solution.' \ 

74 if self.show_reference else '' 

75 

76 mpl_init(fontsize=self.font_size) 

77 cm.create_group_mpl( 

78 self, 

79 self.draw_figures(history, optimiser), 

80 title=u'Jointpar Plot', 

81 section='solution', 

82 feather_icon='crosshair', 

83 description=u''' 

84Source problem parameter's scatter plots, to evaluate the resolution of source 

85parameters and possible trade-offs between pairs of model parameters. 

86 

87A subset of model solutions (from harvest) is shown in two dimensions for all 

88possible parameter pairs as points. The point color indicates the misfit for 

89the model solution with cold colors (blue) for high misfit models and warm 

90colors (red) for low misfit models. The plot ranges are defined by the given 

91parameter bounds and shows the model space of the optimsation. %s''' % sref) 

92 

93 def draw_figures(self, history, optimiser): 

94 

95 color_parameter = self.color_parameter 

96 exclude = self.exclude 

97 include = self.include 

98 nsubplots = self.nsubplots 

99 figsize = self.size_inch 

100 ibootstrap = 0 if self.ibootstrap is None else self.ibootstrap 

101 misfit_cutoff = self.misfit_cutoff 

102 show_ellipses = self.show_ellipses 

103 msize = 1.5 

104 cmap = 'coolwarm' 

105 

106 problem = history.problem 

107 if not problem: 

108 return [] 

109 

110 models = history.models 

111 

112 exclude = list(exclude) 

113 bounds = problem.get_combined_bounds() 

114 for ipar in range(problem.ncombined): 

115 par = problem.combined[ipar] 

116 lo, hi = bounds[ipar] 

117 if lo == hi: 

118 exclude.append(par.name) 

119 

120 xref = problem.get_reference_model() 

121 

122 isort = history.get_sorted_misfits_idx(chain=ibootstrap)[::-1] 

123 models = history.get_sorted_models(chain=ibootstrap)[::-1] 

124 nmodels = history.nmodels 

125 

126 gms = history.get_sorted_misfits(chain=ibootstrap)[::-1] 

127 if misfit_cutoff is not None: 

128 ibest = gms < misfit_cutoff 

129 gms = gms[ibest] 

130 models = models[ibest] 

131 

132 kwargs = {} 

133 

134 if color_parameter == 'dist': 

135 mx = num.mean(models, axis=0) 

136 cov = num.cov(models.T) 

137 mdists = core.mahalanobis_distance(models, mx, cov) 

138 icolor = meta.ordersort(mdists) 

139 

140 elif color_parameter == 'misfit': 

141 iorder = num.arange(nmodels) 

142 icolor = iorder 

143 

144 elif color_parameter in problem.parameter_names: 

145 ind = problem.name_to_index(color_parameter) 

146 icolor = problem.extract(models, ind) 

147 

148 elif color_parameter in history.attribute_names: 

149 icolor = history.get_attribute(color_parameter)[isort] 

150 icolor_need = num.unique(icolor) 

151 

152 colors = [] 

153 for i in range(icolor_need[-1]+1): 

154 colors.append(mpl_graph_color(i)) 

155 

156 cmap = mcolors.ListedColormap(colors) 

157 cmap.set_under(mpl_color('aluminium3')) 

158 kwargs.update(dict(vmin=0, vmax=icolor_need[-1])) 

159 else: 

160 raise meta.GrondError( 

161 'Invalid color_parameter: %s' % color_parameter) 

162 

163 smap = {} 

164 iselected = 0 

165 for ipar in range(problem.ncombined): 

166 par = problem.combined[ipar] 

167 if exclude and par.name in exclude or \ 

168 include and par.name not in include: 

169 continue 

170 

171 smap[iselected] = ipar 

172 iselected += 1 

173 

174 nselected = iselected 

175 

176 if nselected < 2: 

177 logger.warning( 

178 'Cannot draw joinpar figures with less than two parameters ' 

179 'selected.') 

180 return [] 

181 

182 nfig = (nselected - 2) // nsubplots + 1 

183 

184 figs = [] 

185 for ifig in range(nfig): 

186 figs_row = [] 

187 for jfig in range(nfig): 

188 if ifig >= jfig: 

189 item = PlotItem(name='fig_%i_%i' % (ifig, jfig)) 

190 item.attributes['parameters'] = [] 

191 figs_row.append((item, plt.figure(figsize=figsize))) 

192 else: 

193 figs_row.append(None) 

194 

195 figs.append(figs_row) 

196 

197 for iselected in range(nselected): 

198 ipar = smap[iselected] 

199 ypar = problem.combined[ipar] 

200 for jselected in range(iselected): 

201 jpar = smap[jselected] 

202 xpar = problem.combined[jpar] 

203 

204 ixg = (iselected - 1) 

205 iyg = jselected 

206 

207 ix = ixg % nsubplots 

208 iy = iyg % nsubplots 

209 

210 ifig = ixg // nsubplots 

211 jfig = iyg // nsubplots 

212 

213 aind = (nsubplots, nsubplots, (ix * nsubplots) + iy + 1) 

214 

215 item, fig = figs[ifig][jfig] 

216 

217 tlist = item.attributes['parameters'] 

218 if xpar.name not in tlist: 

219 tlist.append(xpar.name) 

220 if ypar.name not in tlist: 

221 tlist.append(ypar.name) 

222 

223 axes = fig.add_subplot(*aind) 

224 

225 axes.axvline(0., color=mpl_color('aluminium3'), lw=0.5) 

226 axes.axhline(0., color=mpl_color('aluminium3'), lw=0.5) 

227 for spine in axes.spines.values(): 

228 spine.set_edgecolor(mpl_color('aluminium5')) 

229 spine.set_linewidth(0.5) 

230 

231 xmin, xmax = fixlim(*xpar.scaled(bounds[jpar])) 

232 ymin, ymax = fixlim(*ypar.scaled(bounds[ipar])) 

233 

234 if not self.show_axes: 

235 if ix == 0 or jselected + 1 == iselected: 

236 for (xpos, xoff, x) in [ 

237 (0.0, 10., xmin), 

238 (1.0, -10., xmax)]: 

239 

240 axes.annotate( 

241 '%.3g%s' % (x, xpar.get_unit_suffix()), 

242 xy=(xpos, 1.05), 

243 xycoords='axes fraction', 

244 xytext=(xoff, 5.), 

245 textcoords='offset points', 

246 verticalalignment='bottom', 

247 horizontalalignment='left', 

248 rotation=45.) 

249 

250 if iy == nsubplots - 1 or jselected + 1 == iselected: 

251 for (ypos, yoff, y) in [ 

252 (0., 10., ymin), 

253 (1.0, -10., ymax)]: 

254 

255 axes.annotate( 

256 '%.3g%s' % (y, ypar.get_unit_suffix()), 

257 xy=(1.0, ypos), 

258 xycoords='axes fraction', 

259 xytext=(5., yoff), 

260 textcoords='offset points', 

261 verticalalignment='bottom', 

262 horizontalalignment='left', 

263 rotation=45.) 

264 

265 axes.set_xlim(xmin, xmax) 

266 axes.set_ylim(ymin, ymax) 

267 

268 if not self.show_axes: 

269 if not self.show_ticks: 

270 axes.get_xaxis().set_ticks([]) 

271 axes.get_yaxis().set_ticks([]) 

272 else: 

273 axes.tick_params(length=4, which='both') 

274 axes.get_yaxis().set_ticklabels([]) 

275 axes.get_xaxis().set_ticklabels([]) 

276 

277 if iselected == nselected - 1 or ix == nsubplots - 1: 

278 axes.annotate( 

279 xpar.get_label(with_unit=False), 

280 xy=(0.5, -0.05), 

281 xycoords='axes fraction', 

282 verticalalignment='top', 

283 horizontalalignment='right', 

284 rotation=45.) 

285 

286 if iy == 0: 

287 axes.annotate( 

288 ypar.get_label(with_unit=False), 

289 xy=(-0.05, 0.5), 

290 xycoords='axes fraction', 

291 verticalalignment='top', 

292 horizontalalignment='right', 

293 rotation=45.) 

294 else: 

295 axes.set_xlabel(xpar.get_label()) 

296 axes.set_ylabel(ypar.get_label()) 

297 

298 fx = problem.extract(models, jpar) 

299 fy = problem.extract(models, ipar) 

300 

301 axes.scatter( 

302 xpar.scaled(fx), 

303 ypar.scaled(fy), 

304 c=icolor, 

305 s=msize, alpha=0.5, cmap=cmap, edgecolors='none', **kwargs) 

306 

307 if show_ellipses: 

308 cov = num.cov((xpar.scaled(fx), ypar.scaled(fy))) 

309 evals, evecs = eigh_sorted(cov) 

310 evals = num.sqrt(evals) 

311 ell = patches.Ellipse( 

312 xy=( 

313 num.mean(xpar.scaled(fx)), 

314 num.mean(ypar.scaled(fy))), 

315 width=evals[0] * 2, 

316 height=evals[1] * 2, 

317 angle=num.rad2deg( 

318 num.arctan2(evecs[1][0], evecs[0][0]))) 

319 

320 ell.set_facecolor('none') 

321 ell.set_edgecolor(mpl_color('aluminium5')) 

322 axes.add_artist(ell) 

323 

324 if self.show_reference: 

325 fx = problem.extract(xref, jpar) 

326 fy = problem.extract(xref, ipar) 

327 

328 ref_color = mpl_color('aluminium6') 

329 ref_color_light = 'none' 

330 axes.plot( 

331 xpar.scaled(fx), ypar.scaled(fy), 's', 

332 mew=1.5, ms=5, mfc=ref_color_light, mec=ref_color) 

333 

334 figs_flat = [] 

335 for figs_row in figs: 

336 figs_flat.extend( 

337 item_fig for item_fig in figs_row if item_fig is not None) 

338 

339 return figs_flat 

340 

341 

342class HistogramPlot(PlotConfig): 

343 ''' 

344 Histograms or Gaussian kernel densities (default) of all parameters 

345 (marginal distributions of model parameters). 

346 

347 The histograms (by default shown as Gaussian kernel densities) show (red 

348 curved solid line) the distributions of the parameters (marginals) along 

349 with some characteristics: The red solid vertical line gives the median of 

350 the distribution and the dashed red vertical line the mean value. Dark gray 

351 vertical lines show reference values (given in the event.txt file). The 

352 overlapping red-shaded areas show the 68% confidence intervals (innermost 

353 area), the 90% confidence intervals (middle area) and the minimum and 

354 maximum values (widest area). The plot ranges are defined by the given 

355 parameter bounds and show the model space. Well resolved model parameters 

356 show peaked distributions. 

357 ''' 

358 

359 name = 'histogram' 

360 size_cm = Tuple.T(2, Float.T(), default=(12.5, 7.5)) 

361 exclude = List.T(String.T()) 

362 include = List.T(String.T()) 

363 method = StringChoice.T( 

364 choices=['gaussian_kde', 'histogram'], 

365 default='gaussian_kde') 

366 show_reference = Bool.T(default=True) 

367 

368 def make(self, environ): 

369 cm = environ.get_plot_collection_manager() 

370 history = environ.get_history(subset='harvest') 

371 

372 mpl_init(fontsize=self.font_size) 

373 cm.create_group_mpl( 

374 self, 

375 self.draw_figures(history), 

376 title=u'Histogram', 

377 section='solution', 

378 feather_icon='bar-chart-2', 

379 description=u''' 

380Distribution of the problem's parameters. 

381 

382The histograms are shown either as Gaussian kernel densities (red curved solid 

383line) or as bar plots the distributions of the parameters (marginals) along 

384with some characteristics: 

385 

386The red solid vertical line gives the median of the distribution and the dashed 

387red vertical line the mean value. Dark gray vertical lines show reference 

388parameter values if given in the event.txt file. The overlapping red-shaded 

389areas show the 68% confidence intervals (innermost area), the 90% confidence 

390intervals (middle area) and the minimum and maximum values (widest area). The 

391plot ranges are defined by the given parameter bounds and show the model 

392space. 

393''') 

394 

395 def draw_figures(self, history): 

396 

397 import scipy.stats 

398 from grond.core import make_stats 

399 

400 exclude = self.exclude 

401 include = self.include 

402 figsize = self.size_inch 

403 fontsize = self.font_size 

404 method = self.method 

405 ref_color = mpl_color('aluminium6') 

406 stats_color = mpl_color('scarletred2') 

407 bar_color = mpl_color('scarletred1') 

408 stats_color3 = mpl_color('scarletred3') 

409 

410 problem = history.problem 

411 

412 models = history.models 

413 

414 bounds = problem.get_combined_bounds() 

415 exclude = list(exclude) 

416 for ipar in range(problem.ncombined): 

417 par = problem.combined[ipar] 

418 vmin, vmax = bounds[ipar] 

419 if vmin == vmax: 

420 exclude.append(par.name) 

421 

422 xref = problem.get_reference_model() 

423 

424 smap = {} 

425 iselected = 0 

426 for ipar in range(problem.ncombined): 

427 par = problem.combined[ipar] 

428 if exclude and par.name in exclude or \ 

429 include and par.name not in include: 

430 continue 

431 

432 smap[iselected] = ipar 

433 iselected += 1 

434 

435 nselected = iselected 

436 del iselected 

437 

438 pnames = [ 

439 problem.combined[smap[iselected]].name 

440 for iselected in range(nselected)] 

441 

442 rstats = make_stats(problem, models, 

443 history.get_primary_chain_misfits(), 

444 pnames=pnames) 

445 

446 for iselected in range(nselected): 

447 ipar = smap[iselected] 

448 par = problem.combined[ipar] 

449 vs = problem.extract(models, ipar) 

450 vmin, vmax = bounds[ipar] 

451 

452 fig = plt.figure(figsize=figsize) 

453 labelpos = mpl_margins( 

454 fig, nw=1, nh=1, w=7., bottom=5., top=1, units=fontsize) 

455 

456 axes = fig.add_subplot(1, 1, 1) 

457 labelpos(axes, 2.5, 2.0) 

458 axes.set_xlabel(par.get_label()) 

459 axes.set_ylabel('PDF') 

460 axes.set_xlim(*fixlim(*par.scaled((vmin, vmax)))) 

461 

462 if method == 'gaussian_kde': 

463 try: 

464 kde = scipy.stats.gaussian_kde(vs) 

465 except Exception: 

466 logger.warning( 

467 'Cannot create plot histogram with gaussian_kde: ' 

468 'possibly all samples have the same value.') 

469 continue 

470 

471 vps = num.linspace(vmin, vmax, 600) 

472 pps = kde(vps) 

473 

474 axes.plot( 

475 par.scaled(vps), par.inv_scaled(pps), color=stats_color) 

476 

477 elif method == 'histogram': 

478 pps, edges = num.histogram( 

479 vs, 

480 bins=num.linspace(vmin, vmax, num=40), 

481 density=True) 

482 vps = 0.5 * (edges[:-1] + edges[1:]) 

483 

484 axes.bar(par.scaled(vps), par.inv_scaled(pps), 

485 par.scaled(2.*(vps - edges[:-1])), 

486 color=bar_color) 

487 

488 pstats = rstats.parameter_stats_list[iselected] 

489 

490 axes.axvspan( 

491 par.scaled(pstats.minimum), 

492 par.scaled(pstats.maximum), 

493 color=stats_color, alpha=0.1) 

494 axes.axvspan( 

495 par.scaled(pstats.percentile16), 

496 par.scaled(pstats.percentile84), 

497 color=stats_color, alpha=0.1) 

498 axes.axvspan( 

499 par.scaled(pstats.percentile5), 

500 par.scaled(pstats.percentile95), 

501 color=stats_color, alpha=0.1) 

502 

503 axes.axvline( 

504 par.scaled(pstats.median), 

505 color=stats_color3, alpha=0.5) 

506 axes.axvline( 

507 par.scaled(pstats.mean), 

508 color=stats_color3, ls=':', alpha=0.5) 

509 

510 if self.show_reference: 

511 axes.axvline( 

512 par.scaled(problem.extract(xref, ipar)), 

513 color=ref_color) 

514 

515 item = PlotItem(name=par.name) 

516 item.attributes['parameters'] = [par.name] 

517 yield item, fig 

518 

519 

520class MTDecompositionPlot(PlotConfig): 

521 ''' 

522 Moment tensor decomposition plot. 

523 ''' 

524 

525 name = 'mt_decomposition' 

526 size_cm = Tuple.T(2, Float.T(), default=(15., 5.)) 

527 cluster_attribute = meta.StringID.T( 

528 optional=True, 

529 help='name of attribute to use as cluster IDs') 

530 show_reference = Bool.T(default=True) 

531 

532 def make(self, environ): 

533 cm = environ.get_plot_collection_manager() 

534 history = environ.get_history(subset='harvest') 

535 mpl_init(fontsize=self.font_size) 

536 cm.create_group_mpl( 

537 self, 

538 self.draw_figures(history), 

539 title=u'MT Decomposition', 

540 section='solution', 

541 feather_icon='sun', 

542 description=u''' 

543Moment tensor decomposition of the best-fitting solution into isotropic, 

544deviatoric and best double couple components. 

545 

546Shown are the ensemble best, the ensemble mean%s and, if available, a reference 

547mechanism. The symbol size indicates the relative strength of the components. 

548The inversion result is consistent and stable if ensemble mean and ensemble 

549best have similar symbol size and patterns. 

550''' % (', cluster results' if self.cluster_attribute else '')) 

551 

552 def draw_figures(self, history): 

553 

554 fontsize = self.font_size 

555 

556 fig = plt.figure(figsize=self.size_inch) 

557 axes = fig.add_subplot(1, 1, 1, aspect=1.0) 

558 fig.subplots_adjust(left=0., right=1., bottom=0., top=1.) 

559 

560 problem = history.problem 

561 models = history.models 

562 

563 if models.size == 0: 

564 logger.warning('Empty models vector.') 

565 return [] 

566 

567 # gms = problem.combine_misfits(history.misfits) 

568 # isort = num.argsort(gms) 

569 # iorder = num.empty_like(isort) 

570 # iorder[isort] = num.arange(iorder.size)[::-1] 

571 

572 ref_source = problem.base_source 

573 

574 mean_source = stats.get_mean_source( 

575 problem, history.models) 

576 

577 best_source = history.get_best_source() 

578 

579 nlines_max = int(round(self.size_cm[1] / 5. * 4. - 1.0)) 

580 

581 if self.cluster_attribute: 

582 cluster_sources = history.mean_sources_by_cluster( 

583 self.cluster_attribute) 

584 else: 

585 cluster_sources = [] 

586 

587 def get_deco(source): 

588 mt = source.pyrocko_moment_tensor() 

589 return mt.standard_decomposition() 

590 

591 lines = [] 

592 lines.append( 

593 ('Ensemble best', get_deco(best_source), mpl_color('aluminium5'))) 

594 

595 lines.append( 

596 ('Ensemble mean', get_deco(mean_source), mpl_color('aluminium5'))) 

597 

598 for (icluster, perc, source) in cluster_sources: 

599 if len(lines) < nlines_max - int(self.show_reference): 

600 lines.append( 

601 (cluster_label(icluster, perc), 

602 get_deco(source), 

603 cluster_color(icluster))) 

604 else: 

605 logger.warning( 

606 'Skipping display of cluster %i because figure height is ' 

607 'too small. Figure height should be at least %g cm.' % ( 

608 icluster, (3 + len(cluster_sources) 

609 + int(self.show_reference)) * 5/4.)) 

610 

611 if self.show_reference: 

612 lines.append( 

613 ('Reference', get_deco(ref_source), mpl_color('aluminium3'))) 

614 

615 moment_full_max = max(deco[-1][0] for (_, deco, _) in lines) 

616 

617 for xpos, label in [ 

618 (0., 'Full'), 

619 (2., 'Isotropic'), 

620 (4., 'Deviatoric'), 

621 (6., 'CLVD'), 

622 (8., 'DC')]: 

623 

624 axes.annotate( 

625 label, 

626 xy=(1 + xpos, nlines_max), 

627 xycoords='data', 

628 xytext=(0., 0.), 

629 textcoords='offset points', 

630 ha='center', 

631 va='center', 

632 color='black', 

633 fontsize=fontsize) 

634 

635 for i, (label, deco, color_t) in enumerate(lines): 

636 ypos = nlines_max - i - 1.0 

637 

638 [(moment_iso, ratio_iso, m_iso), 

639 (moment_dc, ratio_dc, m_dc), 

640 (moment_clvd, ratio_clvd, m_clvd), 

641 (moment_devi, ratio_devi, m_devi), 

642 (moment_full, ratio_full, m_full)] = deco 

643 

644 size0 = moment_full / moment_full_max 

645 

646 axes.annotate( 

647 label, 

648 xy=(-2., ypos), 

649 xycoords='data', 

650 xytext=(0., 0.), 

651 textcoords='offset points', 

652 ha='left', 

653 va='center', 

654 color='black', 

655 fontsize=fontsize) 

656 

657 for xpos, mt_part, ratio, ops in [ 

658 (0., m_full, ratio_full, '-'), 

659 (2., m_iso, ratio_iso, '='), 

660 (4., m_devi, ratio_devi, '='), 

661 (6., m_clvd, ratio_clvd, '+'), 

662 (8., m_dc, ratio_dc, None)]: 

663 

664 if ratio > 1e-4: 

665 try: 

666 beachball.plot_beachball_mpl( 

667 mt_part, axes, 

668 beachball_type='full', 

669 position=(1. + xpos, ypos), 

670 size=0.9 * size0 * math.sqrt(ratio), 

671 size_units='data', 

672 color_t=color_t, 

673 linewidth=1.0) 

674 

675 except beachball.BeachballError as e: 

676 logger.warning(str(e)) 

677 

678 axes.annotate( 

679 'ERROR', 

680 xy=(1. + xpos, ypos), 

681 ha='center', 

682 va='center', 

683 color='red', 

684 fontsize=fontsize) 

685 

686 else: 

687 axes.annotate( 

688 'N/A', 

689 xy=(1. + xpos, ypos), 

690 ha='center', 

691 va='center', 

692 color='black', 

693 fontsize=fontsize) 

694 

695 if ops is not None: 

696 axes.annotate( 

697 ops, 

698 xy=(2. + xpos, ypos), 

699 ha='center', 

700 va='center', 

701 color='black', 

702 fontsize=fontsize) 

703 

704 axes.axison = False 

705 axes.set_xlim(-2.25, 9.75) 

706 axes.set_ylim(-0.5, nlines_max+0.5) 

707 

708 item = PlotItem(name='main') 

709 return [[item, fig]] 

710 

711 

712class MTLocationPlot(SectionPlotConfig): 

713 ''' MT location plot of the best solutions in three cross-sections. ''' 

714 name = 'location_mt' 

715 beachball_type = StringChoice.T( 

716 choices=['full', 'deviatoric', 'dc'], 

717 default='dc') 

718 normalisation_gamma = Float.T( 

719 default=3., 

720 help='Normalisation of colors and alpha as :math:`x^\\gamma`.' 

721 'A linear colormap/alpha with :math:`\\gamma=1`.') 

722 

723 def make(self, environ): 

724 environ.setup_modelling() 

725 cm = environ.get_plot_collection_manager() 

726 history = environ.get_history(subset='harvest') 

727 mpl_init(fontsize=self.font_size) 

728 self._to_be_closed = [] 

729 cm.create_group_mpl( 

730 self, 

731 self.draw_figures(history), 

732 title=u'MT Location', 

733 section='solution', 

734 feather_icon='target', 

735 description=u''' 

736Location plot of the ensemble of best solutions in three cross-sections. 

737 

738The coordinate range is defined by the search space given in the config file. 

739Symbols show best double-couple mechanisms, and colors indicate low (red) and 

740high (blue) misfit. 

741''') 

742 for obj in self._to_be_closed: 

743 obj.close() 

744 

745 def draw_figures(self, history, color_p_axis=False): 

746 from matplotlib import colors 

747 

748 color = 'black' 

749 fontsize = self.font_size 

750 markersize = fontsize * 1.5 

751 beachballsize_small = markersize * 0.5 

752 beachball_type = self.beachball_type 

753 

754 problem = history.problem 

755 sp = SectionPlot(config=self) 

756 self._to_be_closed.append(sp) 

757 

758 fig = sp.fig 

759 axes_en = sp.axes_xy 

760 axes_dn = sp.axes_zy 

761 axes_ed = sp.axes_xz 

762 

763 bounds = problem.get_combined_bounds() 

764 

765 models = history.get_sorted_primary_models()[::-1] 

766 

767 iorder = num.arange(history.nmodels) 

768 

769 for parname, set_label, set_lim in [ 

770 ['east_shift', sp.set_xlabel, sp.set_xlim], 

771 ['north_shift', sp.set_ylabel, sp.set_ylim], 

772 ['depth', sp.set_zlabel, sp.set_zlim]]: 

773 

774 ipar = problem.name_to_index(parname) 

775 par = problem.combined[ipar]