Skip to content
This repository

New colormap normalizations: sqrt, arcsinh #1780

Closed
wants to merge 5 commits into from

6 participants

Adam Ginsburg Damon McDougall Thomas A Caswell Till Stensitzki Phil Elson Varoquaux
Adam Ginsburg

(this may be best to regard as a feature request with sample code included; I doubt my code conforms to mpl standards)

Related to #1355, #1204, and #632, I think there should be more default normalizations available in matplotlib. This pull request includes two that I've made use of: a sqrt/n'th root normalization and an arcsinh normalization. The former is good for emphasizing values near 1 (e.g., ratios of like values) and the latter is useful for high dynamic range images where the low scales are to be emphasized (e.g., images of stars with nebulosity around them); it largely serves the same purposes as SymLogNorm but with a different shape.

In part, I often want to reproduce figures made by ds9, which includes a wider range of default normalizations.

If there's an easier way to accomplish what I'm doing - e.g. just wrapping a function in some Normalize subclass - please let me know; I'm not up to date on new features of matplotlib unless they're thoroughly exampled in the gallery.

Damon McDougall
Collaborator

A question to the other developers: In the interest of preventing a combinatorial explosion by creating a class for each elementary mathematical function, would it be a better idea to have a class that will take a callable in its constructor, and have it normalise the data with respect to the underlying callable? For example, the user could call something like CustomNorm(np.sqrt, vmin=1, vmax=10).

Adam Ginsburg
Damon McDougall
Collaborator

@keflavich I am in agreement with you. The general 'take a callable' approach is more robust to changing requirements.

Thomas A Caswell
Collaborator

+1 for this (this being the callable method).

Till Stensitzki

While preferable, It is not simple: the normalization classes also include the inverse function, else features like colorbar won't work right.

Varoquaux NelleV commented on the diff March 01, 2013
lib/matplotlib/colors.py
((4 lines not shown))
  1003
+
  1004
+class SqrtNorm(Normalize):
  1005
+    """
  1006
+    Normalize a given value to the 0-1 range on a square (or n'th) root scale
  1007
+    """
  1008
+    def __init__(self, vmin=None, vmax=None, clip=False, nthroot=2):
  1009
+        """
  1010
+        nthroot allows cube roots, fourth roots, etc.
  1011
+        """
  1012
+        self.vmin = vmin
  1013
+        self.vmax = vmax
  1014
+        self.clip = clip
  1015
+        self.nthroot = nthroot
  1016
+
  1017
+    def __call__(self, value, clip=None, midpoint=None):
  1018
+
1
Varoquaux Collaborator
NelleV added a note March 01, 2013

PEP8 compliance: there should only be one blank line here.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
lib/matplotlib/colors.py
((21 lines not shown))
  1020
+        if clip is None:
  1021
+            clip = self.clip
  1022
+
  1023
+        if cbook.iterable(value):
  1024
+            vtype = 'array'
  1025
+            val = ma.asarray(value).astype(np.float)
  1026
+        else:
  1027
+            vtype = 'scalar'
  1028
+            val = ma.array([value]).astype(np.float)
  1029
+
  1030
+        self.autoscale_None(val)
  1031
+        vmin, vmax = self.vmin, self.vmax
  1032
+
  1033
+        if vmin > vmax:
  1034
+            raise ValueError("minvalue must be less than or equal to maxvalue")
  1035
+        elif vmin==vmax:
1
Varoquaux Collaborator
NelleV added a note March 01, 2013

PEP8 compliance: there should be whitespaces around '=='.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Varoquaux NelleV commented on the diff March 01, 2013
lib/matplotlib/colors.py
((100 lines not shown))
  1099
+            #result = (ma.arcsinh(val)-np.arcsinh(vmin))/(np.arcsinh(vmax)-np.arcsinh(vmin))
  1100
+            result = ma.arcsinh(result/midpoint) / ma.arcsinh(1./midpoint)
  1101
+        if vtype == 'scalar':
  1102
+            result = result[0]
  1103
+        return result
  1104
+
  1105
+    def autoscale_None(self, A):
  1106
