~ubuntu-branches/ubuntu/saucy/python-scipy/saucy

« back to all changes in this revision

Viewing changes to Lib/sandbox/image/color.py

  • Committer: Bazaar Package Importer
  • Author(s): Ondrej Certik
  • Date: 2008-06-16 22:58:01 UTC
  • mfrom: (2.1.24 intrepid)
  • Revision ID: james.westby@ubuntu.com-20080616225801-irdhrpcwiocfbcmt
Tags: 0.6.0-12
* The description updated to match the current SciPy (Closes: #489149).
* Standards-Version bumped to 3.8.0 (no action needed)
* Build-Depends: netcdf-dev changed to libnetcdf-dev

Show diffs side-by-side

added added

removed removed

Lines of Context:
1
 
 
2
 
import numpy as sb
3
 
import scipy
4
 
import os
5
 
 
6
 
# Various utilities and data for color processing
7
 
#  See http://www.cvrl.org/ for text file of the data
8
 
 
9
 
# rgb  the linear sRGB color space using D65 as a white-point
10
 
#      (IEC 61966-2-1).  Represents standard monitor (w/o gamma correction).
11
 
# rgbp is the non-linear color-space (w/ gamma correction)
12
 
# rgbntsc is the NTSC receiver primary color coordinate system
13
 
# rgbcie is the CIE, monochromatic RGB primary system
14
 
# rgbsb is the Stiles & Burch (1955) 2-deg color coordinate system
15
 
 
16
 
# Primaries for the coordinate systems
17
 
cie_primaries = [700, 546.1, 435.8]
18
 
 
19
 
sb_primaries = [1./155 * 1e5, 1./190 * 1e5, 1./225 * 1e5]
20
 
 
21
 
 
22
 
# Matrices from Jain
23
 
 
24
 
xyz_from_rgbcie = [[0.490, 0.310, 0.200],
25
 
                   [0.177, 0.813, 0.011],
26
 
                   [0.000, 0.010, 0.990]]
27
 
 
28
 
rgbcie_from_xyz = scipy.linalg.inv(xyz_from_rgbcie)
29
 
 
30
 
rgbntsc_from_xyz = [[1.910, -0.533, -0.288],
31
 
                    [-0.985, 2.000, -0.028],
32
 
                    [0.058, -0.118, 0.896]]
33
 
 
34
 
yiq_from_rgbntsc = [[0.299, 0.587, 0.114],
35
 
                    [0.596, -0.274, -0.322],
36
 
                    [0.211, -0.523, 0.312]]
37
 
 
38
 
uvw_from_xyz = [[2.0/3.0, 0, 0],
39
 
                [0,1,0],
40
 
                [-0.5,1.5,0.5]]
41
 
 
42
 
# From sRGB specification
43
 
 
44
 
xyz_from_rgb =  [[0.412453, 0.357580, 0.180423],
45
 
                 [0.212671, 0.715160, 0.072169],
46
 
                 [0.019334, 0.119193, 0.950227]]
47
 
 
48
 
rgb_from_xyz = scipy.linalg.inv(xyz_from_rgb)
49
 
 
50
 
# From http://www.mir.com/DMG/ycbcr.html
51
 
 
52
 
ycbcr_from_rgbp = [[0.299, 0.587, 0.114],
53
 
                   [-0.168736, -0.331264, 0.5],
54
 
                   [0.5, -0.418688, -0.081312]]
55
 
 
56
 
rgbp_from_ycbcr = scipy.linalg.inv(ycbcr_from_rgbp)
57
 
 
58
 
 
59
 
# LMS color space spectral matching curves provide the
60
 
#  spectral response curves of three types of cones.
61
 
#
62
 
#
63
 
# Vos, Estevez, and Walraven (1990)
64
 
# with alteration in S-cone sensitivity from
65
 
#  Stockman and Sharpe (2000)
66
 
# scaled so that sum(LMS,axis=0) has a peak of 1
67
 
#  just like LMS_from_XYZ
68
 
 
69
 
lms_from_rgbsb = [[0.14266235473644004, 0.49009667755566039,
70
 
                   0.028959576047175539],
71
 
                  [0.013676614570405768, 0.35465861798651171,
72
 
                   0.029062883056895625],
73
 
                  [0.0, 0.00029864360424843419, 0.01837806004659253]]
74
 
 
75
 
# use Judd, Vos CIE color matching curves XYZJV
76
 
#  with Stockman and Sharpe(2000) S-cone alteration.
77
 
# scaled so that sum(LMS,axis=0) has a peak of 1
78
 
#  based on CIE standard observer
79
 
 
80
 
lms_from_xyz = [[0.15513920309034629, 0.54298741130344153,
81
 
                 -0.037010041369525896],
82
 
                [-0.15513920309034629, 0.45684891207177714,
83
 
                 0.029689739651154123],
84
 
                [0.0, 6.3686624249879016e-05, 0.0073203016383768691]]
85
 
 
86
 
# Read spectral matching curves from file
87
 
# XYZJV and RGBsb55 are most modern curves to use
88
 
# LMScvrl are the cone response curves from www.cvrl.org
89
 
#   (normalized to peak at 1 for all three cones)
90
 
 
91
 
varnames = ['xyz31','xyz64','xyzjv','rgbsb55','lmscvrl']
92
 
k=-1
93
 
thisdict = globals()
94
 
for name in ['ciexyz31_1.txt','ciexyz64_1.txt','ciexyzjv.txt',
95
 
             'sbrgb2.txt','linss2_10e_1.txt']:
96
 
    k = k + 1
97
 
    name = os.path.join(os.path.split(__file__)[0],'colordata',name)
98
 
    afile = open(name)
99
 
    lines = afile.readlines()
100
 
    afile.close()
101
 
    wlen = []
102
 
    xl = []
103
 
    yl = []
104
 
    zl = []
105
 
    for line in lines:
106
 
        this = line.split(',')
107
 
        if this[0].strip()[0] not in '0123456789':
108
 
            break
109
 
        wlen.append(int(this[0].strip()))
110
 
        xl.append(float(this[1].strip()))
111
 
        yl.append(float(this[2].strip()))
112
 
        try:
113
 
            zl.append(float(this[3].strip()))
114
 
        except ValueError, inst:
115
 
            msg = inst.args[0]
116
 
            if msg.startswith("empty string"):
117
 
                zl.append(0.0)
118
 
            else:
119
 
                raise inst
120
 
 
121
 
    thisdict[varnames[k]] = (wlen,xl,yl,zl)
122
 
 
123
 
del thisdict, wlen, xl, yl, zl, afile, lines, this, line, k, msg, name
124
 
del varnames, inst
125
 
 
126
 
# XYZ white-point coordinates
127
 
#  from http://www.aim-dtp.net/aim/technology/cie_xyz/cie_xyz.htm
128
 
 
129
 
whitepoints = {'CIE A': ['Normal incandescent', 0.4476, 0.4074],
130
 
               'CIE B': ['Direct sunlight', 0.3457, 0.3585],
131
 
               'CIE C': ['Average sunlight', 0.3101, 0.3162],
132
 
               'CIE E': ['Normalized reference', 1.0/3, 1.0/3],
133
 
               'D50' : ['Bright tungsten', 0.3457, 0.3585],
134
 
               'D55' : ['Cloudy daylight', 0.3324, 0.3474],
135
 
               'D65' : ['Daylight', 0.312713, 0.329016],
136
 
               'D75' : ['?', 0.299, 0.3149],
137
 
               'D93' : ['low-quality old CRT', 0.2848, 0.2932]
138
 
               }
139
 
# convert to X,Y,Z white-point
140
 
 
141
 
def triwhite(chrwhite):
142
 
    x,y = chrwhite
143
 
    X = x / y
144
 
    Y = 1.0
145
 
    Z = (1-x-y)/y
146
 
    return X,Y,Z
147
 
 
148
 
for key in whitepoints.keys():
149
 
    whitepoints[key].append(triwhite(whitepoints[key][1:]))
150
 
del key
151
 
 
152
 
 
153
 
def tri2chr(tri,axis=None):
154
 
    """Convert tristimulus values to chromoticity values"""
155
 
    tri = sb.asarray(tri)
156
 
    n = len(tri.shape)
157
 
    if axis is None:
158
 
        axis = coloraxis(tri.shape)
159
 
    slices = []
160
 
    for k in range(n):
161
 
        slices.append(slice(None))
162
 
    slices[axis] = sb.newaxis
163
 
    norm = sb.sum(tri,axis=axis)[slices]
164
 
    slices[axis] = slice(None,2)
165
 
    out = tri[slices]/norm
166
 
    return out
167
 
 
168
 
 
169
 
# find the lowest dimension of size 3
170
 
def coloraxis(shape):
171
 
    for k, val in enumerate(shape):
172
 
        if val == 3:
173
 
            return k
174
 
    raise ValueError, "No Color axis found."
175
 
 
176
 
def convert(matrix,TTT,axis=None):
177
 
    TTT = sb.asarray(TTT)
178
 
    if axis is None:
179
 
        axis = coloraxis(TTT.shape)
180
 
    if (axis != 0):
181
 
        TTT = sb.swapaxes(TTT,0,axis)
182
 
    oldshape = TTT.shape
183
 
    TTT = sb.reshape(TTT,(3,-1))
184
 
    OUT = sb.dot(matrix, TTT)
185
 
    OUT.shape = oldshape
186
 
    if (axis != 0):
187
 
        OUT = sb.swapaxes(OUT,axis,0)
188
 
    return OUT
189
 
 
190
 
def xyz2rgbcie(xyz,axis=None):
191
 
    return convert(rgbcie_from_xyz, xyz, axis)
192
 
 
193
 
def xyz2rgb(xyz,axis=None):
194
 
    return convert(rgb_from_xyz, xyz, axis)
195
 
 
196
 
def rgb2xyz(rgb, axis=None):
197
 
    return convert(xyz_from_rgb, rgb, axis)
198
 
 
199
 
def makeslices(n):
200
 
    slices = []
201
 
    for k in range(n):
202
 
        slices.append(slice(None))
203
 
    return slices
204
 
 
205
 
def separate_colors(xyz,axis=None):
206
 
    if axis is None:
207
 
        axis = coloraxis(xyz.shape)
208
 
    n = len(xyz.shape)
209
 
    slices = makeslices(n)
210
 
    slices[axis] = 0
211
 
    x = xyz[slices]
212
 
    slices[axis] = 1
213
 
    y = xyz[slices]
214
 
    slices[axis] = 2
215
 
    z = xyz[slices]
216
 
    return x, y, z, axis
217
 
 
218
 
def join_colors(c1,c2,c3,axis):
219
 
    c1,c2,c3 = sb.asarray(c1),sb.asarray(c2),sb.asarray(c3)
220
 
    newshape = c1.shape[:axis] + (1,) + c1.shape[axis:]
221
 
    c1.shape = newshape
222
 
    c2.shape = newshape
223
 
    c3.shape = newshape
224
 
    return sb.concatenate((c1,c2,c3),axis=axis)
225
 
 
226
 
def xyz2lab(xyz, axis=None, wp=whitepoints['D65'][-1], doclip=1):
227
 
    x,y,z,axis = separate_colors(xyz,axis)
228
 
    xn,yn,zn = x/wp[0], y/wp[1], z/wp[2]
229
 
    def f(t):
230
 
        eps = 216/24389.
231
 
        kap = 24389/27.
232
 
        return sb.where(t > eps,
233
 
                        sb.power(t, 1.0/3),
234
 
                        (kap*t + 16.0)/116)
235
 
    fx,fy,fz = f(xn), f(yn), f(zn)
236
 
    L = 116*fy - 16
237
 
    a = 500*(fx - fy)
238
 
    b = 200*(fy - fz)
239
 
    if doclip:
240
 
        L = sb.clip(L, 0.0, 100.0)
241
 
        a = sb.clip(a, -500.0, 500.0)
242
 
        b = sb.clip(b, -200.0, 200.0)
243
 
    return join_colors(L,a,b,axis)
244
 
 
245
 
def lab2xyz(lab, axis=None, wp=whitepoints['D65'][-1]):
246
 
    lab = sb.asarray(lab)
247
 
    L,a,b,axis = separate_colors(lab,axis)
248
 
    fy = (L+16)/116.0
249
 
    fz = fy - b / 200.
250
 
    fx = a/500.0 + fy
251
 
    def finv(y):
252
 
        eps3 = (216/24389.)**3
253
 
        kap = 24389/27.
254
 
        return sb.where(y > eps3,
255
 
                        sb.power(y,3),
256
 
                        (116*y-16)/kap)
257
 
    xr, yr, zr = finv(fx), finv(fy), finv(fz)
258
 
    return join_colors(xr*wp[0],yr*wp[1],zr*wp[2],axis)
259
 
 
260
 
def rgb2lab(rgb):
261
 
    return xyz2lab(rgb2xyz(rgb))
262
 
 
263
 
def lab2rgb(lab):
264
 
    return xyz2rgb(lab2xyz(lab))
265
 
 
266
 
#  RGB values that will be displayed on a screen are always
267
 
#  R'G'B' values.  To get the XYZ value of the color that will be
268
 
#  displayed you need a calibrated monitor with a profile
269
 
#  -- someday we should support reading and writing such profiles and
270
 
#     doing color conversion with them.
271
 
#  But, for quick-and-dirty calculation you can often assume the sR'G'B'
272
 
#   coordinate system for your computer, and so the rgbp2rgb will
273
 
#   put you in the linear coordinate system (assuming normalized to [0,1]
274
 
#   sR'G'B' coordiates)
275
 
#
276
 
 
277
 
# sRGB <-> sR'G'B'  equations from
278
 
#   http://www.w3.org/Graphics/Color/sRGB
279
 
#   http://www.srgb.com/basicsofsrgb.htm
280
 
 
281
 
 
282
 
# Macintosh displays are usually gamma = 1.8
283
 
 
284
 
# These transformations are done with normalized [0,1.0] coordinates
285
 
 
286
 
# when gamma is None:
287
 
#    rgb2rgbp gives the nonlinear (gamma corrected) sR'G'B' from
288
 
#       linear sRGB values
289
 
#    approximately the same as rgb**(1.0/2.2)
290
 
# otherwise do a simple gamma calculation
291
 
#    rgbp = rgb**(1.0/gamma)
292
 
 
293
 
def rgb2rgbp(rgb,gamma=None):
294
 
    rgb = sb.asarray(rgb)
295
 
    if gamma is None:
296
 
        eps = 0.0031308
297
 
        return sb.where(rgb < eps, 12.92*rgb,
298
 
                     1.055*rgb**(1.0/2.4) - 0.055)
299
 
    else:
300
 
        return rgb**(1.0/gamma)
301
 
 
302
 
# when gamma is None:
303
 
#     rgbp2rgb gives linear sRGB values from nonlinear sR'G'B' values
304
 
#     approximately the same as rgbp**2.2
305
 
# otherwise do a simple gamma coorection
306
 
#     rgb = rgbp**gamma
307
 
#
308
 
def rgbp2rgb(rgbp,gamma=None):
309
 
    rgbp = sb.asarray(rgbp)
310
 
    if gamma is None:
311
 
        eps = 0.04045
312
 
        return sb.where(rgbp <= eps, rgbp / 12.92,
313
 
                     sb.power((rgbp + 0.055)/1.055,2.4))
314
 
    else:
315
 
        return rgbp**gamma
316
 
 
317
 
# The Y'CbCr coordinate system is useful because
318
 
#
319
 
# Y'CbCr information from here
320
 
#  http://www.mir.com/DMG/ycbcr.html
321
 
# This transforms from rgbp coordinates to normalized
322
 
#  y' cb cr coordinates  y' in [0,1], cb and cr in [-0.5,0.5]
323
 
#
324
 
# To convert to 8-bit use (according to the web-page cited)
325
 
#  Y' = y'*219 + 16   => [16,235]
326
 
#  Cb = cb*224 + 128  => [16,240]
327
 
#  Cr = cr*224 + 128  => [16,240]
328
 
 
329
 
def rgbp2ycbcr(rgbp,axis=None):
330
 
    return convert(ycbcr_from_rgbp, rgbp, axis)
331
 
 
332
 
def ycbcr2rgbp(ycbcr,axis=None):
333
 
    return convert(rgbp_from_ycbcr, ycbcr, axis)
334
 
 
335
 
def rgb2ycbcr(rgb,gamma=None,axis=None):
336
 
    return rgbp2ycbcr(rgb2rgbp(rgb,gamma),axis)
337
 
 
338
 
def ycbcr2rgb(ycbcr,gamma=None,axis=None):
339
 
    return rgbp2rgb(ycbcr2rgbp(ycbcr,axis),gamma)
340
 
 
341
 
def ycbcr_8bit(ycbcr,axis=None):
342
 
    y,cb,cr,axis = separate_colors(ycbcr,axis)
343
 
    Y = sb.asarray((y*219 + 16),sb.UInt8)
344
 
    Cb = sb.asarray((cb*224 + 128),sb.UInt8)
345
 
    Cr = sb.asarray((cr*224 + 128),sb.UInt8)
346
 
    return join_colors(Y,Cb,Cr,axis)
347
 
 
348
 
def ycbcr_norm(YCbCr,axis=None):
349
 
    Y,Cb,Cr,axis = separate_colors(YCbCr,axis)
350
 
    y = (Y-16.)/219
351
 
    cb = (Cb-128.)/224
352
 
    cr = (Cr-128.)/224
353
 
    return join_colors(y,cb,cr,axis)