This commit is contained in:
Noe Brucy
2020-04-21 10:40:29 +02:00
committed by Noe Brucy
parent 77676ec9a4
commit 9e5fc8b4ad
4 changed files with 141 additions and 170 deletions
+1 -1
View File
@@ -18,7 +18,7 @@ from pymses.filters import CellsToPoints
from pymses.analysis import ScalarOperator, FractionOperator, MaxLevelOperator from pymses.analysis import ScalarOperator, FractionOperator, MaxLevelOperator
import subprocess import subprocess
import module_extract as me # import module_extract as me
from mypool import MyPool as Pool from mypool import MyPool as Pool
from functools import partial from functools import partial
+109
View File
@@ -0,0 +1,109 @@
import matplotlib.pyplot as plt
from postprocessor import *
class InteractiveGUI:
"""
This is a matplotlib interactive session to restrain analysis to a specific area
"""
def onbuttonrelease(self, event):
"""Deal with click events"""
button = ["left", "middle", "right"]
toolbar = plt.get_current_fig_manager().toolbar
if toolbar.mode == "zoom rect" and event.inaxes == self.ax_col:
print("zooming ")
xlim = self.ax_col.get_xlim()
ylim = self.ax_col.get_ylim()
self.reset_mask()
elif self.add_mask and event.inaxes == self.ax_col:
self.plot_side()
plt.draw()
def onbuttonpress(self, event):
"""Deal with click events"""
button = ["left", "middle", "right"]
toolbar = plt.get_current_fig_manager().toolbar
if toolbar.mode != "":
print(
"You clicked on something, but toolbar is in mode {:s}.".format(
toolbar.mode
)
)
print(self.add_mask)
if self.add_mask and toolbar.mode == "" and event.inaxes == self.ax_col:
ix, iy = event.xdata, event.ydata
print("Add patch {}, {}".format(ix, iy))
xlim = self.ax_col.get_xlim()
ylim = self.ax_col.get_ylim()
radius = 0.05 * min(abs(xlim[1] - xlim[0]), abs(ylim[1] - ylim[0]))
circle = mpatches.Circle(
[ix, iy], radius, color="black", alpha=0.1, ec="none"
)
self.circles.append(circle)
self.ax_col.add_artist(circle)
self.ax_col.draw_artist(circle)
self.patch_mask = self.patch_mask | (
(self.xx - ix) ** 2 + (self.yy - iy) ** 2 < radius ** 2
)
# self.plot_side()
def onkeypress(self, event):
"""whenever a key is pressed"""
if not event.inaxes:
return
if event.key == "t":
self.add_mask = not self.add_mask
print("Add mode is {}".format(self.add_mask))
elif event.key == "r":
self.reset_mask()
def plot_side(self):
if self.add_mask:
mask = (self.patch_mask & self.mask).flatten()
else:
mask = self.mask.flatten()
self.ax_gamma.clear()
plt.sca(self.ax_gamma)
plot_dcsdrho(self.fluct_maps, mask, tag=self.tag)
self.ax_pdf.clear()
plt.sca(self.ax_pdf)
sigma_pdf(self.fluct_maps, mask, tag=self.tag, nb_bin_hist=self.args.pdf_nb_bin)
def reset_mask(self):
xlim = self.ax_col.get_xlim()
ylim = self.ax_col.get_ylim()
self.mask = (
(self.xx >= xlim[0])
& (self.xx <= xlim[1])
& (self.yy >= ylim[0])
& (self.yy <= ylim[1])
)
self.patch_mask = np.full(self.mask.shape, False)
for circle in self.circles:
circle.remove()
self.circles = []
self.plot_side()
plt.draw()
def __init__(self, num, path="./", pp_params=None):
"""
Interactive plotting
"""
pp = PostProcessor(path, num, pp_params=pp_params, tag="interactive")
pp.pdf_coldens("z")
fluct_map = pp.get_value("/maps/avg_map_coldens_z")
rr = pp.get_value("/maps/rr_z")
fig, ax = plt.subplots(2)
im = ax[0].imshow(fluct_map, origin="lower", cmap="RdBu_r")
cbar = plt.colorbar()
fig.canvas.mpl_connect("button_release_event", self.onbuttonrelease)
fig.canvas.mpl_connect("button_press_event", self.onbuttonpress)
fig.canvas.mpl_connect("key_press_event", self.onkeypress)
plt.tight_layout()
plt.show()
+14 -157
View File
@@ -2,26 +2,22 @@
import sys import sys
import os import os
from functools import partial
import tables import tables
import numpy as np import numpy as np
from scipy.stats import linregress
from numpy.polynomial.polynomial import polyfit
from scipy.ndimage.filters import gaussian_filter1d
from scipy import optimize
import matplotlib as mpl import matplotlib as mpl
if os.environ.get("DISPLAY", "") == "": if os.environ.get("DISPLAY", "") == "":
print("No display found. Using non-interactive Agg backend") print("No display found. Using non-interactive Agg backend")
mpl.use("Agg") mpl.use("Agg")
from matplotlib.patches import Polygon
import pylab as P import pylab as P
from scipy.stats import linregress
import matplotlib.patches as mpatches
from matplotlib.collections import PatchCollection
from functools import partial
from numpy.polynomial.polynomial import polyfit
from scipy.ndimage.filters import gaussian_filter1d
from scipy import optimize
from comparator import * from comparator import *
P.rcParams["image.cmap"] = "plasma" P.rcParams["image.cmap"] = "plasma"
P.rcParams["savefig.dpi"] = 400 P.rcParams["savefig.dpi"] = 400
@@ -105,6 +101,8 @@ class Plotter(Aggregator, BaseProcessor):
# Define rules # Define rules
self.def_rules() self.def_rules()
self.save = None
def _not_self_dep(self, name, dep, dep_arg, overwrite, **kwargs): def _not_self_dep(self, name, dep, dep_arg, overwrite, **kwargs):
if dep in self.comp.rules: if dep in self.comp.rules:
done = self.comp.process(dep, dep_arg, overwrite, overwrite) done = self.comp.process(dep, dep_arg, overwrite, overwrite)
@@ -460,7 +458,7 @@ class Plotter(Aggregator, BaseProcessor):
coordinates="figure", coordinates="figure",
) )
def _plot_radial(self, name, label=None, xlog=False, ylog=False): def _plot_radial(self, name, ax_los, label=None, xlog=False, ylog=False):
radial_bins = self.save.get_node("/radial/radial_bins_" + ax_los).read() radial_bins = self.save.get_node("/radial/radial_bins_" + ax_los).read()
bin_centers = 0.5 * (radial_bins[1:] + radial_bins[:-1]) bin_centers = 0.5 * (radial_bins[1:] + radial_bins[:-1])
@@ -481,8 +479,8 @@ class Plotter(Aggregator, BaseProcessor):
def _plot_hist( def _plot_hist(
self, self,
name, name,
ax_los, ax_los=None,
run, run="",
label=None, label=None,
unit=None, unit=None,
unit_coeff=1.0, unit_coeff=1.0,
@@ -500,6 +498,9 @@ class Plotter(Aggregator, BaseProcessor):
**kwargs **kwargs
): ):
if not ax_los is None:
name = name + "_" + ax_los
node = self.save.get_node("/hist/" + name) node = self.save.get_node("/hist/" + name)
if xlog is None: if xlog is None:
@@ -1045,147 +1046,3 @@ class Plotter(Aggregator, BaseProcessor):
) )
super(Plotter, self).def_rules() super(Plotter, self).def_rules()
class InteractiveGUI:
"""
This is a matplotlib interactive session to restrain analysis to a specific area
"""
def onbuttonrelease(self, event):
"""Deal with click events"""
button = ["left", "middle", "right"]
toolbar = P.get_current_fig_manager().toolbar
if toolbar.mode == "zoom rect" and event.inaxes == self.ax_col:
print("zooming ")
xlim = self.ax_col.get_xlim()
ylim = self.ax_col.get_ylim()
self.reset_mask()
elif self.add_mask and event.inaxes == self.ax_col:
self.plot_side()
P.draw()
def onbuttonpress(self, event):
"""Deal with click events"""
button = ["left", "middle", "right"]
toolbar = P.get_current_fig_manager().toolbar
if toolbar.mode != "":
print(
"You clicked on something, but toolbar is in mode {:s}.".format(
toolbar.mode
)
)
print(self.add_mask)
if self.add_mask and toolbar.mode == "" and event.inaxes == self.ax_col:
ix, iy = event.xdata, event.ydata
print("Add patch {}, {}".format(ix, iy))
xlim = self.ax_col.get_xlim()
ylim = self.ax_col.get_ylim()
radius = 0.05 * min(abs(xlim[1] - xlim[0]), abs(ylim[1] - ylim[0]))
circle = mpatches.Circle(
[ix, iy], radius, color="black", alpha=0.1, ec="none"
)
self.circles.append(circle)
self.ax_col.add_artist(circle)
self.ax_col.draw_artist(circle)
self.patch_mask = self.patch_mask | (
(self.xx - ix) ** 2 + (self.yy - iy) ** 2 < radius ** 2
)
# self.plot_side()
def onkeypress(self, event):
"""whenever a key is pressed"""
if not event.inaxes:
return
if event.key == "t":
self.add_mask = not self.add_mask
print("Add mode is {}".format(self.add_mask))
elif event.key == "r":
self.reset_mask()
def plot_side(self):
if self.add_mask:
mask = (self.patch_mask & self.mask).flatten()
else:
mask = self.mask.flatten()
self.ax_gamma.clear()
P.sca(self.ax_gamma)
plot_dcsdrho(self.fluct_maps, mask, tag=self.tag)
self.ax_pdf.clear()
P.sca(self.ax_pdf)
sigma_pdf(self.fluct_maps, mask, tag=self.tag, nb_bin_hist=self.args.pdf_nb_bin)
def reset_mask(self):
xlim = self.ax_col.get_xlim()
ylim = self.ax_col.get_ylim()
self.mask = (
(self.xx >= xlim[0])
& (self.xx <= xlim[1])
& (self.yy >= ylim[0])
& (self.yy <= ylim[1])
)
self.patch_mask = np.full(self.mask.shape, False)
for circle in self.circles:
circle.remove()
self.circles = []
self.plot_side()
P.draw()
def __init__(self, path, run, num, path_out=None, pp_params=default_params()):
"""
Interactive plotting
"""
self.add_mask = False
self.circles = []
self.tag = tag
# Get plotter object
self.plot = Plotter(path, [run], [num], path_out, pp_params)
im_extent = maps_disk["im_extent"]
fig = P.figure()
self.ax_col = P.subplot(1, 2, 1)
coldens = maps_disk["coldens_z"]
im = self.ax_col.imshow(
coldens, extent=im_extent, origin="lower", norm=mpl.colors.LogNorm()
)
if set_lim:
im.set_clim(0.01, 100)
self.ax_col.set_xlabel(r"$x$")
self.ax_col.set_ylabel(r"$y$")
self.xx, self.yy, self.fluct_maps = disk_pdf(
path,
num,
maps_disk,
tag=self.tag,
force=True,
put_title=False,
interactive=True,
)
coord_flat = zip(self.xx.flatten(), self.yy.flatten())
self.ax_gamma = P.subplot(2, 2, 2)
self.ax_pdf = P.subplot(2, 2, 4)
xlim = self.ax_col.get_xlim()
ylim = self.ax_col.get_ylim()
self.mask = (
(self.xx >= xlim[0])
& (self.xx <= xlim[1])
& (self.yy >= ylim[0])
& (self.yy <= ylim[1])
)
self.patch_mask = np.full(self.mask.shape, False)
self.plot_side()
fig.canvas.mpl_connect("button_release_event", self.onbuttonrelease)
fig.canvas.mpl_connect("button_press_event", self.onbuttonpress)
fig.canvas.mpl_connect("key_press_event", self.onkeypress)
P.tight_layout()
P.show()
+17 -12
View File
@@ -34,14 +34,16 @@ class PostProcessor(HDF5Container):
else: else:
tag_name = "" tag_name = ""
self.filename = path_out + "/postproc_" + tag_name + format(num, "05") + ".h5" self.filename = (
self.path_out + "/postproc_" + tag_name + format(num, "05") + ".h5"
self.cells_filename = (
path_out + "/cells_" + tag_name + format(num, "05") + ".h5"
) )
if not os.path.exists(path_out): self.cells_filename = (
os.makedirs(path_out) self.path_out + "/cells_" + tag_name + format(num, "05") + ".h5"
)
if not os.path.exists(self.path_out):
os.makedirs(self.path_out)
self.open() self.open()
@@ -50,7 +52,10 @@ class PostProcessor(HDF5Container):
self.run = os.path.basename(path) self.run = os.path.basename(path)
self.num = num self.num = num
self._ro = pymses.RamsesOutput( self._ro = pymses.RamsesOutput(
path, num, order=pp_params.pymses.order, verbose=pp_params.pymses.verbose path,
num,
order=self.pp_params.pymses.order,
verbose=self.pp_params.pymses.verbose,
) )
self._amr = self._ro.amr_source(self.pp_params.pymses.variables) self._amr = self._ro.amr_source(self.pp_params.pymses.variables)
self.info = self._ro.info.copy() self.info = self._ro.info.copy()
@@ -61,7 +66,7 @@ class PostProcessor(HDF5Container):
) )
# Density ray tracer # Density ray tracer
if pp_params.pymses.fft: if self.pp_params.pymses.fft:
self._rt = splatting.SplatterProcessor( self._rt = splatting.SplatterProcessor(
self._amr, self._ro.info, self._rho_op self._amr, self._ro.info, self._rho_op
) )
@@ -69,9 +74,9 @@ class PostProcessor(HDF5Container):
self._rt = raytracing.RayTracer(self._amr, self._ro.info, self._rho_op) self._rt = raytracing.RayTracer(self._amr, self._ro.info, self._rho_op)
# Set the extend of the image # Set the extend of the image
self._radius = 0.5 / pp_params.pymses.zoom self._radius = 0.5 / self.pp_params.pymses.zoom
self._lbox = self.info["boxlen"] self._lbox = self.info["boxlen"]
center = pp_params.pymses.center center = self.pp_params.pymses.center
im_extent = [ im_extent = [
(-self._radius + center[0]), (-self._radius + center[0]),
(self._radius + center[0]), (self._radius + center[0]),
@@ -103,13 +108,13 @@ class PostProcessor(HDF5Container):
ax_v = self._axes_v[ax_los] ax_v = self._axes_v[ax_los]
self._cam[ax_los] = Camera( self._cam[ax_los] = Camera(
center=pp_params.pymses.center, center=self.pp_params.pymses.center,
line_of_sight_axis=ax_los, line_of_sight_axis=ax_los,
region_size=[2.0 * self._radius, 2.0 * self._radius], region_size=[2.0 * self._radius, 2.0 * self._radius],
distance=self._radius, distance=self._radius,
far_cut_depth=self._radius, far_cut_depth=self._radius,
up_vector=ax_v, up_vector=ax_v,
map_max_size=pp_params.pymses.map_size, map_max_size=self.pp_params.pymses.map_size,
) )
self.close() self.close()