+        ' autoscale only None-valued vmin or vmax'
  1107
+        if self.vmin is None:
  1108
+            self.vmin = ma.min(A)
  1109
+        if self.vmax is None:
  1110
+            self.vmax = ma.max(A)
  1111
+        if self.vmid is None:
  1112
+            self.vmid = (self.vmax+self.vmin)/2.0
  1113
+
  1114
+
1
Varoquaux Collaborator
NelleV added a note March 01, 2013

PEP compliance: there should be only 2 blank lines here.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Varoquaux NelleV commented on the diff March 01, 2013
lib/matplotlib/colors.py
((36 lines not shown))
  1035
+        elif vmin==vmax:
  1036
+            return 0.0 * val
  1037
+        else:
  1038
+            if clip:
  1039
+                mask = ma.getmask(val)
  1040
+                val = ma.array(np.clip(val.filled(vmax), vmin, vmax),
  1041
+                                mask=mask)
  1042
+            result = (val-vmin) * (1.0/(vmax-vmin))
  1043
+            #result = (ma.arcsinh(val)-np.arcsinh(vmin))/(np.arcsinh(vmax)-np.arcsinh(vmin))
  1044
+            result = result**(1./self.nthroot)
  1045
+        if vtype == 'scalar':
  1046
+            result = result[0]
  1047
+        return result
  1048
+
  1049
+    def autoscale_None(self, A):
  1050
+        ' autoscale only None-valued vmin or vmax'
1
Varoquaux Collaborator
NelleV added a note March 01, 2013

It would be nice to have real docstrings here, ie triple double-quotes.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Phil Elson
Collaborator

In principle, I think I'm in favour of this change. It definitely needs some tests though (they don't need to be pictorial).

@keflavich - do you fancy trying to get this into v1.3?

Adam Ginsburg

Yes, sounds like a good idea. I'll work on it. I guess I should start with autopep8 and then work on tests. When's the 1.3 deadline

Varoquaux
Collaborator

@keflavich autopep8 really doesn't work well: you might have to recheck everything by hand.
On the other hand, the pep8 tool is very strict, but very useful.

Phil Elson
Collaborator

Closing as the PR seems to have stalled. I love the idea of a callable in a Norm constructor though - if somebody wants to take that forwards I'd be more than happy to review.

Cheers,

Phil Elson pelson closed this January 09, 2014
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
This page is out of date. Refresh to see the latest.
148  lib/matplotlib/colors.py
@@ -1000,7 +1000,140 @@ def autoscale_None(self, A):
1000 1000
             self.vmin = ma.min(A)
1001 1001
         if self.vmax is None:
1002 1002
             self.vmax = ma.max(A)
1003  
-            
  1003
+
  1004
+
  1005
+class SqrtNorm(Normalize):
  1006
+    """
  1007
+    Normalize a given value to the 0-1 range on a square (or n'th) root scale
  1008
+    """
  1009
+    def __init__(self, vmin=None, vmax=None, clip=False, nthroot=2):
  1010
+        """
  1011
+        nthroot allows cube roots, fourth roots, etc.
  1012
+        """
  1013
+        self.vmin = vmin
  1014
+        self.vmax = vmax
  1015
+        self.clip = clip
  1016
+        self.nthroot = nthroot
  1017
+
  1018
+    def __call__(self, value, clip=None, midpoint=None):
  1019
+
  1020
+        if clip is None:
  1021
+            clip = self.clip
  1022
+
  1023
+        if cbook.iterable(value):
  1024
+            vtype = 'array'
  1025
+            val = ma.asarray(value).astype(np.float)
  1026
+        else:
  1027
+            vtype = 'scalar'
  1028
+            val = ma.array([value]).astype(np.float)
  1029
+
  1030
+        self.autoscale_None(val)
  1031
+        vmin, vmax = self.vmin, self.vmax
  1032
+
  1033
+        if vmin > vmax:
  1034
+            raise ValueError("minvalue must be less than or equal to maxvalue")
  1035
+        elif vmin == vmax:
  1036
+            return 0.0 * val
  1037
+        else:
  1038
+            if clip:
  1039
