In [ ]:
%matplotlib inline
import matplotlib.pyplot as plt
import numpy as np
import openmhd
# dummy index
vx=0;vy=1;vz=2;pr=3;ro=4;bx=5;by=6;bz=7;ps=8
# reading the data ...
x,y,t,data = openmhd.data_read("data/field-%05d.dat" % 20)
# reading the data (partial domain: [ix1,ix2] x [jx1,jx2])
# x,y,t,data = openmhd.data_read("data/field-%05d.dat" % 20,ix1=0,ix2=100,jx1=11)
# clearing the current figure, if any
plt.clf()
# extent: [left, right, bottom, top]
extent=[x[0],x[-1],y[0],y[-1]]
# 2D plot (vmin/mymin: minimum value, vmax/mymax: max value)
# Note: ().T is necessary for 2-D plot routines (imshow/pcolormesh...)
tmp = np.ndarray((x.size,y.size),np.double)
tmp[:,:] = data[:,:,pr]
mymax = max(tmp.max(), -tmp.min()) if( tmp.max() > 0.0 ) else 0.0
mymin = min(tmp.min(), -tmp.max()) if( tmp.min() < 0.0 ) else 0.0
myimg = plt.imshow(tmp.T,origin='lower',vmin=mymin,vmax=mymax,cmap='jet',extent=extent)
# image operations (e.g. colormaps)
# myimg.set_cmap('jet')
# myimg.set_cmap('RdBu_r') # colortable(70,/reverse) in IDL
# myimg.set_cmap('seismic')
# myimg.set_cmap('gnuplot2')
# useful options
# plt.grid()
plt.xlabel("X",size=16)
plt.ylabel("Y",size=16)
plt.title('Pressure (t = %6.1f)' % t, size=20)
# colorbar
plt.colorbar()
# preparing Vector potential (Az)
az = np.ndarray((x.size,y.size),np.double)
fx = 0.5*(x[1]-x[0])
fy = 0.5*(y[1]-y[0])
az[0,0] = (fy*data[0,0,bx] - fx*data[0,0,by])
for j in range(1,y.size):
az[0,j] = az[0,j-1] + fy*(data[0,j-1,bx]+data[0,j,bx])
for i in range(1,x.size):
az[i,:] = az[i-1,:] - fx*(data[i-1,:,by]+data[i,:,by])
# contour of Az = magnetic field lines
plt.contour(az.T,extent=extent,colors='w',linestyles='solid')
# plot
plt.show()
# image file
# plt.savefig('output.png', dpi=144)
# end
In [ ]:
%matplotlib inline
import matplotlib.pyplot as plt
import numpy as np
import openmhd
import os # for movies
# dummy index
vx=0;vy=1;vz=2;pr=3;ro=4;bx=5;by=6;bz=7;ps=8
# ---- for movies --------------------
# NOTE: "import os" is necessary
if not os.path.exists('movie'):
print('This routine outputs image files in "./movie/".')
os.mkdir('./movie/')
# ---- loop for movies ---------------
for ii in range(0,41):
# reading the data ...
x,y,t,data = openmhd.data_read("data/field-%05d.dat" % ii)
# clearing the current figure, if any
plt.clf()
# extent: [left, right, bottom, top]
extent=[x[0],x[-1],y[0],y[-1]]
# 2D plot (vmin/mymin: minimum value, vmax/mymax: max value)
# Note: ().T is necessary for 2-D plot routines (imshow/pcolormesh...)
tmp = np.ndarray((x.size,y.size),np.double)
tmp[:,:] = data[:,:,pr]
mymax = max(tmp.max(), -tmp.min()) if( tmp.max() > 0.0 ) else 0.0
mymin = min(tmp.min(), -tmp.max()) if( tmp.min() < 0.0 ) else 0.0
myimg = plt.imshow(tmp.T,origin='lower',vmin=mymin,vmax=mymax,cmap='jet',extent=extent)
plt.xlabel("X",size=16)
plt.ylabel("Y",size=16)
plt.title('Pressure (t = %6.1f)' % t, size=20)
# colorbar
plt.colorbar()
# preparing Vector potential (Az)
az = np.ndarray((x.size,y.size),np.double)
fx = 0.5*(x[1]-x[0])
fy = 0.5*(y[1]-y[0])
az[0,0] = (fy*data[0,0,bx] - fx*data[0,0,by])
for j in range(1,y.size):
az[0,j] = az[0,j-1] + fy*(data[0,j-1,bx]+data[0,j,bx])
for i in range(1,x.size):
az[i,:] = az[i-1,:] - fx*(data[i-1,:,by]+data[i,:,by])
# contour of Az = magnetic field lines
plt.contour(az.T,extent=extent,colors='w',linestyles='solid')
# plot
# plt.show()
# image file
plt.savefig('movie/output-%05d.png' % ii, dpi=144)
# ---- for movies --------------------
plt.show()
print()
print('The image files should be found in "./movie/".')
# end
In [ ]:
# Interactive version by jupyter-notebook / ipywidgets
# To use it, please install jupyter and then activate widgetsnbextension.
# $ pip3 install jupyter
# $ [ pip3 install ipywidgets ]
# $ jupyter nbextension enable --py widgetsnbextension
# Then one can run this sample
# $ jupyter-notebook plot.ipynb
%matplotlib inline
import numpy as np
import matplotlib.pyplot as plt
import openmhd
import glob
from ipywidgets import interact
# dummy index
vx=0;vy=1;vz=2;pr=3;ro=4;bx=5;by=6;bz=7;ps=8
labels={"Vx":0,"Vy":1,"Vz":2,"Pressure":3,"Density":4,"Bx":5,"By":6,"Bz":7,"Psi":8}
cmaps=['jet','RdBu_r','seismic','gnuplot2']
datalist = sorted(glob.glob('data/field-?????.dat')) # filelist (regular expression)
@interact (inputdata=datalist,mylabel=labels,mycmap=cmaps)
def plot(inputdata='data/field-00020.dat',mylabel=3,mycmap='jet'):
# reading the data ...
x,y,t,data = openmhd.data_read(inputdata)
# reading the data (partial domain: [ix1,ix2] x [jx1,jx2])
# x,y,t,data = openmhd.data_read("data/field-%05d.dat" % 20,ix1=0,ix2=100,jx1=11)
print(inputdata)
# extent: [left, right, bottom, top]
extent=[x[0],x[-1],y[0],y[-1]]
# 2D plot (vmin/mymin: minimum value, vmax/mymax: max value)
# Note: ().T is necessary for 2-D plot routines (imshow/pcolormesh...)
tmp = np.ndarray((x.size,y.size),np.double)
tmp[:,:] = data[:,:,mylabel]
mymax = max(tmp.max(), -tmp.min()) if( tmp.max() > 0.0 ) else 0.0
mymin = min(tmp.min(), -tmp.max()) if( tmp.min() < 0.0 ) else 0.0
myimg = plt.imshow(tmp.T,origin='lower',vmin=mymin,vmax=mymax,cmap=mycmap,extent=extent)
# useful options
# plt.grid()
plt.xlabel("X",size=16)
plt.ylabel("Y",size=16)
plt.title('%s (t = %6.1f)' % ( (list(labels.keys()))[mylabel], t), size=20)
# colorbar
plt.colorbar()
# preparing Vector potential (Az)
az = np.ndarray((x.size,y.size),np.double)
fx = 0.5*(x[1]-x[0])
fy = 0.5*(y[1]-y[0])
az[0,0] = (fy*data[0,0,bx] - fx*data[0,0,by])
for j in range(1,y.size):
az[0,j] = az[0,j-1] + fy*(data[0,j-1,bx]+data[0,j,bx])
for i in range(1,x.size):
az[i,:] = az[i-1,:] - fx*(data[i-1,:,by]+data[i,:,by])
plt.contour(az.T,extent=extent,colors='w',linestyles='solid')
In [ ]:
# This sample generates an animated GIF file (animation_pillow.gif) with help from Pillow.
# It may be necessary to install pillow:
# $ pip3 install pillow
import glob
from PIL import Image
images = []
for png in sorted(glob.glob("movie/output-*.png")):
im = Image.open(png)
images.append(im)
if len(images) > 1:
images[0].save('animation_pillow.gif', save_all=True, append_images=images[1:], optimize=False, duration=100, loop=0)
else:
print("No output.")
In [ ]:
# This sample generates a MP4 file (animation_cv2.mp4) with help from OpenCV.
# It may be necessary to install opencv-python:
# $ pip3 install opencv-python
import glob
import cv2
from IPython.display import Video
images = []
for png in sorted(glob.glob("movie/output-*.png")):
im = cv2.imread(png)
height, width, layers = im.shape
images.append(im)
out = cv2.VideoWriter("animation_cv2.mp4", cv2.VideoWriter_fourcc(*'MP4V'), 10.0, (width, height))
for i in range(len(images)):
out.write(images[i])
out.release()
Video("animation_cv2.mp4")
In [ ]: