DAMASK_EICMD/python/damask/_colormap.py

806 lines
24 KiB
Python
Raw Permalink Blame History

This file contains ambiguous Unicode characters

This file contains Unicode characters that might be confused with other characters. If you think that this is intentional, you can safely ignore this warning. Use the Escape button to reveal them.

import os
import json
import functools
import colorsys
from typing import Union, TextIO
import numpy as np
import scipy.interpolate as interp
import matplotlib as mpl
if os.name == 'posix' and 'DISPLAY' not in os.environ:
mpl.use('Agg')
import matplotlib.pyplot as plt
from matplotlib import cm
from PIL import Image
from ._typehints import FloatSequence, FileHandle
from . import util
from . import Table
_EPS = 216./24389.
_KAPPA = 24389./27.
_REF_WHITE = np.array([.95047, 1.00000, 1.08883]) # Observer = 2, Illuminant = D65
# ToDo (if needed)
# - support alpha channel (paraview/ASCII/input)
# - support NaN color (paraview)
class Colormap(mpl.colors.ListedColormap):
"""
Enhance matplotlib colormap functionality to be used within DAMASK.
Colors are internally stored as R(ed) G(green) B(lue) values.
The colormap can be used in matplotlib, seaborn, etc., or can
exported to file for external use.
References
----------
K. Moreland, Proceedings of the 5th International Symposium on Advances in Visual Computing, 2009
https://doi.org/10.1007/978-3-642-10520-3_9
P. Eisenlohr et al., International Journal of Plasticity 46:3753, 2013
https://doi.org/10.1016/j.ijplas.2012.09.012
Matplotlib colormaps overview
https://matplotlib.org/stable/tutorials/colors/colormaps.html
"""
def __eq__(self,
other: object) -> bool:
"""
Return self==other.
Test equality of other.
"""
if not isinstance(other, Colormap):
return NotImplemented
return len(self.colors) == len(other.colors) \
and bool(np.all(self.colors == other.colors))
def __add__(self,
other: 'Colormap') -> 'Colormap':
"""
Return self+other.
Concatenate.
"""
return Colormap(np.vstack((self.colors,other.colors)),
f'{self.name}+{other.name}')
def __iadd__(self,
other: 'Colormap') -> 'Colormap':
"""
Return self+=other.
Concatenate (in-place).
"""
return self.__add__(other)
def __mul__(self,
factor: int) -> 'Colormap':
"""
Return self*other.
Repeat.
"""
return Colormap(np.vstack([self.colors]*factor),f'{self.name}*{factor}')
def __imul__(self,
factor: int) -> 'Colormap':
"""
Return self*=other.
Repeat (in-place).
"""
return self.__mul__(factor)
def __invert__(self) -> 'Colormap':
"""
Return ~self.
Reverse.
"""
return self.reversed()
def __repr__(self) -> str:
"""
Return repr(self).
Show as matplotlib figure.
"""
fig = plt.figure(self.name,figsize=(5,.5))
ax1 = fig.add_axes([0, 0, 1, 1])
ax1.set_axis_off()
ax1.imshow(np.linspace(0,1,self.N).reshape(1,-1),
aspect='auto', cmap=self, interpolation='nearest')
plt.show(block=False)
return f'Colormap: {self.name}'
@staticmethod
def from_range(low: FloatSequence,
high: FloatSequence,
name: str = 'DAMASK colormap',
N: int = 256,
model: str = 'rgb') -> 'Colormap':
"""
Create a perceptually uniform colormap between given (inclusive) bounds.
Parameters
----------
low : sequence of float, len (3)
Color definition for minimum value.
high : sequence of float, len (3)
Color definition for maximum value.
name : str, optional
Name of the colormap. Defaults to 'DAMASK colormap'.
N : int, optional
Number of color quantization levels. Defaults to 256.
model : {'rgb', 'hsv', 'hsl', 'xyz', 'lab', 'msh'}
Color model used for input color definitions. Defaults to 'rgb'.
The available color models are:
- 'rgb': Red Green Blue.
- 'hsv': Hue Saturation Value.
- 'hsl': Hue Saturation Luminance.
- 'xyz': CIE Xyz.
- 'lab': CIE Lab.
- 'msh': Msh (for perceptual uniform interpolation).
Returns
-------
new : damask.Colormap
Colormap within given bounds.
Examples
--------
>>> import damask
>>> damask.Colormap.from_range((0,0,1),(0,0,0),'blue_to_black')
"""
toMsh = dict(
rgb=Colormap._rgb2msh,
hsv=Colormap._hsv2msh,
hsl=Colormap._hsl2msh,
xyz=Colormap._xyz2msh,
lab=Colormap._lab2msh,
msh=lambda x:x,
)
if model.lower() not in toMsh:
raise ValueError(f'invalid color model "{model}"')
low_high = np.vstack((low,high)).astype(float)
out_of_bounds = np.bool_(False)
if model.lower() == 'rgb':
out_of_bounds = np.any(low_high<0) or np.any(low_high>1)
elif model.lower() == 'hsv':
out_of_bounds = np.any(low_high<0) or np.any(low_high>[360,1,1])
elif model.lower() == 'hsl':
out_of_bounds = np.any(low_high<0) or np.any(low_high>[360,1,1])
elif model.lower() == 'lab':
out_of_bounds = np.any(low_high[:,0]<0)
if out_of_bounds:
raise ValueError(f'{model.upper()} colors {low_high[0]} | {low_high[1]} are out of bounds')
low_,high_ = map(toMsh[model.lower()],low_high)
msh = map(functools.partial(Colormap._interpolate_msh,low=low_,high=high_),np.linspace(0,1,N))
rgb = np.array(list(map(Colormap._msh2rgb,msh)))
return Colormap(rgb,name=name)
@staticmethod
def from_predefined(name: str,
N: int = 256) -> 'Colormap':
"""
Select from a set of predefined colormaps.
Predefined colormaps (Colormap.predefined) include
native matplotlib colormaps and common DAMASK colormaps.
Parameters
----------
name : str
Name of the colormap.
N : int, optional
Number of color quantization levels. Defaults to 256.
This parameter is not used for matplotlib colormaps
that are of type `ListedColormap`.
Returns
-------
new : damask.Colormap
Predefined colormap.
Examples
--------
>>> import damask
>>> damask.Colormap.from_predefined('strain')
"""
try:
# matplotlib presets
colormap = cm.__dict__[name]
return Colormap(np.array(list(map(colormap,np.linspace(0,1,N)))
if isinstance(colormap,mpl.colors.LinearSegmentedColormap) else
colormap.colors),
name=name)
except KeyError:
# DAMASK presets
definition = Colormap._predefined_DAMASK[name]
return Colormap.from_range(definition['low'],definition['high'],name,N)
def at(self,
fraction : Union[float,FloatSequence]) -> np.ndarray:
"""
Interpolate color at fraction.
Parameters
----------
fraction : float or sequence of float
Fractional coordinate(s) to evaluate Colormap at.
Returns
-------
color : numpy.ndarray, shape(...,4)
RGBA values of interpolated color(s).
Examples
--------
>>> import damask
>>> cmap = damask.Colormap.from_predefined('gray')
>>> cmap.at(0.5)
array([0.5, 0.5, 0.5, 1. ])
>>> 'rgb({},{},{})'.format(*cmap.at(0.5))
'rgb(0.5,0.5,0.5)'
"""
return interp.interp1d(np.linspace(0,1,self.N),
self.colors,
axis=0,
assume_sorted=True)(fraction)
def shade(self,
field: np.ndarray,
bounds: FloatSequence = None,
gap: float = None) -> Image:
"""
Generate PIL image of 2D field using colormap.
Parameters
----------
field : numpy.ndarray, shape (:,:)
Data to be shaded.
bounds : sequence of float, len (2), optional
Value range (left,right) spanned by colormap.
gap : field.dtype, optional
Transparent value. NaN will always be rendered transparent.
Returns
-------
PIL.Image
RGBA image of shaded data.
"""
mask = np.logical_not(np.isnan(field) if gap is None else \
np.logical_or (np.isnan(field), field == gap)) # mask NaN (and gap if present)
l,r = (field[mask].min(),field[mask].max()) if bounds is None else \
(bounds[0],bounds[1])
if abs(delta := r-l) * 1e8 <= (avg := 0.5*abs(r+l)): # delta is similar to numerical noise
l,r = (l-0.5*avg*np.sign(delta),r+0.5*avg*np.sign(delta)) # extend range to have actual data centered within
return Image.fromarray(
(np.dstack((
self.colors[(np.round(np.clip((field-l)/delta,0.0,1.0)*(self.N-1))).astype(np.uint16),:3],
mask.astype(float)
)
)*255
).astype(np.uint8),
mode='RGBA')
def reversed(self,
name: str = None) -> 'Colormap':
"""
Reverse.
Parameters
----------
name : str, optional
Name of the reversed colormap.
Defaults to parent colormap name + '_r'.
Returns
-------
damask.Colormap
Reversed colormap.
Examples
--------
>>> import damask
>>> damask.Colormap.from_predefined('stress').reversed()
"""
rev = super().reversed(name)
return Colormap(np.array(rev.colors),rev.name[:-4] if rev.name.endswith('_r_r') else rev.name)
def _get_file_handle(self,
fname: Union[FileHandle, None],
suffix: str = '') -> TextIO:
"""
Provide file handle.
Parameters
----------
fname : file, str, pathlib.Path, or None
Name or handle of file.
If None, colormap name + suffix.
suffix: str, optional
Extension to use for colormap file.
Returns
-------
f : file object
File handle with write access.
"""
return util.open_text(self.name.replace(' ','_')+suffix if fname is None else fname, 'w')
def save_paraview(self,
fname: FileHandle = None):
"""
Save as JSON file for use in Paraview.
Parameters
----------
fname : file, str, or pathlib.Path, optional
File to store results. Defaults to colormap name + '.json'.
"""
colors = []
for i,c in enumerate(np.round(self.colors,6).tolist()):
colors+=[i]+c
out = [{
'Creator':util.execution_stamp('Colormap'),
'ColorSpace':'RGB',
'Name':self.name,
'DefaultMap':True,
'RGBPoints':colors
}]
fhandle = self._get_file_handle(fname,'.json')
json.dump(out,fhandle,indent=4)
fhandle.write('\n')
def save_ASCII(self,
fname: FileHandle = None):
"""
Save as ASCII file.
Parameters
----------
fname : file, str, or pathlib.Path, optional
File to store results. Defaults to colormap name + '.txt'.
"""
labels = {'RGBA':4} if self.colors.shape[1] == 4 else {'RGB': 3}
t = Table(labels,self.colors,f'Creator: {util.execution_stamp("Colormap")}')
t.save(self._get_file_handle(fname,'.txt'))
def save_GOM(self, fname: FileHandle = None):
"""
Save as ASCII file for use in GOM Aramis.
Parameters
----------
fname : file, str, or pathlib.Path, optional
File to store results. Defaults to colormap name + '.legend'.
"""
# ToDo: test in GOM
GOM_str = '1 1 {name} 9 {name} '.format(name=self.name.replace(" ","_")) \
+ '0 1 0 3 0 0 -1 9 \\ 0 0 0 255 255 255 0 0 255 ' \
+ f'30 NO_UNIT 1 1 64 64 64 255 1 0 0 0 0 0 0 3 0 {self.N}' \
+ ' '.join([f' 0 {c[0]} {c[1]} {c[2]} 255 1' for c in reversed((self.colors*255).astype(np.int64))]) \
+ '\n'
self._get_file_handle(fname,'.legend').write(GOM_str)
def save_gmsh(self,
fname: FileHandle = None):
"""
Save as ASCII file for use in gmsh.
Parameters
----------
fname : file, str, or pathlib.Path, optional
File to store results. Defaults to colormap name + '.msh'.
"""
# ToDo: test in gmsh
gmsh_str = 'View.ColorTable = {\n' \
+'\n'.join([f'{c[0]},{c[1]},{c[2]},' for c in self.colors[:,:3]*255]) \
+'\n}\n'
self._get_file_handle(fname,'.msh').write(gmsh_str)
@staticmethod
def _interpolate_msh(frac: float,
low: np.ndarray,
high: np.ndarray) -> np.ndarray:
"""
Interpolate in Msh color space.
This interpolation gives a perceptually uniform colormap.
References
----------
https://www.kennethmoreland.com/color-maps/ColorMapsExpanded.pdf
https://www.kennethmoreland.com/color-maps/diverging_map.py
"""
def rad_diff(a,b):
return abs(a[2]-b[2])
def adjust_hue(msh_sat, msh_unsat):
"""If saturation of one of the two colors is much less than the other, hue of the less."""
if msh_sat[0] >= msh_unsat[0]:
return msh_sat[2]
hSpin = msh_sat[1]/np.sin(msh_sat[1])*np.sqrt(msh_unsat[0]**2.0-msh_sat[0]**2)/msh_sat[0]
if msh_sat[2] < - np.pi/3.0: hSpin *= -1.0
return msh_sat[2] + hSpin
lo = np.array(low)
hi = np.array(high)
if (lo[1] > 0.05 and hi[1] > 0.05 and rad_diff(lo,hi) > np.pi/3.0):
M_mid = max(lo[0],hi[0],88.0)
if frac < 0.5:
hi = np.array([M_mid,0.0,0.0])
frac *= 2.0
else:
lo = np.array([M_mid,0.0,0.0])
frac = 2.0*frac - 1.0
if lo[1] < 0.05 < hi[1]:
lo[2] = adjust_hue(hi,lo)
elif hi[1] < 0.05 < lo[1]:
hi[2] = adjust_hue(lo,hi)
return (1.0 - frac) * lo + frac * hi
_predefined_mpl= {'Perceptually Uniform Sequential': [
'viridis', 'plasma', 'inferno', 'magma', 'cividis'],
'Sequential': [
'Greys', 'Purples', 'Blues', 'Greens', 'Oranges', 'Reds',
'YlOrBr', 'YlOrRd', 'OrRd', 'PuRd', 'RdPu', 'BuPu',
'GnBu', 'PuBu', 'YlGnBu', 'PuBuGn', 'BuGn', 'YlGn'],
'Sequential (2)': [
'binary', 'gist_yarg', 'gist_gray', 'gray', 'bone', 'pink',
'spring', 'summer', 'autumn', 'winter', 'cool', 'Wistia',
'hot', 'afmhot', 'gist_heat', 'copper'],
'Diverging': [
'PiYG', 'PRGn', 'BrBG', 'PuOr', 'RdGy', 'RdBu',
'RdYlBu', 'RdYlGn', 'Spectral', 'coolwarm', 'bwr', 'seismic'],
'Cyclic': ['twilight', 'twilight_shifted', 'hsv'],
'Qualitative': [
'Pastel1', 'Pastel2', 'Paired', 'Accent',
'Dark2', 'Set1', 'Set2', 'Set3',
'tab10', 'tab20', 'tab20b', 'tab20c'],
'Miscellaneous': [
'flag', 'prism', 'ocean', 'gist_earth', 'terrain', 'gist_stern',
'gnuplot', 'gnuplot2', 'CMRmap', 'cubehelix', 'brg',
'gist_rainbow', 'rainbow', 'jet', 'nipy_spectral', 'gist_ncar']}
_predefined_DAMASK = {'orientation': {'low': [0.933334,0.878432,0.878431], # noqa
'high': [0.250980,0.007843,0.000000]},
'strain': {'low': [0.941177,0.941177,0.870588],
'high': [0.266667,0.266667,0.000000]},
'stress': {'low': [0.878432,0.874511,0.949019],
'high': [0.000002,0.000000,0.286275]}}
predefined = dict(**{'DAMASK':list(_predefined_DAMASK)},**_predefined_mpl)
@staticmethod
def _hsv2rgb(hsv: np.ndarray) -> np.ndarray:
"""
Hue Saturation Value to Red Green Blue.
Parameters
----------
hsv : numpy.ndarray, shape (3)
HSV values.
Returns
-------
rgb : numpy.ndarray, shape (3)
RGB values.
"""
return np.array(colorsys.hsv_to_rgb(hsv[0]/360.,hsv[1],hsv[2]))
@staticmethod
def _rgb2hsv(rgb: np.ndarray) -> np.ndarray:
"""
Red Green Blue to Hue Saturation Value.
Parameters
----------
rgb : numpy.ndarray, shape (3)
RGB values.
Returns
-------
hsv : numpy.ndarray, shape (3)
HSV values.
"""
h,s,v = colorsys.rgb_to_hsv(rgb[0],rgb[1],rgb[2])
return np.array([h*360,s,v])
@staticmethod
def _hsl2rgb(hsl: np.ndarray) -> np.ndarray:
"""
Hue Saturation Luminance to Red Green Blue.
Parameters
----------
hsl : numpy.ndarray, shape (3)
HSL values.
Returns
-------
rgb : numpy.ndarray, shape (3)
RGB values.
"""
return np.array(colorsys.hls_to_rgb(hsl[0]/360.,hsl[2],hsl[1]))
@staticmethod
def _rgb2hsl(rgb: np.ndarray) -> np.ndarray:
"""
Red Green Blue to Hue Saturation Luminance.
Parameters
----------
rgb : numpy.ndarray, shape (3)
RGB values.
Returns
-------
hsl : numpy.ndarray, shape (3)
HSL values.
"""
h,l,s = colorsys.rgb_to_hls(rgb[0],rgb[1],rgb[2])
return np.array([h*360,s,l])
@staticmethod
def _xyz2rgb(xyz: np.ndarray) -> np.ndarray:
"""
CIE Xyz to Red Green Blue.
Parameters
----------
xyz : numpy.ndarray, shape (3)
CIE Xyz values.
Returns
-------
rgb : numpy.ndarray, shape (3)
RGB values.
References
----------
https://www.easyrgb.com/en/math.php
"""
rgb_lin = np.dot(np.array([
[ 3.240969942,-1.537383178,-0.498610760],
[-0.969243636, 1.875967502, 0.041555057],
[ 0.055630080,-0.203976959, 1.056971514]
]),xyz)
with np.errstate(invalid='ignore'):
rgb = np.where(rgb_lin>0.0031308,rgb_lin**(1.0/2.4)*1.0555-0.0555,rgb_lin*12.92)
return np.clip(rgb,0.,1.)
@staticmethod
def _rgb2xyz(rgb: np.ndarray) -> np.ndarray:
"""
Red Green Blue to CIE Xyz.
Parameters
----------
rgb : numpy.ndarray, shape (3)
RGB values.
Returns
-------
xyz : numpy.ndarray, shape (3)
CIE Xyz values.
References
----------
https://www.easyrgb.com/en/math.php
"""
rgb_lin = np.where(rgb>0.04045,((rgb+0.0555)/1.0555)**2.4,rgb/12.92)
return np.dot(np.array([
[0.412390799,0.357584339,0.180480788],
[0.212639006,0.715168679,0.072192315],
[0.019330819,0.119194780,0.950532152]
]),rgb_lin)
@staticmethod
def _lab2xyz(lab: np.ndarray,
ref_white: np.ndarray = _REF_WHITE) -> np.ndarray:
"""
CIE Lab to CIE Xyz.
Parameters
----------
lab : numpy.ndarray, shape (3)
CIE lab values.
ref_white : numpy.ndarray, shape (3)
Reference white, default value is the standard 2° observer for D65.
Returns
-------
xyz : numpy.ndarray, shape (3)
CIE Xyz values.
References
----------
http://www.brucelindbloom.com/index.html?Eqn_Lab_to_XYZ.html
"""
f_x = (lab[0]+16.)/116. + lab[1]/500.
f_z = (lab[0]+16.)/116. - lab[2]/200.
return np.array([
f_x**3. if f_x**3. > _EPS else (116.*f_x-16.)/_KAPPA,
((lab[0]+16.)/116.)**3 if lab[0]>_KAPPA*_EPS else lab[0]/_KAPPA,
f_z**3. if f_z**3. > _EPS else (116.*f_z-16.)/_KAPPA
])*ref_white
@staticmethod
def _xyz2lab(xyz: np.ndarray,
ref_white: np.ndarray = _REF_WHITE) -> np.ndarray:
"""
CIE Xyz to CIE Lab.
Parameters
----------
xyz : numpy.ndarray, shape (3)
CIE Xyz values.
ref_white : numpy.ndarray, shape (3)
Reference white, default value is the standard 2° observer for D65.
Returns
-------
lab : numpy.ndarray, shape (3)
CIE lab values.
References
----------
http://www.brucelindbloom.com/index.html?Eqn_Lab_to_XYZ.html
"""
f = np.where(xyz/ref_white > _EPS,(xyz/ref_white)**(1./3.),(_KAPPA*xyz/ref_white+16.)/116.)
return np.array([
116.0 * f[1] - 16.0,
500.0 * (f[0] - f[1]),
200.0 * (f[1] - f[2])
])
@staticmethod
def _lab2msh(lab: np.ndarray) -> np.ndarray:
"""
CIE Lab to Msh.
Parameters
----------
lab : numpy.ndarray, shape (3)
CIE lab values.
Returns
-------
msh : numpy.ndarray, shape (3)
Msh values.
References
----------
https://www.kennethmoreland.com/color-maps/ColorMapsExpanded.pdf
https://www.kennethmoreland.com/color-maps/diverging_map.py
"""
M = np.linalg.norm(lab)
return np.array([
M,
np.arccos(lab[0]/M) if M>1e-8 else 0.,
np.arctan2(lab[2],lab[1]) if M>1e-8 else 0.,
])
@staticmethod
def _msh2lab(msh: np.ndarray) -> np.ndarray:
"""
Msh to CIE Lab.
Parameters
----------
msh : numpy.ndarray, shape (3)
Msh values.
Returns
-------
lab : numpy.ndarray, shape (3)
CIE lab values.
References
----------
https://www.kennethmoreland.com/color-maps/ColorMapsExpanded.pdf
https://www.kennethmoreland.com/color-maps/diverging_map.py
"""
return np.array([
msh[0] * np.cos(msh[1]),
msh[0] * np.sin(msh[1]) * np.cos(msh[2]),
msh[0] * np.sin(msh[1]) * np.sin(msh[2])
])
@staticmethod
def _lab2rgb(lab: np.ndarray) -> np.ndarray:
return Colormap._xyz2rgb(Colormap._lab2xyz(lab))
@staticmethod
def _rgb2lab(rgb: np.ndarray) -> np.ndarray:
return Colormap._xyz2lab(Colormap._rgb2xyz(rgb))
@staticmethod
def _msh2rgb(msh: np.ndarray) -> np.ndarray:
return Colormap._lab2rgb(Colormap._msh2lab(msh))
@staticmethod
def _rgb2msh(rgb: np.ndarray) -> np.ndarray:
return Colormap._lab2msh(Colormap._rgb2lab(rgb))
@staticmethod
def _hsv2msh(hsv: np.ndarray) -> np.ndarray:
return Colormap._rgb2msh(Colormap._hsv2rgb(hsv))
@staticmethod
def _hsl2msh(hsl: np.ndarray) -> np.ndarray:
return Colormap._rgb2msh(Colormap._hsl2rgb(hsl))
@staticmethod
def _xyz2msh(xyz: np.ndarray) -> np.ndarray:
return Colormap._lab2msh(Colormap._xyz2lab(xyz))