+                mask = ma.getmask(val)
  1040
+                val = ma.array(np.clip(val.filled(vmax), vmin, vmax),
  1041
+                               mask=mask)
  1042
+            result = (val-vmin) * (1.0/(vmax-vmin))
  1043
+            #(arcsinh(val)-arcsinh(vmin))/(arcsinh(vmax)-arcsinh(vmin))
  1044
+            result = result**(1./self.nthroot)
  1045
+        if vtype == 'scalar':
  1046
+            result = result[0]
  1047
+        return result
  1048
+
  1049
+    def autoscale_None(self, A):
  1050
+        """ autoscale only None-valued vmin or vmax """
  1051
+        if self.vmin is None:
  1052
+            self.vmin = ma.min(A)
  1053
+        if self.vmax is None:
  1054
+            self.vmax = ma.max(A)
  1055
+
  1056
+    def inverse(self, value):
  1057
+        if not self.scaled():
  1058
+            raise ValueError("Not invertible until scaled")
  1059
+        vmin, vmax = self.vmin, self.vmax
  1060
+
  1061
+        if cbook.iterable(value):
  1062
+            val = ma.asarray(value)
  1063
+            # r = sqrt(v-vmin/(vmax-vmin))
  1064
+            # v = r**2 * (vmax-vmin) + vmin
  1065
+            return val**self.nthroot * (vmax-vmin) + vmin
  1066
+        else:
  1067
+            return value**self.nthroot * (vmax-vmin) + vmin
  1068
+
  1069
+
  1070
+class AsinhNorm(Normalize):
  1071
+    """
  1072
+    Normalize a range of values to 0-1 on an arcsinh scale (nice alternative to
  1073
+    SymLogNormalize)
  1074
+    """
  1075
+    def __init__(self, vmin=None, vmax=None, clip=False, vmid=None):
  1076
+        self.vmid = vmid
  1077
+        self.vmin = vmin
  1078
+        self.vmax = vmax
  1079
+        self.clip = clip
  1080
+
  1081
+    def __call__(self, value, clip=None, midpoint=None):
  1082
+
  1083
+        if clip is None:
  1084
+            clip = self.clip
  1085
+
  1086
+        if cbook.iterable(value):
  1087
+            vtype = 'array'
  1088
+            val = ma.asarray(value).astype(np.float)
  1089
+        else:
  1090
+            vtype = 'scalar'
  1091
+            val = ma.array([value]).astype(np.float)
  1092
+
  1093
+        self.autoscale_None(val)
  1094
+        vmin, vmax = self.vmin, self.vmax
  1095
+
  1096
+        vmid = self.vmid if self.vmid is not None else (vmax+vmin)/2.0
  1097
+
  1098
+        if midpoint is None:
  1099
+            midpoint = (vmid - vmin) / (vmax - vmin)
  1100
+
  1101
+        if vmin > vmax:
  1102
+            raise ValueError("minvalue must be less than or equal to maxvalue")
  1103
+        elif vmin == vmax:
  1104
+            return 0.0 * val
  1105
+        else:
  1106
+            if clip:
  1107
+                mask = ma.getmask(val)
  1108
+                val = ma.array(np.clip(val.filled(vmax), vmin, vmax),
  1109
+                               mask=mask)
  1110
+            result = (val-vmin) * (1.0/(vmax-vmin))
  1111
+            result = ma.arcsinh(result/midpoint) / ma.arcsinh(1./midpoint)
  1112
+        if vtype == 'scalar':
  1113
+            result = result[0]
  1114
+        return result
  1115
+
  1116
+    def autoscale_None(self, A):
  1117
+        ' autoscale only None-valued vmin or vmax'
  1118
+        if self.vmin is None:
  1119
+            self.vmin = ma.min(A)
  1120
+        if self.vmax is None:
  1121
+            self.vmax = ma.max(A)
  1122
+        if self.vmid is None:
  1123
+            self.vmid = (self.vmax+self.vmin)/2.0
  1124
+
  1125
+    def inverse(self, value):
  1126
+        if not self.scaled():
  1127
+            raise ValueError("Not invertible until scaled")
  1128
