~ubuntu-branches/ubuntu/raring/python-scipy/raring-proposed

« back to all changes in this revision

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

  • Committer: Bazaar Package Importer
  • Author(s): Matthias Klose
  • Date: 2007-01-07 14:12:12 UTC
  • mfrom: (1.1.1 upstream)
  • Revision ID: james.westby@ubuntu.com-20070107141212-mm0ebkh5b37hcpzn
* Remove build dependency on python-numpy-dev.
* python-scipy: Depend on python-numpy instead of python-numpy-dev.
* Package builds on other archs than i386. Closes: #402783.

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) 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)