Commit 280d4224 authored by Stephan Kuschel's avatar Stephan Kuschel
Browse files

Merge remote-tracking branch 'gh/master' into docs

parents 45f76e52 6efda939
name: run-tests
on: [push, pull_request]
jobs:
latest:
runs-on: ubuntu-latest
strategy:
matrix:
python-version: [3.5, 3.6, 3.7, 3.8]
steps:
- uses: actions/checkout@v2
- name: Set up Python ${{ matrix.python-version }}
uses: actions/setup-python@v1
with:
python-version: ${{ matrix.python-version }}
- name: Install dependencies
run: |
python -m pip install --upgrade pip
pip install -r pip-requirements.txt
- name: run tests
run: |
./run-tests.py
other-os:
runs-on: ${{ matrix.os }}
strategy:
matrix:
os: [macos-latest, windows-latest]
python-version: [3.x]
# have to copy steps from above, as anchors are currently
# not supported by github workflow (Jan 2020).
steps:
- uses: actions/checkout@v2
- name: Set up Python ${{ matrix.python-version }}
uses: actions/setup-python@v1
with:
python-version: ${{ matrix.python-version }}
- name: Install dependencies
run: |
python -m pip install --upgrade pip
pip install -r pip-requirements.txt
- name: run tests
run: |
./run-tests.py
# auto generated by setup.py
/examples/_openPMDdata
/_examplepictures
# Byte-compiled / optimized / DLL files
......
# Check this file before committing a new version:
# http://lint.travis-ci.org/
language: python
sudo: false
dist: xenial
notifications:
email:
on_success: change
on_failure: change
env:
- MATPLOTLIB_V=* NUMPY_V=* SCIPY_V=* CYTHON_V=* PYGMENTS_V="=*"
# - MATPLOTLIB_V=2.0.2 NUMPY_V=* SCIPY_V=* CYTHON_V=*
jobs:
fast_finish: true
# allow_failures:
# have this allowed as there is a bug in matplotlib 2.1.0
# - env: MATPLOTLIB_V=* NUMPY_V=* SCIPY_V=* CYTHON_V=*
include:
- stage: Fast test
python: 3.5
virtualenv:
system_site_packages: true
before_install: skip
install:
- pip install 'matplotlib<3.1'
- pip install -r pip-requirements.txt
before_script: skip
script: ./run-tests.py --fast
- stage: Tests
python: 2.7
env: MATPLOTLIB_V=2.2.3
- python: 3.5
env: MATPLOTLIB_V=3.0
- python: 3.6
- python: 3.7
# Ubuntu 14.04.5 LTS
# python 2.6 is explicitly NOT supported!
- name: "Ubuntu 14.04.5 LTS - Python 3.4"
python: 3.4
env: MATPLOTLIB_V=1.3.1 NUMPY_V=1.8.1 SCIPY_V=0.13.3 CYTHON_V=0.20.1 EXTRAPACK="libgfortran<3" PYGMENTS_V="<2.4"
# Ubuntu 16.04.3 LTS
- name: "Ubuntu 16.04.3 LTS - Python 2.7"
python: 2.7
env: MATPLOTLIB_V=1.5.1 NUMPY_V=1.11.0 SCIPY_V=0.17.0 CYTHON_V=0.23.4
- name: "Ubuntu 16.04.3 LTS - Python 3.5"
python: 3.5
env: MATPLOTLIB_V=1.5.1 NUMPY_V=1.11.0 SCIPY_V=0.17.0 CYTHON_V=0.23.4
# Ubuntu 18.04 LTS
- name: "Ubuntu 18.04 LTS - Python 2.7"
python: 2.7
env: MATPLOTLIB_V=2.1.1 NUMPY_V=1.13.3 SCIPY_V=0.19.1 CYTHON_V=0.26.1
- name: "Ubuntu 18.04 LTS - Python 3.6"
python: 3.6
env: MATPLOTLIB_V=2.1.1 NUMPY_V=1.13.3 SCIPY_V=0.19.1 CYTHON_V=0.26.1
# Install conda
before_install:
- wget https://repo.continuum.io/miniconda/Miniconda3-latest-Linux-x86_64.sh -O miniconda.sh;
- chmod +x miniconda.sh
- ./miniconda.sh -b -p $HOME/miniconda
- export PATH=$HOME/miniconda/bin:$PATH
- conda update --yes conda;
- conda create -yn pyenv python=$TRAVIS_PYTHON_VERSION atlas numpy=$NUMPY_V scipy=$SCIPY_V matplotlib=$MATPLOTLIB_V nose cython=$CYTHON_V h5py numexpr sphinx pygments$PYGMENTS_V $EXTRAPACK
- source activate pyenv
install:
- pip install -r pip-requirements.txt
- python setup.py install
before_script:
- uname -a
- free -m
- df -h
- ulimit -a
- python -V
- cython --version
# run tests
script:
- ./run-tests.py
......@@ -2,9 +2,33 @@
## current master
* New convenience method `Field.copy`
**Highlights**
* Parallelized implementation of `Field.map_coordinates`
* Added support for reading files written by the fbpic v0.13.0 code ( https://github.com/fbpic ). The fields can be accessed by `Er` and `Etheta`, which have been introduced to the fbpic data reader. Particles are saved in Cartesian coordinates, hence the interface does not change there.
* Reimplementation of `Field.fft`, see below.
**Incompatible adjustments to previous version**
* Reimplementation of `Field.fft`. This changes the phases of Fourier transforms in a way to make it more consistent. However, if your code depends on the phases, `Field.fft()` now has a parameter `old_behaviour` that can be used to switch back to the old behaviour.
* Indexing a field by a number (integer or float) will now remove the according axis altogether, instead of leaving behind a length-1 axis.
A new class `KeepDim` was introduced through which the old behaviour can still be used.
Behaviour of PostPic before this change:
```
field.shape == (x,y,z)
field[:, 0.0, :].shape == (x,1,z)
```
Using the new class `KeepDim`, it is possible to retain that behaviour in the new version:
```
field.shape == (x,y,z)
field[:, 0.0, :].shape == (x, z)
field[:, KeepDim(0.0), :].shape == (x,1,z)
```
**Other improvements and new features**
* New convenience method `Field.copy`
## v0.4
2019-01-14
......
......@@ -27,19 +27,26 @@ Additionally postpic can plot and label your plot automatically. This makes it e
Installation
------------
Postpic can be used with python2 and python3. However the usage of python3 is recommended.
Postpic can be used with `python2` and `python3`. However the usage of `python3` is recommended.
**Users** should install the latest version directly from github using pip (python package manager):
**Users** should install the latest version directly from github using `pip` (python package manager):
`pip install --user git+https://github.com/skuschel/postpic.git`
The latest *release* is also available ib the python package index [pypi](https://pypi.python.org/pypi/postpic/), thus it can be installed by using the python package manager pip:
Please note that, depending on your system python setup, `pip` may default to `python2`.
In that case you will need use `pip3` instead to make sure that postpic is installed for `python3`:
`pip install postpic`
`pip3 install --user git+https://github.com/skuschel/postpic.git`
The latest *release* is also available in the python package index [pypi](https://pypi.python.org/pypi/postpic/), thus it can be installed by using the python package manager pip:
`pip install --user postpic`
In both cases the `--user` option makes sure that the package can be installed without `root` privileges by installing it into your user profile.
**Developers** should clone the git repository (or their fork of it) and install it using
`python setup.py develop --user`
`pip install --user -e .`
This command will link the current folder to global python scope, such that changing the code will immediately update the installed package.
......
kspace-test-2d/
_openPMDdata/
#!/usr/bin/env python
# coding: utf-8
# In[1]:
import sys
def download(url, file):
import urllib3
import shutil
import os
if os.path.isfile(file):
return True
try:
urllib3.disable_warnings()
http = urllib3.PoolManager()
print('downloading {:} ...'.format(file))
with http.request('GET', url, preload_content=False) as r, open(file, 'wb') as out_file:
shutil.copyfileobj(r, out_file)
success = True
except urllib3.exceptions.MaxRetryError:
success = False
return success
def main():
import os
import os.path as osp
datadir = 'examples/kspace-test-2d'
tarball = osp.join(datadir, 'kspace-test-2d-resources.tar.gz')
tarball_url = 'https://blinne.net/files/postpic/kspace-test-2d-resources.tar.gz'
if not os.path.exists(datadir):
os.mkdir(datadir)
s = download(tarball_url, tarball)
if s:
import hashlib
chksum = hashlib.sha256(open(tarball, 'rb').read()).hexdigest()
if chksum != "54bbeae19e1f412ddd475f565c2193e92922ed8a22e7eb3ecb4d73b5cf193b24":
os.remove(tarball)
s = False
if not s:
print('Failed to Download example data. Skipping this example.')
return
import tarfile
tar = tarfile.open(tarball)
tar.extractall(datadir)
import matplotlib
matplotlib.use('Agg')
font = {'size' : 12}
matplotlib.rc('font', **font)
import copy
import postpic as pp
import numpy as np
import pickle
try:
import sdf
sdfavail = True
except ImportError:
sdfavail = False
# basic constants
micro = 1e-6
femto = 1e-15
c = pp.PhysicalConstants.c
# known parameters from simulation
lam = 0.5 * micro
k0 = 2*np.pi/lam
f0 = pp.PhysicalConstants.c/lam
pymajorver = str(sys.version_info[0])
if sdfavail:
pp.chooseCode("EPOCH")
dump = pp.readDump(osp.join(datadir, '0002.sdf'))
plotter = pp.plotting.plottercls(dump, autosave=False)
fields = dict()
for fc in ['Ey', 'Bz']:
fields[fc] = getattr(dump, fc)()
fields[fc].saveto(osp.join(datadir, '0002_'+fc+pymajorver), compressed=False)
t = dump.time()
dt = t/dump.timestep()
pickle.dump(dict(t=t, dt=dt), open(osp.join(datadir,'0002_meta'+pymajorver+'.pickle'), 'wb'))
else:
fields = dict()
for fc in ['Ey', 'Bz']:
fields[fc] = pp.Field.loadfrom(osp.join(datadir,'0002_{}.npz').format(fc+pymajorver))
meta = pickle.load(open(osp.join(datadir,'0002_meta'+pymajorver+'.pickle'), 'rb'))
t = meta['t']
dt = meta['dt']
plotter = pp.plotting.plottercls(None, autosave=False)
Ey = fields['Ey']
# grid spacing from field
dx = [ax.grid[1] - ax.grid[0] for ax in Ey.axes]
#w0 = 2*np.pi*f0
#w0 = pp.PhysicalConstants.c * k0
wn = np.pi/dt
omega_yee = pp.helper.omega_yee_factory(dx=dx, dt=dt)
wyee = omega_yee([k0,0])
w0 = pp.helper.omega_free([k0,0])
print('dx', dx)
print('t', t)
print('dt', dt)
print('wn', wn)
print('w0', w0)
print('wyee', wyee)
print('wyee/w0', wyee/w0)
print('wyee/wn', wyee/wn)
print('lam/dx[0]', lam/dx[0])
print('cos(1/2 wyee dt)', np.cos(1/2 * wyee * dt))
vg_yee = c*np.cos(k0*dx[0]/2.0)/np.sqrt(1-(c*dt/dx[0]*np.sin(k0*dx[0]/2.0))**2)
print('vg/c', vg_yee/c)
r = np.sqrt(1.0 - (pp.PhysicalConstants.c * dt)**2 * (1/dx[0]*np.sin(1/2.0*k0*dx[0]))**2)
print('r', r)
# In[2]:
omega_yee = pp.helper.omega_yee_factory(dx=Ey.spacing, dt=dt)
lin_int_response_omega = pp.helper._linear_interpolation_frequency_response(dt)
lin_int_response_k = pp.helper._linear_interpolation_frequency_response_on_k(lin_int_response_omega,
Ey.fft().axes, omega_yee)
lin_int_response_k_vac = pp.helper._linear_interpolation_frequency_response_on_k(lin_int_response_omega,
Ey.fft().axes, pp.helper.omega_free)
# In[3]:
_=plotter.plotField(Ey[:,0.0])
# In[4]:
kspace = dict()
component = "Ey"
# In[5]:
#key = 'default epoch map=yee, omega=vac'
key = 'linresponse map=yee, omega=vac'
if sdfavail:
kspace[key] = abs(dump.kspace_Ey(solver='yee'))
else:
kspace[key] = abs(pp.helper.kspace(component, fields=dict(Ey=fields['Ey'],
Bz=fields['Bz'].fft()/lin_int_response_k),
interpolation='fourier'))
# using the helper function `kspace_epoch_like` would yield same result:
# kspace[key] = abs(pp.helper.kspace_epoch_like(component, fields=fields, dt=dt, omega_func=omega_yee))
normalisation = 1.0/np.max(kspace[key].matrix)
kspace[key] *= normalisation
kspace[key].name = r'corrected $\vec{k}$-space, $\omega_0=c|\vec{k}|$'
kspace[key].unit = ''
# In[6]:
key = 'simple fft'
kspace[key] = abs(Ey.fft()) * normalisation
kspace[key].name = r'plain fft'
kspace[key].unit = ''
# In[7]:
key = 'fourier'
kspace[key] = abs(pp.helper.kspace(component,
fields=fields,
interpolation='fourier')
) * normalisation
kspace[key].name = u'naïve $\\vec{k}$-space, $\\omega_0=c|\\vec{k}|$'
kspace[key].unit = ''
# In[8]:
key = 'fourier yee'
kspace[key] = abs(pp.helper.kspace(component,
fields=fields, interpolation='fourier',
omega_func=omega_yee)
) * normalisation
kspace[key].name = u'naïve $\\vec{k}$-space, $\\omega_0=\\omega_\\mathrm{grid}$'
kspace[key].unit = ''
# In[9]:
key = 'linresponse map=yee, omega=yee'
kspace[key] = abs(pp.helper.kspace(component, fields=dict(Ey=fields['Ey'],
Bz=fields['Bz'].fft()/lin_int_response_k),
interpolation='fourier', omega_func=omega_yee)
) * normalisation
kspace[key].name = r'corrected $\vec{k}$-space, $\omega_0=\omega_\mathrm{grid}$'
kspace[key].unit = ''
# In[10]:
slices = [slice(360-120, 360+120), slice(120, 121)]
# In[11]:
keys = ['simple fft',
'fourier yee',
'fourier',
'linresponse map=yee, omega=yee',
'linresponse map=yee, omega=vac'
]
figure2 = plotter.plotFields1d(*[kspace[k][slices] for k in keys],
log10plot=True, ylim=(5e-17, 5))
figure2.set_figwidth(8)
figure2.set_figheight(6)
while figure2.axes[0].texts:
figure2.axes[0].texts[-1].remove()
figure2.axes[0].set_title('')
figure2.axes[0].set_ylabel(r'$|E_y(k_x,0,0)|\, [a. u.]$')
figure2.axes[0].set_xlabel(r'$k_x\,[m^{-1}]$')
figure2.tight_layout()
figure2.savefig(osp.join(datadir, 'gaussian-kspace.pdf'))
print("Integrated ghost peaks")
for k in keys:
I = kspace[k][:0.0,:].integrate().matrix
print(k, I)
if k == 'linresponse map=yee, omega=vac':
if I < 30000000.:
print('linresponse map=yee, omega=vac value is low enough: YES' )
else:
print('linresponse map=yee, omega=vac value is low enough: NO' )
print('Something is WRONG' )
sys.exit(1)
# In[13]:
if sdfavail:
kspace_ey = dump.kspace_Ey(solver='yee')
else:
kspace_ey = pp.helper.kspace_epoch_like(component, fields=fields, dt=dt, omega_func=omega_yee)
complex_ey = kspace_ey.fft()
envelope_ey_2d = abs(complex_ey)[:,:].squeeze()
try:
from skimage.restoration import unwrap_phase
phase_ey = complex_ey.replace_data( unwrap_phase(np.angle(complex_ey)) )
except ImportError:
phase_ey = complex_ey.replace_data( np.angle(complex_ey) )
phase_ey_2d = phase_ey[:,:].squeeze()
# In[14]:
ey = complex_ey.real[-0.1e-5:0.2e-5, :]
#ey = Ey
ey.name = r'$E_y$'
ey.unit = r'$\frac{\mathrm{V}}{\mathrm{m}}$'
figure = plotter.plotField2d(ey)
figure.set_figwidth(6)
figure.axes[0].set_title(r'')#$\Re E_y\ [\frac{\mathrm{V}}{\mathrm{m}}]$')
figure.axes[0].set_xlabel(u'$x\\,[µm]$')
figure.axes[0].set_ylabel(u'$y\\,[µm]$')
import matplotlib.ticker as ticker
ticks_x = ticker.FuncFormatter(lambda x, pos: '{0:g}'.format(x/1e-6))
figure.axes[0].xaxis.set_major_formatter(ticks_x)
figure.axes[0].yaxis.set_major_formatter(ticks_x)
try:
figure.axes[0].images[0].colorbar.remove()
except AttributeError:
pass
figure.colorbar(figure.axes[0].images[0], format='%6.0e', pad=0.15, label=r'$\Re E_y\ [\frac{\mathrm{V}}{\mathrm{m}}]$')
axes2 = figure.axes[0].twinx()
axes2.set_ylabel(r'$\operatorname{Arg}(E_y)\, [\pi],\;|E_y|\, [a. u.]$')
env = abs(complex_ey[-0.1e-5:0.2e-5,0.0])
m = np.max(env)
env = (env/m*40/np.pi).squeeze()
p = phase_ey[-0.1e-5:0.2e-5,0.0].squeeze()/np.pi
_ = axes2.plot(env.grid, env.matrix, label=r'$|E_y|\, [a. u.]$')
_ = axes2.plot(p.grid, p.matrix, label=r'$\operatorname{Arg}(E_y)\, [\pi]$')
handles, labels = axes2.get_legend_handles_labels()
axes2.legend(handles, labels)
#figure.axes[0].set_title(r'$E_y\,[V/m]$')
figure.set_figwidth(6)
figure.set_figheight(6)
figure.tight_layout()
figure.savefig(osp.join(datadir, 'gaussian-env-arg.pdf'))
# In[ ]:
if __name__=='__main__':
main()
......@@ -27,6 +27,11 @@ import scipy as sp
import scipy.signal as sps
import collections
try:
from collections.abc import Iterable
except ImportError:
from collections import Iterable
def np_meshgrid(*args, **kwargs):
if len(args) == 0:
......@@ -50,9 +55,9 @@ def np_moveaxis(*args, **kwargs):
a, source, destination = args
# twice a quick implementation of numpy.numeric.normalize_axis_tuple
if not isinstance(source, collections.Iterable):
if not isinstance(source, Iterable):
source = (source,)
if not isinstance(destination, collections.Iterable):
if not isinstance(destination, Iterable):
destination = (destination,)
source = [s % a.ndim for s in source]
destination = [d % a.ndim for d in destination]
......
......@@ -14,7 +14,7 @@
# You should have received a copy of the GNU General Public License
# along with postpic. If not, see <http://www.gnu.org/licenses/>.
#
# Stephan Kuschel 2014
# Stephan Kuschel 2014-2019
# Alexander Blinne, 2017
"""
Field related routines.
......@@ -51,12 +51,12 @@ class FieldAnalyzer(object):
pass
# General interface for everything
def _createfieldfromdata(self, data, gridkey):
def _createfieldfromdata(self, data, gridkey, **kwargs):
ret = Field(np.float64(data))
self.setgridtofield(ret, gridkey)
self.setgridtofield(ret, gridkey, **kwargs)
return ret
def createfieldfromkey(self, key, gridkey=None):
def createfieldfromkey(self, key, gridkey=None, **kwargs):
'''
This method creates a Field object from the data identified by "key".
The Grid is also inferred from that key unless an alternate "gridkey"
......@@ -64,27 +64,39 @@ class FieldAnalyzer(object):
'''
if gridkey is None:
gridkey = key
ret = self._createfieldfromdata(self.data(key), gridkey)
ret = self._createfieldfromdata(self.data(key, **kwargs), gridkey, **kwargs)
ret.name = key
return ret
def getaxisobj(self, gridkey, axis):
def getaxisobj(self, gridkey, axis, **kwargs):
'''
returns an Axis object for the "axis" and the grid defined by "gridkey".