+        vmin, vmax = self.vmin, self.vmax
  1129
+
  1130
+        if cbook.iterable(value):
  1131
+            val = ma.asarray(value)
  1132
+            # r = arcsinh(v-vmin/(vmax-vmin) / midpoint) / arcsinh(1/midpoint)
  1133
+            # v = sinh(r * arcsinh(1/midpoint)) * midpoint * (vmax-vmin) + vmin
  1134
+            return np.sinh(val * np.arcsinh(1./midpoint)) * midpoint * (vmax-vmin) + vmin
  1135
+        else:
  1136
+            return np.sinh(value * np.arcsinh(1./midpoint)) * midpoint * (vmax-vmin) + vmin
1004 1137
 
1005 1138
 class SymLogNorm(Normalize):
1006 1139
     """
@@ -1039,7 +1172,7 @@ def __call__(self, value, clip=None):
1039 1172
         result, is_scalar = self.process_value(value)
1040 1173
         self.autoscale_None(result)
1041 1174
         vmin, vmax = self.vmin, self.vmax
1042  
-        
  1175
+
1043 1176
         if vmin > vmax:
1044 1177
             raise ValueError("minvalue must be less than or equal to maxvalue")
1045 1178
         elif vmin == vmax:
@@ -1048,7 +1181,7 @@ def __call__(self, value, clip=None):
1048 1181
             if clip:
1049 1182
                 mask = ma.getmask(result)
1050 1183
                 result = ma.array(np.clip(result.filled(vmax), vmin, vmax),
1051  
-                                mask=mask)
  1184
+                                  mask=mask)
1052 1185
             # in-place equivalent of above can be much faster
1053 1186
             resdat = self._transform(result.data)
1054 1187
             resdat -= self._lower
@@ -1057,8 +1190,8 @@ def __call__(self, value, clip=None):
1057 1190
         if is_scalar:
1058 1191
             result = result[0]
1059 1192
         return result
1060  
-    
1061  
-    def _transform(self, a):        
  1193
+
  1194
+    def _transform(self, a):
1062 1195
         """
1063 1196
         Inplace transformation.
1064 1197
         """
@@ -1069,7 +1202,7 @@ def _transform(self, a):
1069 1202
         a[masked] = log
1070 1203
         a[~masked] *= self._linscale_adj
1071 1204
         return a
1072  
-        
  1205
+
1073 1206
     def _inv_transform(self, a):
1074 1207
         """
1075 1208
         Inverse inplace Transformation.
@@ -1081,7 +1214,7 @@ def _inv_transform(self, a):
1081 1214
         a[masked] = exp
1082 1215
         a[~masked] /= self._linscale_adj
1083 1216
         return a
1084  
-        
  1217
+
1085 1218
     def _transform_vmin_vmax(self):
1086 1219
         """
1087 1220
         Calculates vmin and vmax in the transformed system.
@@ -1090,7 +1223,6 @@ def _transform_vmin_vmax(self):
1090 1223
         arr = np.array([vmax, vmin])
1091 1224
         self._upper, self._lower  = self._transform(arr)
1092 1225
 
1093  
-
1094 1226
     def inverse(self, value):
1095 1227
         if not self.scaled():
1096 1228
             raise ValueError("Not invertible until scaled")
44  lib/matplotlib/tests/test_colors.py
@@ -8,6 +8,7 @@
8 8
 import matplotlib.colors as mcolors
9 9
 import matplotlib.cm as cm
10 10
 
  11
+
11 12
 def test_colormap_endian():
12 13
     """
13 14
     Github issue #1005: a bug in putmask caused erroneous
@@ -23,6 +24,7 @@ def test_colormap_endian():
23 24
         #print(anative.dtype.isnative, aforeign.dtype.isnative)
24 25
         assert_array_equal(cmap(anative), cmap(aforeign))
25 26
 
  27
+
26 28
 def test_BoundaryNorm():
