Imported Upstream version 6.1
[platform/upstream/ffmpeg.git] / libavfilter / vf_ciescope.c
1 /*
2  * Copyright (c) 2000 John Walker
3  * Copyright (c) 2016 Paul B Mahol
4  *
5  * This file is part of FFmpeg.
6  *
7  * FFmpeg is free software; you can redistribute it and/or
8  * modify it under the terms of the GNU Lesser General Public
9  * License as published by the Free Software Foundation; either
10  * version 2.1 of the License, or (at your option) any later version.
11  *
12  * FFmpeg is distributed in the hope that it will be useful,
13  * but WITHOUT ANY WARRANTY; without even the implied warranty of
14  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
15  * Lesser General Public License for more details.
16  *
17  * You should have received a copy of the GNU Lesser General Public
18  * License along with FFmpeg; if not, write to the Free Software
19  * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
20  */
21
22 #include "libavutil/avassert.h"
23 #include "libavutil/intreadwrite.h"
24 #include "libavutil/opt.h"
25 #include "libavutil/parseutils.h"
26 #include "libavutil/pixdesc.h"
27 #include "avfilter.h"
28 #include "formats.h"
29 #include "internal.h"
30 #include "video.h"
31
32 enum CieSystem {
33     XYY,
34     UCS,
35     LUV,
36     NB_CIE
37 };
38
39 enum ColorsSystems {
40     NTSCsystem,
41     EBUsystem,
42     SMPTEsystem,
43     SMPTE240Msystem,
44     APPLEsystem,
45     wRGBsystem,
46     CIE1931system,
47     Rec709system,
48     Rec2020system,
49     DCIP3,
50     NB_CS
51 };
52
53 typedef struct CiescopeContext {
54     const AVClass *class;
55     int color_system;
56     unsigned gamuts;
57     int size;
58     int show_white;
59     int correct_gamma;
60     int cie;
61     float intensity;
62     float contrast;
63     int background;
64     int fill;
65
66     float log2lin[65536];
67     float igamma;
68     float i[3][3];
69     float m[3][3];
70     AVFrame *f;
71     void (*filter)(AVFilterContext *ctx, const uint8_t *ptr,
72                    ptrdiff_t linesize,
73                    float *cx, float *cy, int x, int y);
74 } CiescopeContext;
75
76 #define OFFSET(x) offsetof(CiescopeContext, x)
77 #define FLAGS AV_OPT_FLAG_FILTERING_PARAM|AV_OPT_FLAG_VIDEO_PARAM
78
79 static const AVOption ciescope_options[] = {
80     { "system",     "set color system", OFFSET(color_system), AV_OPT_TYPE_INT, {.i64=Rec709system}, 0, NB_CS-1, FLAGS, "system" },
81     {   "ntsc",       "NTSC 1953 Y'I'O' (ITU-R BT.470 System M)", 0, AV_OPT_TYPE_CONST, {.i64=NTSCsystem},     0, 0, FLAGS, "system" },
82     {   "470m",       "NTSC 1953 Y'I'O' (ITU-R BT.470 System M)", 0, AV_OPT_TYPE_CONST, {.i64=NTSCsystem},     0, 0, FLAGS, "system" },
83     {   "ebu",        "EBU Y'U'V' (PAL/SECAM) (ITU-R BT.470 System B, G)", 0, AV_OPT_TYPE_CONST, {.i64=EBUsystem},      0, 0, FLAGS, "system" },
84     {   "470bg",      "EBU Y'U'V' (PAL/SECAM) (ITU-R BT.470 System B, G)", 0, AV_OPT_TYPE_CONST, {.i64=EBUsystem},      0, 0, FLAGS, "system" },
85     {   "smpte",      "SMPTE-C RGB",            0, AV_OPT_TYPE_CONST, {.i64=SMPTEsystem},    0, 0, FLAGS, "system" },
86     {   "240m",       "SMPTE-240M Y'PbPr",      0, AV_OPT_TYPE_CONST, {.i64=SMPTE240Msystem},0, 0, FLAGS, "system" },
87     {   "apple",      "Apple RGB",              0, AV_OPT_TYPE_CONST, {.i64=APPLEsystem},    0, 0, FLAGS, "system" },
88     {   "widergb",    "Adobe Wide Gamut RGB",   0, AV_OPT_TYPE_CONST, {.i64=wRGBsystem},     0, 0, FLAGS, "system" },
89     {   "cie1931",    "CIE 1931 RGB",           0, AV_OPT_TYPE_CONST, {.i64=CIE1931system},  0, 0, FLAGS, "system" },
90     {   "hdtv",       "ITU.BT-709 Y'CbCr",      0, AV_OPT_TYPE_CONST, {.i64=Rec709system},   0, 0, FLAGS, "system" },
91     {   "rec709",     "ITU.BT-709 Y'CbCr",      0, AV_OPT_TYPE_CONST, {.i64=Rec709system},   0, 0, FLAGS, "system" },
92     {   "uhdtv",      "ITU-R.BT-2020",          0, AV_OPT_TYPE_CONST, {.i64=Rec2020system},  0, 0, FLAGS, "system" },
93     {   "rec2020",    "ITU-R.BT-2020",          0, AV_OPT_TYPE_CONST, {.i64=Rec2020system},  0, 0, FLAGS, "system" },
94     {   "dcip3",      "DCI-P3",                 0, AV_OPT_TYPE_CONST, {.i64=DCIP3},          0, 0, FLAGS, "system" },
95     { "cie",        "set cie system", OFFSET(cie), AV_OPT_TYPE_INT,   {.i64=XYY}, 0, NB_CIE-1, FLAGS, "cie" },
96     {   "xyy",      "CIE 1931 xyY", 0, AV_OPT_TYPE_CONST, {.i64=XYY}, 0, 0, FLAGS, "cie" },
97     {   "ucs",      "CIE 1960 UCS", 0, AV_OPT_TYPE_CONST, {.i64=UCS}, 0, 0, FLAGS, "cie" },
98     {   "luv",      "CIE 1976 Luv", 0, AV_OPT_TYPE_CONST, {.i64=LUV}, 0, 0, FLAGS, "cie" },
99     { "gamuts",     "set what gamuts to draw", OFFSET(gamuts), AV_OPT_TYPE_FLAGS, {.i64=0}, 0, 0xFFF, FLAGS, "gamuts" },
100     {   "ntsc",     NULL, 0, AV_OPT_TYPE_CONST, {.i64=1<<NTSCsystem},      0, 0, FLAGS, "gamuts" },
101     {   "470m",     NULL, 0, AV_OPT_TYPE_CONST, {.i64=1<<NTSCsystem},      0, 0, FLAGS, "gamuts" },
102     {   "ebu",      NULL, 0, AV_OPT_TYPE_CONST, {.i64=1<<EBUsystem},       0, 0, FLAGS, "gamuts" },
103     {   "470bg",    NULL, 0, AV_OPT_TYPE_CONST, {.i64=1<<EBUsystem},       0, 0, FLAGS, "gamuts" },
104     {   "smpte",    NULL, 0, AV_OPT_TYPE_CONST, {.i64=1<<SMPTEsystem},     0, 0, FLAGS, "gamuts" },
105     {   "240m",     NULL, 0, AV_OPT_TYPE_CONST, {.i64=1<<SMPTE240Msystem}, 0, 0, FLAGS, "gamuts" },
106     {   "apple",    NULL, 0, AV_OPT_TYPE_CONST, {.i64=1<<APPLEsystem},     0, 0, FLAGS, "gamuts" },
107     {   "widergb",  NULL, 0, AV_OPT_TYPE_CONST, {.i64=1<<wRGBsystem},      0, 0, FLAGS, "gamuts" },
108     {   "cie1931",  NULL, 0, AV_OPT_TYPE_CONST, {.i64=1<<CIE1931system},   0, 0, FLAGS, "gamuts" },
109     {   "hdtv",     NULL, 0, AV_OPT_TYPE_CONST, {.i64=1<<Rec709system},    0, 0, FLAGS, "gamuts" },
110     {   "rec709",   NULL, 0, AV_OPT_TYPE_CONST, {.i64=1<<Rec709system},    0, 0, FLAGS, "gamuts" },
111     {   "uhdtv",    NULL, 0, AV_OPT_TYPE_CONST, {.i64=1<<Rec2020system},   0, 0, FLAGS, "gamuts" },
112     {   "rec2020",  NULL, 0, AV_OPT_TYPE_CONST, {.i64=1<<Rec2020system},   0, 0, FLAGS, "gamuts" },
113     {   "dcip3",    NULL, 0, AV_OPT_TYPE_CONST, {.i64=1<<DCIP3},           0, 0, FLAGS, "gamuts" },
114     { "size",       "set ciescope size", OFFSET(size), AV_OPT_TYPE_INT, {.i64=512}, 256, 8192, FLAGS },
115     { "s",          "set ciescope size", OFFSET(size), AV_OPT_TYPE_INT, {.i64=512}, 256, 8192, FLAGS },
116     { "intensity",  "set ciescope intensity", OFFSET(intensity), AV_OPT_TYPE_FLOAT, {.dbl=0.001}, 0, 1, FLAGS },
117     { "i",          "set ciescope intensity", OFFSET(intensity), AV_OPT_TYPE_FLOAT, {.dbl=0.001}, 0, 1, FLAGS },
118     { "contrast",   NULL, OFFSET(contrast), AV_OPT_TYPE_FLOAT, {.dbl=0.75},  0, 1, FLAGS },
119     { "corrgamma",  NULL, OFFSET(correct_gamma), AV_OPT_TYPE_BOOL, {.i64=1}, 0, 1, FLAGS },
120     { "showwhite",  NULL, OFFSET(show_white), AV_OPT_TYPE_BOOL, {.i64=0}, 0, 1, FLAGS },
121     { "gamma",      NULL, OFFSET(igamma), AV_OPT_TYPE_DOUBLE, {.dbl=2.6}, 0.1, 6, FLAGS },
122     { "fill",       "fill with CIE colors", OFFSET(fill), AV_OPT_TYPE_BOOL, {.i64=1}, 0, 1, FLAGS },
123     { NULL }
124 };
125
126 AVFILTER_DEFINE_CLASS(ciescope);
127
128 static const enum AVPixelFormat in_pix_fmts[] = {
129     AV_PIX_FMT_RGB24,
130     AV_PIX_FMT_RGBA,
131     AV_PIX_FMT_RGB48,
132     AV_PIX_FMT_RGBA64,
133     AV_PIX_FMT_XYZ12,
134     AV_PIX_FMT_NONE
135 };
136
137 static const enum AVPixelFormat out_pix_fmts[] = {
138     AV_PIX_FMT_RGBA64,
139     AV_PIX_FMT_NONE
140 };
141
142 static int query_formats(AVFilterContext *ctx)
143 {
144     int ret;
145
146     if ((ret = ff_formats_ref(ff_make_format_list(in_pix_fmts), &ctx->inputs[0]->outcfg.formats)) < 0)
147         return ret;
148
149     if ((ret = ff_formats_ref(ff_make_format_list(out_pix_fmts), &ctx->outputs[0]->incfg.formats)) < 0)
150         return ret;
151
152     return 0;
153 }
154
155 static int config_output(AVFilterLink *outlink)
156 {
157     CiescopeContext *s = outlink->src->priv;
158
159     outlink->h = outlink->w = s->size;
160     outlink->sample_aspect_ratio = (AVRational){1,1};
161
162     return 0;
163 }
164
165 /* A  color  system is defined by the CIE x and y  coordinates of its
166    three primary illuminants and the x and y coordinates of the  white
167    point. */
168
169 struct ColorSystem {
170     float xRed, yRed,        /* Red primary illuminant */
171           xGreen, yGreen,    /* Green primary illuminant */
172           xBlue, yBlue,      /* Blue primary illuminant */
173           xWhite, yWhite,    /* White point */
174           gamma;             /* gamma of nonlinear correction */
175 };
176
177 static float const spectral_chromaticity[][3] = {
178     { 0.175560, 0.005294, 0.819146 },
179     { 0.175483, 0.005286, 0.819231 },
180     { 0.175400, 0.005279, 0.819321 },
181     { 0.175317, 0.005271, 0.819412 },
182     { 0.175237, 0.005263, 0.819500 },
183     { 0.175161, 0.005256, 0.819582 },
184     { 0.175088, 0.005247, 0.819665 },
185     { 0.175015, 0.005236, 0.819749 },
186     { 0.174945, 0.005226, 0.819829 },
187     { 0.174880, 0.005221, 0.819899 },
188     { 0.174821, 0.005221, 0.819959 },
189     { 0.174770, 0.005229, 0.820001 },
190     { 0.174722, 0.005238, 0.820040 },
191     { 0.174665, 0.005236, 0.820098 },
192     { 0.174595, 0.005218, 0.820187 },
193     { 0.174510, 0.005182, 0.820309 },
194     { 0.174409, 0.005127, 0.820464 },
195     { 0.174308, 0.005068, 0.820624 },
196     { 0.174222, 0.005017, 0.820761 },
197     { 0.174156, 0.004981, 0.820863 },
198     { 0.174112, 0.004964, 0.820924 },
199     { 0.174088, 0.004964, 0.820948 },
200     { 0.174073, 0.004973, 0.820955 },
201     { 0.174057, 0.004982, 0.820961 },
202     { 0.174036, 0.004986, 0.820978 },
203     { 0.174008, 0.004981, 0.821012 },
204     { 0.173972, 0.004964, 0.821064 },
205     { 0.173932, 0.004943, 0.821125 },
206     { 0.173889, 0.004926, 0.821185 },
207     { 0.173845, 0.004916, 0.821239 },
208     { 0.173801, 0.004915, 0.821284 },
209     { 0.173754, 0.004925, 0.821321 },
210     { 0.173705, 0.004937, 0.821358 },
211     { 0.173655, 0.004944, 0.821401 },
212     { 0.173606, 0.004940, 0.821454 },
213     { 0.173560, 0.004923, 0.821517 },
214     { 0.173514, 0.004895, 0.821590 },
215     { 0.173468, 0.004865, 0.821667 },
216     { 0.173424, 0.004836, 0.821740 },
217     { 0.173380, 0.004813, 0.821807 },
218     { 0.173337, 0.004797, 0.821866 },
219     { 0.173291, 0.004786, 0.821923 },
220     { 0.173238, 0.004779, 0.821983 },
221     { 0.173174, 0.004775, 0.822051 },
222     { 0.173101, 0.004774, 0.822125 },
223     { 0.173021, 0.004775, 0.822204 },
224     { 0.172934, 0.004781, 0.822285 },
225     { 0.172843, 0.004791, 0.822366 },
226     { 0.172751, 0.004799, 0.822450 },
227     { 0.172662, 0.004802, 0.822536 },
228     { 0.172577, 0.004799, 0.822624 },
229     { 0.172489, 0.004795, 0.822715 },
230     { 0.172396, 0.004796, 0.822808 },
231     { 0.172296, 0.004803, 0.822901 },
232     { 0.172192, 0.004815, 0.822993 },
233     { 0.172087, 0.004833, 0.823081 },
234     { 0.171982, 0.004855, 0.823163 },
235     { 0.171871, 0.004889, 0.823240 },
236     { 0.171741, 0.004939, 0.823319 },
237     { 0.171587, 0.005010, 0.823402 },
238     { 0.171407, 0.005102, 0.823490 },
239     { 0.171206, 0.005211, 0.823583 },
240     { 0.170993, 0.005334, 0.823674 },
241     { 0.170771, 0.005470, 0.823759 },
242     { 0.170541, 0.005621, 0.823838 },
243     { 0.170301, 0.005789, 0.823911 },
244     { 0.170050, 0.005974, 0.823976 },
245     { 0.169786, 0.006177, 0.824037 },
246     { 0.169505, 0.006398, 0.824097 },
247     { 0.169203, 0.006639, 0.824158 },
248     { 0.168878, 0.006900, 0.824222 },
249     { 0.168525, 0.007184, 0.824291 },
250     { 0.168146, 0.007491, 0.824363 },
251     { 0.167746, 0.007821, 0.824433 },
252     { 0.167328, 0.008175, 0.824496 },
253     { 0.166895, 0.008556, 0.824549 },
254     { 0.166446, 0.008964, 0.824589 },
255     { 0.165977, 0.009402, 0.824622 },
256     { 0.165483, 0.009865, 0.824652 },
257     { 0.164963, 0.010351, 0.824687 },
258     { 0.164412, 0.010858, 0.824731 },
259     { 0.163828, 0.011385, 0.824787 },
260     { 0.163210, 0.011937, 0.824853 },
261     { 0.162552, 0.012520, 0.824928 },
262     { 0.161851, 0.013137, 0.825011 },
263     { 0.161105, 0.013793, 0.825102 },
264     { 0.160310, 0.014491, 0.825199 },
265     { 0.159466, 0.015232, 0.825302 },
266     { 0.158573, 0.016015, 0.825412 },
267     { 0.157631, 0.016840, 0.825529 },
268     { 0.156641, 0.017705, 0.825654 },
269     { 0.155605, 0.018609, 0.825786 },
270     { 0.154525, 0.019556, 0.825920 },
271     { 0.153397, 0.020554, 0.826049 },
272     { 0.152219, 0.021612, 0.826169 },
273     { 0.150985, 0.022740, 0.826274 },
274     { 0.149691, 0.023950, 0.826359 },
275     { 0.148337, 0.025247, 0.826416 },
276     { 0.146928, 0.026635, 0.826437 },
277     { 0.145468, 0.028118, 0.826413 },
278     { 0.143960, 0.029703, 0.826337 },
279     { 0.142405, 0.031394, 0.826201 },
280     { 0.140796, 0.033213, 0.825991 },
281     { 0.139121, 0.035201, 0.825679 },
282     { 0.137364, 0.037403, 0.825233 },
283     { 0.135503, 0.039879, 0.824618 },
284     { 0.133509, 0.042692, 0.823798 },
285     { 0.131371, 0.045876, 0.822753 },
286     { 0.129086, 0.049450, 0.821464 },
287     { 0.126662, 0.053426, 0.819912 },
288     { 0.124118, 0.057803, 0.818079 },
289     { 0.121469, 0.062588, 0.815944 },
290     { 0.118701, 0.067830, 0.813468 },
291     { 0.115807, 0.073581, 0.810612 },
292     { 0.112776, 0.079896, 0.807328 },
293     { 0.109594, 0.086843, 0.803563 },
294     { 0.106261, 0.094486, 0.799253 },
295     { 0.102776, 0.102864, 0.794360 },
296     { 0.099128, 0.112007, 0.788865 },
297     { 0.095304, 0.121945, 0.782751 },
298     { 0.091294, 0.132702, 0.776004 },
299     { 0.087082, 0.144317, 0.768601 },
300     { 0.082680, 0.156866, 0.760455 },
301     { 0.078116, 0.170420, 0.751464 },
302     { 0.073437, 0.185032, 0.741531 },
303     { 0.068706, 0.200723, 0.730571 },
304     { 0.063993, 0.217468, 0.718539 },
305     { 0.059316, 0.235254, 0.705430 },
306     { 0.054667, 0.254096, 0.691238 },
307     { 0.050031, 0.274002, 0.675967 },
308     { 0.045391, 0.294976, 0.659633 },
309     { 0.040757, 0.316981, 0.642262 },
310     { 0.036195, 0.339900, 0.623905 },
311     { 0.031756, 0.363598, 0.604646 },
312     { 0.027494, 0.387921, 0.584584 },
313     { 0.023460, 0.412703, 0.563837 },
314     { 0.019705, 0.437756, 0.542539 },
315     { 0.016268, 0.462955, 0.520777 },
316     { 0.013183, 0.488207, 0.498610 },
317     { 0.010476, 0.513404, 0.476120 },
318     { 0.008168, 0.538423, 0.453409 },
319     { 0.006285, 0.563068, 0.430647 },
320     { 0.004875, 0.587116, 0.408008 },
321     { 0.003982, 0.610447, 0.385570 },
322     { 0.003636, 0.633011, 0.363352 },
323     { 0.003859, 0.654823, 0.341318 },
324     { 0.004646, 0.675898, 0.319456 },
325     { 0.006011, 0.696120, 0.297869 },
326     { 0.007988, 0.715342, 0.276670 },
327     { 0.010603, 0.733413, 0.255984 },
328     { 0.013870, 0.750186, 0.235943 },
329     { 0.017766, 0.765612, 0.216622 },
330     { 0.022244, 0.779630, 0.198126 },
331     { 0.027273, 0.792104, 0.180623 },
332     { 0.032820, 0.802926, 0.164254 },
333     { 0.038852, 0.812016, 0.149132 },
334     { 0.045328, 0.819391, 0.135281 },
335     { 0.052177, 0.825164, 0.122660 },
336     { 0.059326, 0.829426, 0.111249 },
337     { 0.066716, 0.832274, 0.101010 },
338     { 0.074302, 0.833803, 0.091894 },
339     { 0.082053, 0.834090, 0.083856 },
340     { 0.089942, 0.833289, 0.076769 },
341     { 0.097940, 0.831593, 0.070468 },
342     { 0.106021, 0.829178, 0.064801 },
343     { 0.114161, 0.826207, 0.059632 },
344     { 0.122347, 0.822770, 0.054882 },
345     { 0.130546, 0.818928, 0.050526 },
346     { 0.138702, 0.814774, 0.046523 },
347     { 0.146773, 0.810395, 0.042832 },
348     { 0.154722, 0.805864, 0.039414 },
349     { 0.162535, 0.801238, 0.036226 },
350     { 0.170237, 0.796519, 0.033244 },
351     { 0.177850, 0.791687, 0.030464 },
352     { 0.185391, 0.786728, 0.027881 },
353     { 0.192876, 0.781629, 0.025495 },
354     { 0.200309, 0.776399, 0.023292 },
355     { 0.207690, 0.771055, 0.021255 },
356     { 0.215030, 0.765595, 0.019375 },
357     { 0.222337, 0.760020, 0.017643 },
358     { 0.229620, 0.754329, 0.016051 },
359     { 0.236885, 0.748524, 0.014591 },
360     { 0.244133, 0.742614, 0.013253 },
361     { 0.251363, 0.736606, 0.012031 },
362     { 0.258578, 0.730507, 0.010916 },
363     { 0.265775, 0.724324, 0.009901 },
364     { 0.272958, 0.718062, 0.008980 },
365     { 0.280129, 0.711725, 0.008146 },
366     { 0.287292, 0.705316, 0.007391 },
367     { 0.294450, 0.698842, 0.006708 },
368     { 0.301604, 0.692308, 0.006088 },
369     { 0.308760, 0.685712, 0.005528 },
370     { 0.315914, 0.679063, 0.005022 },
371     { 0.323066, 0.672367, 0.004566 },
372     { 0.330216, 0.665628, 0.004156 },
373     { 0.337363, 0.658848, 0.003788 },
374     { 0.344513, 0.652028, 0.003459 },
375     { 0.351664, 0.645172, 0.003163 },
376     { 0.358814, 0.638287, 0.002899 },
377     { 0.365959, 0.631379, 0.002662 },
378     { 0.373102, 0.624451, 0.002448 },
379     { 0.380244, 0.617502, 0.002254 },
380     { 0.387379, 0.610542, 0.002079 },
381     { 0.394507, 0.603571, 0.001922 },
382     { 0.401626, 0.596592, 0.001782 },
383     { 0.408736, 0.589607, 0.001657 },
384     { 0.415836, 0.582618, 0.001546 },
385     { 0.422921, 0.575631, 0.001448 },
386     { 0.429989, 0.568649, 0.001362 },
387     { 0.437036, 0.561676, 0.001288 },
388     { 0.444062, 0.554714, 0.001224 },
389     { 0.451065, 0.547766, 0.001169 },
390     { 0.458041, 0.540837, 0.001123 },
391     { 0.464986, 0.533930, 0.001084 },
392     { 0.471899, 0.527051, 0.001051 },
393     { 0.478775, 0.520202, 0.001023 },
394     { 0.485612, 0.513389, 0.001000 },
395     { 0.492405, 0.506615, 0.000980 },
396     { 0.499151, 0.499887, 0.000962 },
397     { 0.505845, 0.493211, 0.000944 },
398     { 0.512486, 0.486591, 0.000923 },
399     { 0.519073, 0.480029, 0.000899 },
400     { 0.525600, 0.473527, 0.000872 },
401     { 0.532066, 0.467091, 0.000843 },
402     { 0.538463, 0.460725, 0.000812 },
403     { 0.544787, 0.454434, 0.000779 },
404     { 0.551031, 0.448225, 0.000744 },
405     { 0.557193, 0.442099, 0.000708 },
406     { 0.563269, 0.436058, 0.000673 },
407     { 0.569257, 0.430102, 0.000641 },
408     { 0.575151, 0.424232, 0.000616 },
409     { 0.580953, 0.418447, 0.000601 },
410     { 0.586650, 0.412758, 0.000591 },
411     { 0.592225, 0.407190, 0.000586 },
412     { 0.597658, 0.401762, 0.000580 },
413     { 0.602933, 0.396497, 0.000571 },
414     { 0.608035, 0.391409, 0.000556 },
415     { 0.612977, 0.386486, 0.000537 },
416     { 0.617779, 0.381706, 0.000516 },
417     { 0.622459, 0.377047, 0.000493 },
418     { 0.627037, 0.372491, 0.000472 },
419     { 0.631521, 0.368026, 0.000453 },
420     { 0.635900, 0.363665, 0.000435 },
421     { 0.640156, 0.359428, 0.000416 },
422     { 0.644273, 0.355331, 0.000396 },
423     { 0.648233, 0.351395, 0.000372 },
424     { 0.652028, 0.347628, 0.000344 },
425     { 0.655669, 0.344018, 0.000313 },
426     { 0.659166, 0.340553, 0.000281 },
427     { 0.662528, 0.337221, 0.000251 },
428     { 0.665764, 0.334011, 0.000226 },
429     { 0.668874, 0.330919, 0.000207 },
430     { 0.671859, 0.327947, 0.000194 },
431     { 0.674720, 0.325095, 0.000185 },
432     { 0.677459, 0.322362, 0.000179 },
433     { 0.680079, 0.319747, 0.000174 },
434     { 0.682582, 0.317249, 0.000170 },
435     { 0.684971, 0.314863, 0.000167 },
436     { 0.687250, 0.312586, 0.000164 },
437     { 0.689426, 0.310414, 0.000160 },
438     { 0.691504, 0.308342, 0.000154 },
439     { 0.693490, 0.306366, 0.000145 },
440     { 0.695389, 0.304479, 0.000133 },
441     { 0.697206, 0.302675, 0.000119 },
442     { 0.698944, 0.300950, 0.000106 },
443     { 0.700606, 0.299301, 0.000093 },
444     { 0.702193, 0.297725, 0.000083 },
445     { 0.703709, 0.296217, 0.000074 },
446     { 0.705163, 0.294770, 0.000067 },
447     { 0.706563, 0.293376, 0.000061 },
448     { 0.707918, 0.292027, 0.000055 },
449     { 0.709231, 0.290719, 0.000050 },
450     { 0.710500, 0.289453, 0.000047 },
451     { 0.711724, 0.288232, 0.000044 },
452     { 0.712901, 0.287057, 0.000041 },
453     { 0.714032, 0.285929, 0.000040 },
454     { 0.715117, 0.284845, 0.000038 },
455     { 0.716159, 0.283804, 0.000036 },
456     { 0.717159, 0.282806, 0.000035 },
457     { 0.718116, 0.281850, 0.000034 },
458     { 0.719033, 0.280935, 0.000032 },
459     { 0.719912, 0.280058, 0.000030 },
460     { 0.720753, 0.279219, 0.000028 },
461     { 0.721555, 0.278420, 0.000026 },
462     { 0.722315, 0.277662, 0.000023 },
463     { 0.723032, 0.276948, 0.000020 },
464     { 0.723702, 0.276282, 0.000016 },
465     { 0.724328, 0.275660, 0.000012 },
466     { 0.724914, 0.275078, 0.000007 },
467     { 0.725467, 0.274530, 0.000003 },
468     { 0.725992, 0.274008, 0.000000 },
469     { 0.726495, 0.273505, 0.000000 },
470     { 0.726975, 0.273025, 0.000000 },
471     { 0.727432, 0.272568, 0.000000 },
472     { 0.727864, 0.272136, 0.000000 },
473     { 0.728272, 0.271728, 0.000000 },
474     { 0.728656, 0.271344, 0.000000 },
475     { 0.729020, 0.270980, 0.000000 },
476     { 0.729361, 0.270639, 0.000000 },
477     { 0.729678, 0.270322, 0.000000 },
478     { 0.729969, 0.270031, 0.000000 },
479     { 0.730234, 0.269766, 0.000000 },
480     { 0.730474, 0.269526, 0.000000 },
481     { 0.730693, 0.269307, 0.000000 },
482     { 0.730896, 0.269104, 0.000000 },
483     { 0.731089, 0.268911, 0.000000 },
484     { 0.731280, 0.268720, 0.000000 },
485     { 0.731467, 0.268533, 0.000000 },
486     { 0.731650, 0.268350, 0.000000 },
487     { 0.731826, 0.268174, 0.000000 },
488     { 0.731993, 0.268007, 0.000000 },
489     { 0.732150, 0.267850, 0.000000 },
490     { 0.732300, 0.267700, 0.000000 },
491     { 0.732443, 0.267557, 0.000000 },
492     { 0.732581, 0.267419, 0.000000 },
493     { 0.732719, 0.267281, 0.000000 },
494     { 0.732859, 0.267141, 0.000000 },
495     { 0.733000, 0.267000, 0.000000 },
496     { 0.733142, 0.266858, 0.000000 },
497     { 0.733281, 0.266719, 0.000000 },
498     { 0.733417, 0.266583, 0.000000 },
499     { 0.733551, 0.266449, 0.000000 },
500     { 0.733683, 0.266317, 0.000000 },
501     { 0.733813, 0.266187, 0.000000 },
502     { 0.733936, 0.266064, 0.000000 },
503     { 0.734047, 0.265953, 0.000000 },
504     { 0.734143, 0.265857, 0.000000 },
505     { 0.734221, 0.265779, 0.000000 },
506     { 0.734286, 0.265714, 0.000000 },
507     { 0.734341, 0.265659, 0.000000 },
508     { 0.734390, 0.265610, 0.000000 },
509     { 0.734438, 0.265562, 0.000000 },
510     { 0.734482, 0.265518, 0.000000 },
511     { 0.734523, 0.265477, 0.000000 },
512     { 0.734560, 0.265440, 0.000000 },
513     { 0.734592, 0.265408, 0.000000 },
514     { 0.734621, 0.265379, 0.000000 },
515     { 0.734649, 0.265351, 0.000000 },
516     { 0.734673, 0.265327, 0.000000 },
517     { 0.734690, 0.265310, 0.000000 },
518     { 0.734690, 0.265310, 0.000000 },
519     { 0.734690, 0.265310, 0.000000 },
520     { 0.734690, 0.265310, 0.000000 },
521     { 0.734690, 0.265310, 0.000000 },
522     { 0.734690, 0.265310, 0.000000 },
523     { 0.734690, 0.265310, 0.000000 },
524     { 0.734690, 0.265310, 0.000000 },
525     { 0.734690, 0.265310, 0.000000 },
526     { 0.734690, 0.265310, 0.000000 },
527     { 0.734690, 0.265310, 0.000000 },
528     { 0.734690, 0.265310, 0.000000 },
529     { 0.734690, 0.265310, 0.000000 },
530     { 0.734690, 0.265310, 0.000000 },
531     { 0.734690, 0.265310, 0.000000 },
532     { 0.734690, 0.265310, 0.000000 },
533     { 0.734690, 0.265310, 0.000000 },
534     { 0.734690, 0.265310, 0.000000 },
535     { 0.734690, 0.265310, 0.000000 },
536     { 0.734690, 0.265310, 0.000000 },
537     { 0.734690, 0.265310, 0.000000 },
538     { 0.734690, 0.265310, 0.000000 },
539     { 0.734690, 0.265310, 0.000000 },
540     { 0.734690, 0.265310, 0.000000 },
541     { 0.734690, 0.265310, 0.000000 },
542     { 0.734690, 0.265310, 0.000000 },
543     { 0.734690, 0.265310, 0.000000 },
544     { 0.734690, 0.265310, 0.000000 },
545     { 0.734690, 0.265310, 0.000000 },
546     { 0.734690, 0.265310, 0.000000 },
547     { 0.734690, 0.265310, 0.000000 },
548     { 0.734690, 0.265310, 0.000000 },
549     { 0.734690, 0.265310, 0.000000 },
550     { 0.734690, 0.265310, 0.000000 },
551     { 0.734690, 0.265310, 0.000000 },
552     { 0.734690, 0.265310, 0.000000 },
553     { 0.734690, 0.265310, 0.000000 },
554     { 0.734690, 0.265310, 0.000000 },
555     { 0.734690, 0.265310, 0.000000 },
556     { 0.734690, 0.265310, 0.000000 },
557     { 0.734690, 0.265310, 0.000000 },
558     { 0.734690, 0.265310, 0.000000 },
559     { 0.734690, 0.265310, 0.000000 },
560     { 0.734690, 0.265310, 0.000000 },
561     { 0.734690, 0.265310, 0.000000 },
562     { 0.734690, 0.265310, 0.000000 },
563     { 0.734690, 0.265310, 0.000000 },
564     { 0.734690, 0.265310, 0.000000 },
565     { 0.734690, 0.265310, 0.000000 },
566     { 0.734690, 0.265310, 0.000000 },
567     { 0.734690, 0.265310, 0.000000 },
568     { 0.734690, 0.265310, 0.000000 },
569     { 0.734690, 0.265310, 0.000000 },
570     { 0.734690, 0.265310, 0.000000 },
571     { 0.734690, 0.265310, 0.000000 },
572     { 0.734690, 0.265310, 0.000000 },
573     { 0.734690, 0.265310, 0.000000 },
574     { 0.734690, 0.265310, 0.000000 },
575     { 0.734690, 0.265310, 0.000000 },
576     { 0.734690, 0.265310, 0.000000 },
577     { 0.734690, 0.265310, 0.000000 },
578     { 0.734690, 0.265310, 0.000000 },
579     { 0.734690, 0.265310, 0.000000 },
580     { 0.734690, 0.265310, 0.000000 },
581     { 0.734690, 0.265310, 0.000000 },
582     { 0.734690, 0.265310, 0.000000 },
583     { 0.734690, 0.265310, 0.000000 },
584     { 0.734690, 0.265310, 0.000000 },
585     { 0.734690, 0.265310, 0.000000 },
586     { 0.734690, 0.265310, 0.000000 },
587     { 0.734690, 0.265310, 0.000000 },
588     { 0.734690, 0.265310, 0.000000 },
589     { 0.734690, 0.265310, 0.000000 },
590     { 0.734690, 0.265310, 0.000000 },
591     { 0.734690, 0.265310, 0.000000 },
592     { 0.734690, 0.265310, 0.000000 },
593     { 0.734690, 0.265310, 0.000000 },
594     { 0.734690, 0.265310, 0.000000 },
595     { 0.734690, 0.265310, 0.000000 },
596     { 0.734690, 0.265310, 0.000000 },
597     { 0.734690, 0.265310, 0.000000 },
598     { 0.734690, 0.265310, 0.000000 },
599     { 0.734690, 0.265310, 0.000000 },
600     { 0.734690, 0.265310, 0.000000 },
601     { 0.734690, 0.265310, 0.000000 },
602     { 0.734690, 0.265310, 0.000000 },
603     { 0.734690, 0.265310, 0.000000 },
604     { 0.734690, 0.265310, 0.000000 },
605     { 0.734690, 0.265310, 0.000000 },
606     { 0.734690, 0.265310, 0.000000 },
607     { 0.734690, 0.265310, 0.000000 },
608     { 0.734690, 0.265310, 0.000000 },
609     { 0.734690, 0.265310, 0.000000 },
610     { 0.734690, 0.265310, 0.000000 },
611     { 0.734690, 0.265310, 0.000000 },
612     { 0.734690, 0.265310, 0.000000 },
613     { 0.734690, 0.265310, 0.000000 },
614     { 0.734690, 0.265310, 0.000000 },
615     { 0.734690, 0.265310, 0.000000 },
616     { 0.734690, 0.265310, 0.000000 },
617     { 0.734690, 0.265310, 0.000000 },
618     { 0.734690, 0.265310, 0.000000 },
619     { 0.734690, 0.265310, 0.000000 },
620     { 0.734690, 0.265310, 0.000000 },
621     { 0.734690, 0.265310, 0.000000 },
622     { 0.734690, 0.265310, 0.000000 },
623     { 0.734690, 0.265310, 0.000000 },
624     { 0.734690, 0.265310, 0.000000 },
625     { 0.734690, 0.265310, 0.000000 },
626     { 0.734690, 0.265310, 0.000000 },
627     { 0.734690, 0.265310, 0.000000 },
628     { 0.734690, 0.265310, 0.000000 },
629     { 0.734690, 0.265310, 0.000000 },
630     { 0.734690, 0.265310, 0.000000 },
631     { 0.734690, 0.265310, 0.000000 },
632     { 0.734690, 0.265310, 0.000000 },
633     { 0.734690, 0.265310, 0.000000 },
634     { 0.734690, 0.265310, 0.000000 },
635     { 0.734690, 0.265310, 0.000000 },
636     { 0.734690, 0.265310, 0.000000 },
637     { 0.734690, 0.265310, 0.000000 },
638     { 0.734690, 0.265310, 0.000000 },
639     { 0.734690, 0.265310, 0.000000 },
640     { 0.734690, 0.265310, 0.000000 },
641     { 0.734690, 0.265310, 0.000000 },
642     { 0.734690, 0.265310, 0.000000 },
643     { 0.734690, 0.265310, 0.000000 },
644     { 0.734690, 0.265310, 0.000000 },
645     { 0.734690, 0.265310, 0.000000 },
646     { 0.734690, 0.265310, 0.000000 },
647     { 0.734690, 0.265310, 0.000000 },
648     { 0.734690, 0.265310, 0.000000 },
649 };
650
651
652 /* Standard white point chromaticities. */
653
654 #define C     0.310063, 0.316158
655 #define E     1.0/3.0, 1.0/3.0
656 #define D50   0.34570, 0.3585
657 #define D65   0.312713, 0.329016
658
659 /* Gamma of nonlinear correction.
660    See Charles Poynton's ColorFAQ Item 45 and GammaFAQ Item 6 at
661    http://www.inforamp.net/~poynton/ColorFAQ.html
662    http://www.inforamp.net/~poynton/GammaFAQ.html
663 */
664
665 #define GAMMA_REC709    0.      /* Rec. 709 */
666
667 static const struct ColorSystem color_systems[] = {
668     [NTSCsystem] = {
669         0.67,  0.33,  0.21,  0.71,  0.14,  0.08,
670         C, GAMMA_REC709
671     },
672     [EBUsystem] = {
673         0.64,  0.33,  0.29,  0.60,  0.15,  0.06,
674         D65, GAMMA_REC709
675     },
676     [SMPTEsystem] = {
677         0.630, 0.340, 0.310, 0.595, 0.155, 0.070,
678         D65, GAMMA_REC709
679     },
680     [SMPTE240Msystem] = {
681         0.670, 0.330, 0.210, 0.710, 0.150, 0.060,
682         D65, GAMMA_REC709
683     },
684     [APPLEsystem] = {
685         0.625, 0.340, 0.280, 0.595, 0.115, 0.070,
686         D65, GAMMA_REC709
687     },
688     [wRGBsystem] = {
689         0.7347, 0.2653, 0.1152, 0.8264, 0.1566, 0.0177,
690         D50, GAMMA_REC709
691     },
692     [CIE1931system] = {
693         0.7347, 0.2653, 0.2738, 0.7174, 0.1666, 0.0089,
694         E, GAMMA_REC709
695     },
696     [Rec709system] = {
697         0.64,  0.33,  0.30,  0.60,  0.15,  0.06,
698         D65, GAMMA_REC709
699     },
700     [Rec2020system] = {
701         0.708,  0.292,  0.170,  0.797,  0.131,  0.046,
702         D65, GAMMA_REC709
703     },
704     [DCIP3] = {
705         0.680,  0.320,  0.265,  0.690,  0.150,  0.060,
706         0.314,  0.351, GAMMA_REC709
707     },
708 };
709
710 /*
711 static struct ColorSystem CustomSystem = {
712     "Custom",
713     0.64,  0.33,  0.30,  0.60,  0.15,  0.06,
714     D65, GAMMA_REC709
715 };
716 */
717
718 static void
719 uv_to_xy(float  const u,
720          float  const v,
721          float *const xc,
722          float *const yc)
723 {
724 /*
725     Given 1970 coordinates u, v, determine 1931 chromaticities x, y
726 */
727     *xc = 3.f*u / (2.f*u - 8.f*v + 4.f);
728     *yc = 2.f*v / (2.f*u - 8.f*v + 4.f);
729 }
730
731 static void
732 upvp_to_xy(float   const up,
733            float   const vp,
734            float * const xc,
735            float * const yc)
736 {
737 /*
738     Given 1976 coordinates u', v', determine 1931 chromaticities x, y
739 */
740     *xc = 9*up / (6*up - 16*vp + 12);
741     *yc = 4*vp / (6*up - 16*vp + 12);
742 }
743
744 static void
745 xy_to_upvp(float xc,
746            float yc,
747            float * const up,
748            float * const vp)
749 {
750 /*
751     Given 1931 chromaticities x, y, determine 1976 coordinates u', v'
752 */
753     const float scale = 1.f / (-2.f*xc + 12.f*yc + 3.f);
754     *up = 4.f*xc * scale;
755     *vp = 9.f*yc * scale;
756 }
757
758 static void
759 xy_to_uv(float xc,
760          float yc,
761          float * const u,
762          float * const v)
763 {
764 /*
765     Given 1931 chromaticities x, y, determine 1960 coordinates u, v
766 */
767     const float scale = 1.f / (-2.f*xc + 12.f*yc + 3.f);
768     *u = 4.f*xc * scale;
769     *v = 6.f*yc * scale;
770 }
771
772 static void
773 xyz_to_rgb(const float m[3][3],
774            float xc, float yc, float zc,
775            float * const r, float * const g, float * const b)
776 {
777     *r = m[0][0]*xc + m[0][1]*yc + m[0][2]*zc;
778     *g = m[1][0]*xc + m[1][1]*yc + m[1][2]*zc;
779     *b = m[2][0]*xc + m[2][1]*yc + m[2][2]*zc;
780 }
781
782 static void invert_matrix3x3(float in[3][3], float out[3][3])
783 {
784     float m00 = in[0][0], m01 = in[0][1], m02 = in[0][2],
785            m10 = in[1][0], m11 = in[1][1], m12 = in[1][2],
786            m20 = in[2][0], m21 = in[2][1], m22 = in[2][2];
787     int i, j;
788     float det;
789
790     out[0][0] =  (m11 * m22 - m21 * m12);
791     out[0][1] = -(m01 * m22 - m21 * m02);
792     out[0][2] =  (m01 * m12 - m11 * m02);
793     out[1][0] = -(m10 * m22 - m20 * m12);
794     out[1][1] =  (m00 * m22 - m20 * m02);
795     out[1][2] = -(m00 * m12 - m10 * m02);
796     out[2][0] =  (m10 * m21 - m20 * m11);
797     out[2][1] = -(m00 * m21 - m20 * m01);
798     out[2][2] =  (m00 * m11 - m10 * m01);
799
800     det = m00 * out[0][0] + m10 * out[0][1] + m20 * out[0][2];
801     det = 1.0 / det;
802
803     for (i = 0; i < 3; i++) {
804         for (j = 0; j < 3; j++)
805             out[i][j] *= det;
806     }
807 }
808
809 static void get_rgb2xyz_matrix(struct ColorSystem system, float m[3][3])
810 {
811     float S[3], X[4], Z[4];
812     int i;
813
814     X[0] = system.xRed   / system.yRed;
815     X[1] = system.xGreen / system.yGreen;
816     X[2] = system.xBlue  / system.yBlue;
817     X[3] = system.xWhite / system.yWhite;
818
819     Z[0] = (1 - system.xRed   - system.yRed)   / system.yRed;
820     Z[1] = (1 - system.xGreen - system.yGreen) / system.yGreen;
821     Z[2] = (1 - system.xBlue  - system.yBlue)  / system.yBlue;
822     Z[3] = (1 - system.xWhite - system.yWhite) / system.yWhite;
823
824     for (i = 0; i < 3; i++) {
825         m[0][i] = X[i];
826         m[1][i] = 1;
827         m[2][i] = Z[i];
828     }
829
830     invert_matrix3x3(m, m);
831
832     for (i = 0; i < 3; i++)
833         S[i] = m[i][0] * X[3] + m[i][1] * 1 + m[i][2] * Z[3];
834
835     for (i = 0; i < 3; i++) {
836         m[0][i] = S[i] * X[i];
837         m[1][i] = S[i] * 1;
838         m[2][i] = S[i] * Z[i];
839     }
840 }
841
842 static void
843 rgb_to_xy(float rc,
844           float gc,
845           float bc,
846           float * const x,
847           float * const y,
848           float * const z,
849           const float m[3][3])
850 {
851     float scale;
852
853     *x = m[0][0] * rc + m[0][1] * gc + m[0][2] * bc;
854     *y = m[1][0] * rc + m[1][1] * gc + m[1][2] * bc;
855     *z = m[2][0] * rc + m[2][1] * gc + m[2][2] * bc;
856
857     scale = *x + *y + *z;
858     scale = 1.f / scale;
859     *x = *x * scale;
860     *y = *y * scale;
861 }
862
863 static int
864 constrain_rgb(float * const r,
865               float * const g,
866               float * const b)
867 {
868 /*----------------------------------------------------------------------------
869     If  the  requested RGB shade contains a negative weight for one of
870     the primaries, it lies outside the color  gamut  accessible  from
871     the  given  triple  of  primaries.  Desaturate it by adding white,
872     equal quantities of R, G, and B, enough to make RGB all positive.
873 -----------------------------------------------------------------------------*/
874     float w;
875
876     /* Amount of white needed is w = - min(0, *r, *g, *b) */
877     w = (0 < *r) ? 0 : *r;
878     w = (w < *g) ? w : *g;
879     w = (w < *b) ? w : *b;
880     w = - w;
881
882     /* Add just enough white to make r, g, b all positive. */
883     if (w > 0) {
884         *r += w;  *g += w; *b += w;
885
886         return 1;                     /* Color modified to fit RGB gamut */
887     }
888
889     return 0;                         /* Color within RGB gamut */
890 }
891
892 static void
893 gamma_correct(const struct ColorSystem * const cs,
894               float *                   const c)
895 {
896 /*----------------------------------------------------------------------------
897     Transform linear RGB values to nonlinear RGB values.
898
899     Rec. 709 is ITU-R Recommendation BT. 709 (1990)
900     ``Basic Parameter Values for the HDTV Standard for the Studio and for
901     International Programme Exchange'', formerly CCIR Rec. 709.
902
903     For details see
904        http://www.inforamp.net/~poynton/ColorFAQ.html
905        http://www.inforamp.net/~poynton/GammaFAQ.html
906 -----------------------------------------------------------------------------*/
907     float gamma;
908     float cc;
909
910     gamma = cs->gamma;
911
912     if (gamma == 0.) {
913         /* Rec. 709 gamma correction. */
914         cc = 0.018;
915         if (*c < cc) {
916             *c *= (1.099 * pow(cc, 0.45) - 0.099) / cc;
917         } else {
918             *c = 1.099 * pow(*c, 0.45) - 0.099;
919         }
920     } else {
921     /* Nonlinear color = (Linear color)^(1/gamma) */
922         *c = pow(*c, 1./gamma);
923     }
924 }
925
926
927
928 static void
929 gamma_correct_rgb(const struct ColorSystem * const cs,
930                   float * const r,
931                   float * const g,
932                   float * const b)
933 {
934     gamma_correct(cs, r);
935     gamma_correct(cs, g);
936     gamma_correct(cs, b);
937 }
938
939 /* Sz(X) is the displacement in pixels of a displacement of X normalized
940    distance units.  (A normalized distance unit is 1/512 of the smaller
941    dimension of the canvas)
942 */
943 #define Sz(x) (((x) * (int)FFMIN(w, h)) / 512)
944
945 static void
946 monochrome_color_location(float waveLength, int w, int h,
947                           int cie, int *xP, int *yP)
948 {
949     const int ix = waveLength - 360;
950     const float pX = spectral_chromaticity[ix][0];
951     const float pY = spectral_chromaticity[ix][1];
952     const float pZ = spectral_chromaticity[ix][2];
953     const float px = pX / (pX + pY + pZ);
954     const float py = pY / (pX + pY + pZ);
955
956     if (cie == LUV) {
957         float up, vp;
958
959         xy_to_upvp(px, py, &up, &vp);
960         *xP = up * (w - 1);
961         *yP = (h - 1) - vp * (h - 1);
962     } else if (cie == UCS) {
963         float u, v;
964
965         xy_to_uv(px, py, &u, &v);
966         *xP = u * (w - 1);
967         *yP = (h - 1) - v * (h - 1);
968     } else if (cie == XYY) {
969         *xP = px * (w - 1);
970         *yP = (h - 1) - py * (h - 1);
971     } else {
972         av_assert0(0);
973     }
974 }
975
976 static void
977 find_tongue(uint16_t* const pixels,
978             int       const w,
979             int       const linesize,
980             int       const row,
981             int *     const presentP,
982             int *     const leftEdgeP,
983             int *     const rightEdgeP)
984 {
985     int i;
986
987     for (i = 0; i < w && pixels[row * linesize + i * 4 + 0] == 0; i++)
988         ;
989
990     if (i >= w) {
991         *presentP = 0;
992     } else {
993         int j;
994         int const leftEdge = i;
995
996         *presentP = 1;
997
998         for (j = w - 1; j >= leftEdge && pixels[row * linesize + j * 4 + 0] == 0; j--)
999             ;
1000
1001         *rightEdgeP = j;
1002         *leftEdgeP = leftEdge;
1003     }
1004 }
1005
1006 static void draw_line(uint16_t *const pixels, int linesize,
1007                       int x0, int y0, int x1, int y1,
1008                       int w, int h,
1009                       const uint16_t *const rgbcolor)
1010 {
1011     int sx = x0 < x1 ? 1 : -1, sy = y0 < y1 ? 1 : -1, x2;
1012     int dx = FFABS(x1-x0), dy = FFABS(y1-y0), err = dx * dx + dy * dy;
1013     int e2 = err == 0 ? 1 : 0xffffff / (dx + dy);
1014
1015     dx *= e2;
1016     dy *= e2;
1017     err = dx - dy;
1018
1019     for (;;) {
1020         pixels[y0 * linesize + x0 * 4 + 0] = rgbcolor[0]-(FFABS(err - dx + dy) >> 8);
1021         pixels[y0 * linesize + x0 * 4 + 1] = rgbcolor[1]-(FFABS(err - dx + dy) >> 8);
1022         pixels[y0 * linesize + x0 * 4 + 2] = rgbcolor[2]-(FFABS(err - dx + dy) >> 8);
1023         pixels[y0 * linesize + x0 * 4 + 3] = rgbcolor[3]-(FFABS(err - dx + dy) >> 8);
1024
1025         e2 = err;
1026         x2 = x0;
1027         if (2 * e2 >= -dx) {
1028             if (x0 == x1)
1029                 break;
1030             if (e2 + dy < 0xff0000) {
1031                 pixels[(y0 + sy) * linesize + x0 * 4 + 0] = rgbcolor[0]-(FFABS(e2 + dy) >> 8);
1032                 pixels[(y0 + sy) * linesize + x0 * 4 + 1] = rgbcolor[1]-(FFABS(e2 + dy) >> 8);
1033                 pixels[(y0 + sy) * linesize + x0 * 4 + 2] = rgbcolor[2]-(FFABS(e2 + dy) >> 8);
1034                 pixels[(y0 + sy) * linesize + x0 * 4 + 3] = rgbcolor[3]-(FFABS(e2 + dy) >> 8);
1035             }
1036             err -= dy;
1037             x0 += sx;
1038         }
1039
1040         if (2 * e2 <= dy) {
1041             if (y0 == y1)
1042                 break;
1043             if (dx - e2 < 0xff0000) {
1044                 pixels[y0 * linesize + (x2 + sx) * 4 + 0] = rgbcolor[0]-(FFABS(dx - e2) >> 8);
1045                 pixels[y0 * linesize + (x2 + sx) * 4 + 1] = rgbcolor[1]-(FFABS(dx - e2) >> 8);
1046                 pixels[y0 * linesize + (x2 + sx) * 4 + 2] = rgbcolor[2]-(FFABS(dx - e2) >> 8);
1047                 pixels[y0 * linesize + (x2 + sx) * 4 + 3] = rgbcolor[3]-(FFABS(dx - e2) >> 8);
1048             }
1049             err += dx;
1050             y0 += sy;
1051         }
1052     }
1053 }
1054
1055 static void draw_rline(uint16_t *const pixels, int linesize,
1056                        int x0, int y0, int x1, int y1,
1057                        int w, int h)
1058 {
1059     int dx  = FFABS(x1 - x0), sx = x0 < x1 ? 1 : -1;
1060     int dy  = FFABS(y1 - y0), sy = y0 < y1 ? 1 : -1;
1061     int err = (dx > dy ? dx : -dy) / 2, e2;
1062
1063     for (;;) {
1064         pixels[y0 * linesize + x0 * 4 + 0] = 65535 - pixels[y0 * linesize + x0 * 4 + 0];
1065         pixels[y0 * linesize + x0 * 4 + 1] = 65535 - pixels[y0 * linesize + x0 * 4 + 1];
1066         pixels[y0 * linesize + x0 * 4 + 2] = 65535 - pixels[y0 * linesize + x0 * 4 + 2];
1067         pixels[y0 * linesize + x0 * 4 + 3] = 65535;
1068
1069         if (x0 == x1 && y0 == y1)
1070             break;
1071
1072         e2 = err;
1073
1074         if (e2 >-dx) {
1075             err -= dy;
1076             x0 += sx;
1077         }
1078
1079         if (e2 < dy) {
1080             err += dx;
1081             y0 += sy;
1082         }
1083     }
1084 }
1085
1086 static void
1087 tongue_outline(uint16_t* const pixels,
1088                int       const linesize,
1089                int       const w,
1090                int       const h,
1091                uint16_t  const maxval,
1092                int       const cie)
1093 {
1094     const uint16_t rgbcolor[4] = { maxval, maxval, maxval, maxval };
1095     int wavelength;
1096     int lx, ly;
1097     int fx, fy;
1098
1099     for (wavelength = 360; wavelength <= 830; wavelength++) {
1100         int icx, icy;
1101
1102         monochrome_color_location(wavelength, w, h, cie,
1103                                   &icx, &icy);
1104
1105         if (wavelength > 360)
1106             draw_line(pixels, linesize, lx, ly, icx, icy, w, h, rgbcolor);
1107         else {
1108             fx = icx;
1109             fy = icy;
1110         }
1111         lx = icx;
1112         ly = icy;
1113     }
1114     draw_line(pixels, linesize, lx, ly, fx, fy, w, h, rgbcolor);
1115 }
1116
1117 static void
1118 fill_in_tongue(uint16_t*                  const pixels,
1119                int                        const linesize,
1120                int                        const w,
1121                int                        const h,
1122                uint16_t                   const maxval,
1123                const struct ColorSystem * const cs,
1124                float                      const m[3][3],
1125                int                        const cie,
1126                int                        const correct_gamma,
1127                float                      const contrast)
1128 {
1129     int y;
1130
1131     /* Scan the image line by line and  fill  the  tongue  outline
1132        with the RGB values determined by the color system for the x-y
1133        co-ordinates within the tongue.
1134     */
1135
1136     for (y = 0; y < h; ++y) {
1137         int  present;  /* There is some tongue on this line */
1138         int leftEdge; /* x position of leftmost pixel in tongue on this line */
1139         int rightEdge; /* same, but rightmost */
1140
1141         find_tongue(pixels, w, linesize, y, &present, &leftEdge, &rightEdge);
1142
1143         if (present) {
1144             int x;
1145
1146             for (x = leftEdge; x <= rightEdge; ++x) {
1147                 float cx, cy, cz, jr, jg, jb, jmax;
1148                 int r, g, b, mx = maxval;
1149
1150                 if (cie == LUV) {
1151                     float up, vp;
1152                     up = ((float) x) / (w - 1);
1153                     vp = 1.0 - ((float) y) / (h - 1);
1154                     upvp_to_xy(up, vp, &cx, &cy);
1155                     cz = 1.0 - (cx + cy);
1156                 } else if (cie == UCS) {
1157                     float u, v;
1158                     u = ((float) x) / (w - 1);
1159                     v = 1.0 - ((float) y) / (h - 1);
1160                     uv_to_xy(u, v, &cx, &cy);
1161                     cz = 1.0 - (cx + cy);
1162                 } else if (cie == XYY) {
1163                     cx = ((float) x) / (w - 1);
1164                     cy = 1.0 - ((float) y) / (h - 1);
1165                     cz = 1.0 - (cx + cy);
1166                 } else {
1167                     av_assert0(0);
1168                 }
1169
1170                 xyz_to_rgb(m, cx, cy, cz, &jr, &jg, &jb);
1171
1172                 /* Check whether the requested color  is  within  the
1173                    gamut  achievable with the given color system.  If
1174                    not, draw it in a reduced  intensity,  interpolated
1175                    by desaturation to the closest within-gamut color. */
1176
1177                 if (constrain_rgb(&jr, &jg, &jb))
1178                     mx *= contrast;
1179
1180                 jmax = FFMAX3(jr, jg, jb);
1181                 if (jmax > 0) {
1182                     jr = jr / jmax;
1183                     jg = jg / jmax;
1184                     jb = jb / jmax;
1185                 }
1186                 /* gamma correct from linear rgb to nonlinear rgb. */
1187                 if (correct_gamma)
1188                     gamma_correct_rgb(cs, &jr, &jg, &jb);
1189                 r = mx * jr;
1190                 g = mx * jg;
1191                 b = mx * jb;
1192                 pixels[y * linesize + x * 4 + 0] = r;
1193                 pixels[y * linesize + x * 4 + 1] = g;
1194                 pixels[y * linesize + x * 4 + 2] = b;
1195                 pixels[y * linesize + x * 4 + 3] = 65535;
1196             }
1197         }
1198     }
1199 }
1200
1201 static void
1202 plot_white_point(uint16_t*      pixels,
1203                  int      const linesize,
1204                  int      const w,
1205                  int      const h,
1206                  int      const maxval,
1207                  int      const color_system,
1208                  int      const cie)
1209 {
1210     const struct ColorSystem *cs = &color_systems[color_system];
1211     int wx, wy;
1212
1213     if (cie == LUV) {
1214         float wup, wvp;
1215         xy_to_upvp(cs->xWhite, cs->yWhite, &wup, &wvp);
1216         wx = (w - 1) * wup;
1217         wy = (h - 1) - ((int) ((h - 1) * wvp));
1218     } else if (cie == UCS) {
1219         float wu, wv;
1220         xy_to_uv(cs->xWhite, cs->yWhite, &wu, &wv);
1221         wx = (w - 1) * wu;
1222         wy = (h - 1) - ((int) ((h - 1) * wv));
1223     } else if (cie == XYY) {
1224         wx = (w - 1) * cs->xWhite;
1225         wy = (h - 1) - ((int) ((h - 1) * cs->yWhite));
1226     } else {
1227         av_assert0(0);
1228     }
1229
1230     draw_rline(pixels, linesize,
1231                wx + Sz(3), wy, wx + Sz(10), wy,
1232                w, h);
1233     draw_rline(pixels, linesize,
1234                wx - Sz(3), wy, wx - Sz(10), wy,
1235                w, h);
1236     draw_rline(pixels, linesize,
1237                wx, wy + Sz(3), wx, wy + Sz(10),
1238                w, h);
1239     draw_rline(pixels, linesize,
1240                wx, wy - Sz(3), wx, wy - Sz(10),
1241                w, h);
1242 }
1243
1244 static int draw_background(AVFilterContext *ctx)
1245 {
1246     CiescopeContext *s = ctx->priv;
1247     const struct ColorSystem *cs = &color_systems[s->color_system];
1248     AVFilterLink *outlink = ctx->outputs[0];
1249     int w = s->size;
1250     int h = s->size;
1251     uint16_t *pixels;
1252
1253     if ((s->f = ff_get_video_buffer(outlink, outlink->w, outlink->h)) == NULL)
1254         return AVERROR(ENOMEM);
1255     pixels = (uint16_t *)s->f->data[0];
1256
1257     tongue_outline(pixels, s->f->linesize[0] / 2, w, h, 65535, s->cie);
1258
1259     if (s->fill)
1260         fill_in_tongue(pixels, s->f->linesize[0] / 2, w, h, 65535, cs, (const float (*)[3])s->i, s->cie,
1261                        s->correct_gamma, s->contrast);
1262
1263     return 0;
1264 }
1265
1266 static void filter_rgb48(AVFilterContext *ctx, const uint8_t *ptr,
1267                          ptrdiff_t linesize,
1268                          float *cx, float *cy, int x, int y)
1269 {
1270     CiescopeContext *s = ctx->priv;
1271     const float scale = 1.f / 65535.f;
1272     const uint16_t *src = (const uint16_t*)(ptr + linesize * y + x * 6);
1273     float r = (src[0] + 0.01f) * scale;
1274     float g = (src[1] + 0.01f) * scale;
1275     float b = (src[2] + 0.01f) * scale;
1276     float cz;
1277
1278     rgb_to_xy(r, g, b, cx, cy, &cz, (const float (*)[3])s->m);
1279 }
1280
1281 static void filter_rgba64(AVFilterContext *ctx, const uint8_t *ptr,
1282                           ptrdiff_t linesize,
1283                           float *cx, float *cy, int x, int y)
1284 {
1285     CiescopeContext *s = ctx->priv;
1286     const float scale = 1.f / 65535.f;
1287     const uint16_t *src = (const uint16_t*)(ptr + linesize * y + x * 8);
1288     float r = (src[0] + 0.01f) * scale;
1289     float g = (src[1] + 0.01f) * scale;
1290     float b = (src[2] + 0.01f) * scale;
1291     float cz;
1292
1293     rgb_to_xy(r, g, b, cx, cy, &cz, (const float (*)[3])s->m);
1294 }
1295
1296 static void filter_rgb24(AVFilterContext *ctx, const uint8_t *ptr,
1297                          ptrdiff_t linesize,
1298                          float *cx, float *cy, int x, int y)
1299 {
1300     CiescopeContext *s = ctx->priv;
1301     const float scale = 1.f / 255.f;
1302     const uint8_t *src = ptr + linesize * y + x * 3;
1303     float r = (src[0] + 0.01f) * scale;
1304     float g = (src[1] + 0.01f) * scale;
1305     float b = (src[2] + 0.01f) * scale;
1306     float cz;
1307
1308     rgb_to_xy(r, g, b, cx, cy, &cz, (const float (*)[3])s->m);
1309 }
1310
1311 static void filter_rgba(AVFilterContext *ctx, const uint8_t *ptr,
1312                         ptrdiff_t linesize,
1313                         float *cx, float *cy, int x, int y)
1314 {
1315     CiescopeContext *s = ctx->priv;
1316     const float scale = 1.f / 255.f;
1317     const uint8_t *src = ptr + linesize * y + x * 4;
1318     float r = (src[0] + 0.01f) * scale;
1319     float g = (src[1] + 0.01f) * scale;
1320     float b = (src[2] + 0.01f) * scale;
1321     float cz;
1322
1323     rgb_to_xy(r, g, b, cx, cy, &cz, (const float (*)[3])s->m);
1324 }
1325
1326 static void filter_xyz(AVFilterContext *ctx, const uint8_t *ptr,
1327                        ptrdiff_t linesize,
1328                        float *cx, float *cy, int x, int y)
1329 {
1330     CiescopeContext *s = ctx->priv;
1331     const uint16_t* src = (uint16_t *)(ptr + linesize * y + x * 6);
1332     float lx = s->log2lin[src[0]];
1333     float ly = s->log2lin[src[1]];
1334     float lz = s->log2lin[src[2]];
1335     float sum = lx + ly + lz;
1336
1337     if (sum == 0)
1338         sum = 1;
1339     *cx = lx / sum;
1340     *cy = ly / sum;
1341 }
1342
1343 static void plot_gamuts(uint16_t *pixels, int linesize, int w, int h,
1344                         int cie, int gamuts)
1345 {
1346     int i;
1347
1348     for (i = 0; i < NB_CS; i++) {
1349         const struct ColorSystem *cs = &color_systems[i];
1350         int rx, ry, gx, gy, bx, by;
1351
1352         if (!((1 << i) & gamuts))
1353             continue;
1354         if (cie == LUV) {
1355             float wup, wvp;
1356             xy_to_upvp(cs->xRed, cs->yRed, &wup, &wvp);
1357             rx = (w - 1) * wup;
1358             ry = (h - 1) - ((int) ((h - 1) * wvp));
1359             xy_to_upvp(cs->xGreen, cs->yGreen, &wup, &wvp);
1360             gx = (w - 1) * wup;
1361             gy = (h - 1) - ((int) ((h - 1) * wvp));
1362             xy_to_upvp(cs->xBlue, cs->yBlue, &wup, &wvp);
1363             bx = (w - 1) * wup;
1364             by = (h - 1) - ((int) ((h - 1) * wvp));
1365         } else if (cie == UCS) {
1366             float wu, wv;
1367             xy_to_uv(cs->xRed, cs->yRed, &wu, &wv);
1368             rx = (w - 1) * wu;
1369             ry = (h - 1) - ((int) ((h - 1) * wv));
1370             xy_to_uv(cs->xGreen, cs->yGreen, &wu, &wv);
1371             gx = (w - 1) * wu;
1372             gy = (h - 1) - ((int) ((h - 1) * wv));
1373             xy_to_uv(cs->xBlue, cs->yBlue, &wu, &wv);
1374             bx = (w - 1) * wu;
1375             by = (h - 1) - ((int) ((h - 1) * wv));
1376         } else if (cie == XYY) {
1377             rx = (w - 1) * cs->xRed;
1378             ry = (h - 1) - ((int) ((h - 1) * cs->yRed));
1379             gx = (w - 1) * cs->xGreen;
1380             gy = (h - 1) - ((int) ((h - 1) * cs->yGreen));
1381             bx = (w - 1) * cs->xBlue;
1382             by = (h - 1) - ((int) ((h - 1) * cs->yBlue));
1383         } else {
1384             av_assert0(0);
1385         }
1386
1387         draw_rline(pixels, linesize, rx, ry, gx, gy, w, h);
1388         draw_rline(pixels, linesize, gx, gy, bx, by, w, h);
1389         draw_rline(pixels, linesize, bx, by, rx, ry, w, h);
1390     }
1391 }
1392
1393 static int filter_frame(AVFilterLink *inlink, AVFrame *in)
1394 {
1395     AVFilterContext *ctx  = inlink->dst;
1396     CiescopeContext *s = ctx->priv;
1397     AVFilterLink *outlink = ctx->outputs[0];
1398     int i = s->intensity * 65535;
1399     int w = outlink->w;
1400     int h = outlink->h;
1401     AVFrame *out;
1402     int ret, x, y;
1403
1404     out = ff_get_video_buffer(outlink, outlink->w, outlink->h);
1405     if (!out) {
1406         av_frame_free(&in);
1407         return AVERROR(ENOMEM);
1408     }
1409     out->pts = in->pts;
1410     out->duration = in->duration;
1411
1412     if (!s->background) {
1413         ret = draw_background(ctx);
1414         if (ret < 0) {
1415             av_frame_free(&out);
1416             return ret;
1417         }
1418         s->background = 1;
1419     }
1420     for (y = 0; y < outlink->h; y++) {
1421         memset(out->data[0] + y * out->linesize[0], 0, outlink->w * 8);
1422     }
1423
1424     for (y = 0; y < in->height; y++) {
1425         const uint8_t *src = in->data[0];
1426         const ptrdiff_t src_linesize = in->linesize[0];
1427         uint16_t *dst = (uint16_t *)out->data[0];
1428         const ptrdiff_t linesize = out->linesize[0] / 2;
1429         const int w_1 = w - 1;
1430         const int h_1 = h - 1;
1431
1432         for (x = 0; x < in->width; x++) {
1433             float cx, cy;
1434             int wx, wy, pos;
1435             int r, g, b;
1436
1437             s->filter(ctx, src, src_linesize, &cx, &cy, x, y);
1438
1439             if (s->cie == LUV) {
1440                 float up, vp;
1441                 xy_to_upvp(cx, cy, &up, &vp);
1442                 cx = up;
1443                 cy = vp;
1444             } else if (s->cie == UCS) {
1445                 float u, v;
1446                 xy_to_uv(cx, cy, &u, &v);
1447                 cx = u;
1448                 cy = v;
1449             }
1450
1451             wx = w_1 * cx;
1452             wy = h_1 - h_1 * cy;
1453
1454             if (wx < 0 || wx >= w ||
1455                 wy < 0 || wy >= h)
1456                 continue;
1457
1458             pos = wy * linesize + wx * 4;
1459             r = dst[pos + 0] + i;
1460             g = dst[pos + 1] + i;
1461             b = dst[pos + 2] + i;
1462
1463             dst[pos + 0] = FFMIN(r, 65535);
1464             dst[pos + 1] = FFMIN(g, 65535);
1465             dst[pos + 2] = FFMIN(b, 65535);
1466             dst[pos + 3] = 65535;
1467         }
1468     }
1469
1470     for (y = 0; y < outlink->h; y++) {
1471         uint16_t *dst = (uint16_t *)(out->data[0] + y * out->linesize[0]);
1472         const uint16_t *src = (const uint16_t *)(s->f->data[0] + y * s->f->linesize[0]);
1473         for (x = 0; x < outlink->w; x++) {
1474             const int xx = x * 4;
1475             if (dst[xx + 3] == 0) {
1476                 dst[xx + 0] = src[xx + 0];
1477                 dst[xx + 1] = src[xx + 1];
1478                 dst[xx + 2] = src[xx + 2];
1479                 dst[xx + 3] = src[xx + 3];
1480             }
1481         }
1482     }
1483
1484     if (s->show_white)
1485         plot_white_point((uint16_t *)out->data[0], out->linesize[0] / 2,
1486                          outlink->w, outlink->h, 65535,
1487                          s->color_system, s->cie);
1488
1489     plot_gamuts((uint16_t *)out->data[0], out->linesize[0] / 2,
1490                 outlink->w, outlink->h,
1491                 s->cie, s->gamuts);
1492
1493     av_frame_free(&in);
1494     return ff_filter_frame(outlink, out);
1495 }
1496
1497 static void av_cold uninit(AVFilterContext *ctx)
1498 {
1499     CiescopeContext *s = ctx->priv;
1500
1501     av_frame_free(&s->f);
1502 }
1503
1504 static int config_input(AVFilterLink *inlink)
1505 {
1506     CiescopeContext *s = inlink->dst->priv;
1507     int i;
1508
1509     get_rgb2xyz_matrix(color_systems[s->color_system], s->m);
1510     invert_matrix3x3(s->m, s->i);
1511
1512     switch (inlink->format) {
1513     case AV_PIX_FMT_RGB24:
1514         s->filter = filter_rgb24;
1515         break;
1516     case AV_PIX_FMT_RGBA:
1517         s->filter = filter_rgba;
1518         break;
1519     case AV_PIX_FMT_RGB48:
1520         s->filter = filter_rgb48;
1521         break;
1522     case AV_PIX_FMT_RGBA64:
1523         s->filter = filter_rgba64;
1524         break;
1525     case AV_PIX_FMT_XYZ12:
1526         s->filter = filter_xyz;
1527         for (i = 0; i < 65536; i++)
1528             s->log2lin[i] = pow(i / 65535., s->igamma) * 65535.;
1529         break;
1530     default:
1531         av_assert0(0);
1532     }
1533
1534     return 0;
1535 }
1536
1537 static const AVFilterPad inputs[] = {
1538     {
1539         .name         = "default",
1540         .type         = AVMEDIA_TYPE_VIDEO,
1541         .filter_frame = filter_frame,
1542         .config_props = config_input,
1543     },
1544 };
1545
1546 static const AVFilterPad outputs[] = {
1547     {
1548         .name         = "default",
1549         .type         = AVMEDIA_TYPE_VIDEO,
1550         .config_props = config_output,
1551     },
1552 };
1553
1554 const AVFilter ff_vf_ciescope = {
1555     .name          = "ciescope",
1556     .description   = NULL_IF_CONFIG_SMALL("Video CIE scope."),
1557     .priv_size     = sizeof(CiescopeContext),
1558     .priv_class    = &ciescope_class,
1559     .uninit        = uninit,
1560     FILTER_INPUTS(inputs),
1561     FILTER_OUTPUTS(outputs),
1562     FILTER_QUERY_FUNC(query_formats),
1563 };