27 29
     """
28 30
     Github issue #1258: interpolation was failing with numpy
@@ -36,16 +38,17 @@ def test_BoundaryNorm():
36 38
     ncolors = len(boundaries)
37 39
     bn = mcolors.BoundaryNorm(boundaries, ncolors)
38 40
     assert_array_equal(bn(vals), expected)
39  
-    
  41
+
  42
+
40 43
 def test_LogNorm():
41 44
     """
42  
-    LogNorm igornoed clip, now it has the same
43  
-    behavior as Normalize, e.g. values > vmax are bigger than 1
44  
-    without clip, with clip they are 1.
  45
+    LogNorm ignored clip, now it has the same behavior as Normalize, e.g.
  46
+    values > vmax are bigger than 1 without clip, with clip they are 1.
45 47
     """
46 48
     ln = mcolors.LogNorm(clip=True, vmax=5)
47 49
     assert_array_equal(ln([1, 6]), [0, 1.0])
48 50
 
  51
+
49 52
 def test_Normalize():
50 53
     norm = mcolors.Normalize()
51 54
     vals = np.arange(-10, 10, 1, dtype=np.float)
@@ -61,19 +64,49 @@ def test_SymLogNorm():
61 64
     norm = mcolors.SymLogNorm(3, vmax=5, linscale=1.2)
62 65
     vals = np.array([-30, -1, 2, 6], dtype=np.float)
63 66
     normed_vals = norm(vals)
64  
-    expected = [ 0., 0.53980074, 0.826991, 1.02758204]
  67
+    expected = [0., 0.53980074, 0.826991, 1.02758204]
65 68
     assert_array_almost_equal(normed_vals, expected)
66 69
     _inverse_tester(norm, vals)
67 70
     _scalar_tester(norm, vals)
68 71
     _mask_tester(norm, vals)
69 72
 
70 73
 
  74
+def test_SqrtNorm():
  75
+    """
  76
+    Test SqrtNorm behavior
  77
+    """
  78
+    norm = mcolors.SqrtNorm(vmin=-3, vmax=5)
  79
+    vals = np.array([-30, -1, 2, 6], dtype=np.float)
  80
+    normed_vals = norm(vals)
  81
+    expected = [np.nan, 0.5, 0.79056942,  1.06066017]
  82
+    assert_array_almost_equal(normed_vals, expected)
  83
+    _inverse_tester(norm, vals) # note that this accepts NaN == -30 because it is masked
  84
+    _scalar_tester(norm, vals)
  85
+    _mask_tester(norm, vals)
  86
+
  87
+
  88
+def test_AsinhNorm():
  89
+    """
  90
+    Test AsinhNorm behavior
  91
+    """
  92
+    norm = mcolors.AsinhNorm(vmin=-3, vmax=5, vmid=0)
  93
+    vals = np.array([-30, -1, 2, 6], dtype=np.float)
  94
+    normed_vals = norm(vals)
  95
+    expected = [-1.694637815353593, 0.36613618958718575, 0.751895902558991,
  96
+            1.0650312050423278]
  97
+    assert_array_almost_equal(normed_vals, expected)
  98
+    _inverse_tester(norm, vals) # note that this accepts NaN == -30 because it is masked
  99
+    _scalar_tester(norm, vals)
  100
+    _mask_tester(norm, vals)
  101
+
  102
+
71 103
 def _inverse_tester(norm_instance, vals):
72 104
     """
73 105
     Checks if the inverse of the given normalization is working.
74 106
     """
75 107
     assert_array_almost_equal(norm_instance.inverse(norm_instance(vals)), vals)
76 108
 
  109
+
77 110
 def _scalar_tester(norm_instance, vals):
78 111
     """
79 112
     Checks if scalars and arrays are handled the same way.
@@ -82,6 +115,7 @@ def _scalar_tester(norm_instance, vals):
82 115
     scalar_result = [norm_instance(float(v)) for v in vals]
83 116
     assert_array_almost_equal(scalar_result, norm_instance(vals))
84 117
 
  118
+
85 119
 def _mask_tester(norm_instance, vals):
86 120
     """
87 121
     Checks mask handling
Commit_comment_tip

Tip: You can add notes to lines in a file. Hover to the left of a line to make a note

Something went wrong with that request. Please try again.