1 /*
2 * Copyright (c) 2016 ReneBrals
3 * Copyright (c) 2021 Paul B Mahol
4 *
5 * This file is part of FFmpeg.
6 *
7 * Permission is hereby granted, free of charge, to any person obtaining a copy
8 * of this software and associated documentation files (the "Software"), to deal
9 * in the Software without restriction, including without limitation the rights
10 * to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
11 * copies of the Software, and to permit persons to whom the Software is
12 * furnished to do so, subject to the following conditions:
13 *
14 * The above copyright notice and this permission notice shall be included in all
15 * copies or substantial portions of the Software.
16 *
17 * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
18 * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
19 * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
20 * AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
21 * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
22 * OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
23 * SOFTWARE.
24 */
25
26 #include "libavutil/avassert.h"
27 #include "libavutil/imgutils.h"
28 #include "libavutil/intreadwrite.h"
29 #include "libavutil/opt.h"
30 #include "libavutil/pixdesc.h"
31 #include "avfilter.h"
32 #include "formats.h"
33 #include "framesync.h"
34 #include "internal.h"
35 #include "video.h"
36
37 enum MorphModes {
38 ERODE,
39 DILATE,
40 OPEN,
41 CLOSE,
42 GRADIENT,
43 TOPHAT,
44 BLACKHAT,
45 NB_MODES
46 };
47
48 typedef struct IPlane {
49 uint8_t **img;
50 int w, h;
51 int range;
52 int depth;
53 int type_size;
54
55 void (*max_out_place)(uint8_t *c, const uint8_t *a, const uint8_t *b, int x);
56 void (*min_out_place)(uint8_t *c, const uint8_t *a, const uint8_t *b, int x);
57 void (*diff_rin_place)(uint8_t *a, const uint8_t *b, int x);
58 void (*max_in_place)(uint8_t *a, const uint8_t *b, int x);
59 void (*min_in_place)(uint8_t *a, const uint8_t *b, int x);
60 void (*diff_in_place)(uint8_t *a, const uint8_t *b, int x);
61 } IPlane;
62
63 typedef struct LUT {
64 /* arr is shifted from base_arr by FFMAX(min_r, 0).
65 * arr != NULL means "lut completely allocated" */
66 uint8_t ***arr;
67 uint8_t ***base_arr;
68 int min_r;
69 int max_r;
70 int I;
71 int X;
72 int pre_pad_x;
73 int type_size;
74 } LUT;
75
76 typedef struct chord {
77 int x;
78 int y;
79 int l;
80 int i;
81 } chord;
82
83 typedef struct chord_set {
84 chord *C;
85 int size;
86 int cap;
87
88 int *R;
89 int Lnum;
90
91 int minX;
92 int maxX;
93 int minY;
94 int maxY;
95 unsigned nb_elements;
96 } chord_set;
97
98 typedef struct MorphoContext {
99 const AVClass *class;
100 FFFrameSync fs;
101
102 chord_set SE[4];
103 IPlane SEimg[4];
104 IPlane g[4], f[4], h[4];
105 LUT Ty[2][4];
106
107 int mode;
108 int planes;
109 int structures;
110
111 int planewidth[4];
112 int planeheight[4];
113 int splanewidth[4];
114 int splaneheight[4];
115 int depth;
116 int type_size;
117 int nb_planes;
118
119 int got_structure[4];
120
121 AVFrame *temp;
122
123 int64_t *plane_f, *plane_g;
124 } MorphoContext;
125
126 #define OFFSET(x) offsetof(MorphoContext, x)
127 #define FLAGS AV_OPT_FLAG_VIDEO_PARAM | AV_OPT_FLAG_FILTERING_PARAM | AV_OPT_FLAG_RUNTIME_PARAM
128
129 static const AVOption morpho_options[] = {
130 { "mode", "set morphological transform", OFFSET(mode), AV_OPT_TYPE_INT, {.i64=0}, 0, NB_MODES-1, FLAGS, "mode" },
131 { "erode", NULL, 0, AV_OPT_TYPE_CONST, {.i64=ERODE}, 0, 0, FLAGS, "mode" },
132 { "dilate", NULL, 0, AV_OPT_TYPE_CONST, {.i64=DILATE}, 0, 0, FLAGS, "mode" },
133 { "open", NULL, 0, AV_OPT_TYPE_CONST, {.i64=OPEN}, 0, 0, FLAGS, "mode" },
134 { "close", NULL, 0, AV_OPT_TYPE_CONST, {.i64=CLOSE}, 0, 0, FLAGS, "mode" },
135 { "gradient",NULL, 0, AV_OPT_TYPE_CONST, {.i64=GRADIENT},0, 0, FLAGS, "mode" },
136 { "tophat",NULL, 0, AV_OPT_TYPE_CONST, {.i64=TOPHAT}, 0, 0, FLAGS, "mode" },
137 { "blackhat",NULL, 0, AV_OPT_TYPE_CONST, {.i64=BLACKHAT},0, 0, FLAGS, "mode" },
138 { "planes", "set planes to filter", OFFSET(planes), AV_OPT_TYPE_INT, {.i64=7}, 0, 15, FLAGS },
139 { "structure", "when to process structures", OFFSET(structures), AV_OPT_TYPE_INT, {.i64=1}, 0, 1, FLAGS, "str" },
140 { "first", "process only first structure, ignore rest", 0, AV_OPT_TYPE_CONST, {.i64=0}, 0, 0, FLAGS, "str" },
141 { "all", "process all structure", 0, AV_OPT_TYPE_CONST, {.i64=1}, 0, 0, FLAGS, "str" },
142 { NULL }
143 };
144
145 FRAMESYNC_DEFINE_CLASS(morpho, MorphoContext, fs);
146
147 static const enum AVPixelFormat pix_fmts[] = {
148 AV_PIX_FMT_YUVA444P, AV_PIX_FMT_YUV444P, AV_PIX_FMT_YUV440P,
149 AV_PIX_FMT_YUVJ444P, AV_PIX_FMT_YUVJ440P,
150 AV_PIX_FMT_YUVA422P, AV_PIX_FMT_YUV422P, AV_PIX_FMT_YUVA420P, AV_PIX_FMT_YUV420P,
151 AV_PIX_FMT_YUVJ422P, AV_PIX_FMT_YUVJ420P,
152 AV_PIX_FMT_YUVJ411P, AV_PIX_FMT_YUV411P, AV_PIX_FMT_YUV410P,
153 AV_PIX_FMT_GBRP, AV_PIX_FMT_GBRAP, AV_PIX_FMT_GRAY8, AV_PIX_FMT_GRAY9,
154 AV_PIX_FMT_YUV420P9, AV_PIX_FMT_YUV422P9, AV_PIX_FMT_YUV444P9, AV_PIX_FMT_GBRP9,
155 AV_PIX_FMT_YUVA420P9, AV_PIX_FMT_YUVA422P9, AV_PIX_FMT_YUVA444P9,
156 AV_PIX_FMT_YUV420P10, AV_PIX_FMT_YUV422P10, AV_PIX_FMT_YUV444P10,
157 AV_PIX_FMT_YUV420P12, AV_PIX_FMT_YUV422P12, AV_PIX_FMT_YUV444P12, AV_PIX_FMT_YUV440P12,
158 AV_PIX_FMT_YUV420P14, AV_PIX_FMT_YUV422P14, AV_PIX_FMT_YUV444P14,
159 AV_PIX_FMT_YUV420P16, AV_PIX_FMT_YUV422P16, AV_PIX_FMT_YUV444P16,
160 AV_PIX_FMT_YUVA420P10, AV_PIX_FMT_YUVA422P10, AV_PIX_FMT_YUVA444P10,
161 AV_PIX_FMT_YUVA422P12, AV_PIX_FMT_YUVA444P12,
162 AV_PIX_FMT_YUVA420P16, AV_PIX_FMT_YUVA422P16, AV_PIX_FMT_YUVA444P16,
163 AV_PIX_FMT_GBRP10, AV_PIX_FMT_GBRP12, AV_PIX_FMT_GBRP14, AV_PIX_FMT_GBRP16,
164 AV_PIX_FMT_GBRAP10, AV_PIX_FMT_GBRAP12, AV_PIX_FMT_GBRAP16,
165 AV_PIX_FMT_GRAY10, AV_PIX_FMT_GRAY12, AV_PIX_FMT_GRAY14, AV_PIX_FMT_GRAY16,
166 AV_PIX_FMT_NONE
167 };
168
min_fun(uint8_t * c,const uint8_t * a,const uint8_t * b,int x)169 static void min_fun(uint8_t *c, const uint8_t *a, const uint8_t *b, int x)
170 {
171 for (int i = 0; i < x; i++)
172 c[i] = FFMIN(b[i], a[i]);
173 }
174
mininplace_fun(uint8_t * a,const uint8_t * b,int x)175 static void mininplace_fun(uint8_t *a, const uint8_t *b, int x)
176 {
177 for (int i = 0; i < x; i++)
178 a[i] = FFMIN(a[i], b[i]);
179 }
180
max_fun(uint8_t * c,const uint8_t * a,const uint8_t * b,int x)181 static void max_fun(uint8_t *c, const uint8_t *a, const uint8_t *b, int x)
182 {
183 for (int i = 0; i < x; i++)
184 c[i] = FFMAX(a[i], b[i]);
185 }
186
maxinplace_fun(uint8_t * a,const uint8_t * b,int x)187 static void maxinplace_fun(uint8_t *a, const uint8_t *b, int x)
188 {
189 for (int i = 0; i < x; i++)
190 a[i] = FFMAX(a[i], b[i]);
191 }
192
diff_fun(uint8_t * a,const uint8_t * b,int x)193 static void diff_fun(uint8_t *a, const uint8_t *b, int x)
194 {
195 for (int i = 0; i < x; i++)
196 a[i] = FFMAX(b[i] - a[i], 0);
197 }
198
diffinplace_fun(uint8_t * a,const uint8_t * b,int x)199 static void diffinplace_fun(uint8_t *a, const uint8_t *b, int x)
200 {
201 for (int i = 0; i < x; i++)
202 a[i] = FFMAX(a[i] - b[i], 0);
203 }
204
min16_fun(uint8_t * cc,const uint8_t * aa,const uint8_t * bb,int x)205 static void min16_fun(uint8_t *cc, const uint8_t *aa, const uint8_t *bb, int x)
206 {
207 const uint16_t *a = (const uint16_t *)aa;
208 const uint16_t *b = (const uint16_t *)bb;
209 uint16_t *c = (uint16_t *)cc;
210
211 for (int i = 0; i < x; i++)
212 c[i] = FFMIN(b[i], a[i]);
213 }
214
mininplace16_fun(uint8_t * aa,const uint8_t * bb,int x)215 static void mininplace16_fun(uint8_t *aa, const uint8_t *bb, int x)
216 {
217 uint16_t *a = (uint16_t *)aa;
218 const uint16_t *b = (const uint16_t *)bb;
219
220 for (int i = 0; i < x; i++)
221 a[i] = FFMIN(a[i], b[i]);
222 }
223
diff16_fun(uint8_t * aa,const uint8_t * bb,int x)224 static void diff16_fun(uint8_t *aa, const uint8_t *bb, int x)
225 {
226 const uint16_t *b = (const uint16_t *)bb;
227 uint16_t *a = (uint16_t *)aa;
228
229 for (int i = 0; i < x; i++)
230 a[i] = FFMAX(b[i] - a[i], 0);
231 }
232
diffinplace16_fun(uint8_t * aa,const uint8_t * bb,int x)233 static void diffinplace16_fun(uint8_t *aa, const uint8_t *bb, int x)
234 {
235 uint16_t *a = (uint16_t *)aa;
236 const uint16_t *b = (const uint16_t *)bb;
237
238 for (int i = 0; i < x; i++)
239 a[i] = FFMAX(a[i] - b[i], 0);
240 }
241
max16_fun(uint8_t * cc,const uint8_t * aa,const uint8_t * bb,int x)242 static void max16_fun(uint8_t *cc, const uint8_t *aa, const uint8_t *bb, int x)
243 {
244 const uint16_t *a = (const uint16_t *)aa;
245 const uint16_t *b = (const uint16_t *)bb;
246 uint16_t *c = (uint16_t *)cc;
247
248 for (int i = 0; i < x; i++)
249 c[i] = FFMAX(a[i], b[i]);
250 }
251
maxinplace16_fun(uint8_t * aa,const uint8_t * bb,int x)252 static void maxinplace16_fun(uint8_t *aa, const uint8_t *bb, int x)
253 {
254 uint16_t *a = (uint16_t *)aa;
255 const uint16_t *b = (const uint16_t *)bb;
256
257 for (int i = 0; i < x; i++)
258 a[i] = FFMAX(a[i], b[i]);
259 }
260
alloc_lut(LUT * Ty,chord_set * SE,int type_size,int mode)261 static int alloc_lut(LUT *Ty, chord_set *SE, int type_size, int mode)
262 {
263 const int min = FFMAX(Ty->min_r, 0);
264 const int max = min + (Ty->max_r - Ty->min_r);
265 int pre_pad_x = 0;
266
267 if (SE->minX < 0)
268 pre_pad_x = 0 - SE->minX;
269 Ty->pre_pad_x = pre_pad_x;
270 Ty->type_size = type_size;
271
272 Ty->base_arr = av_calloc(max + 1, sizeof(*Ty->base_arr));
273 if (!Ty->base_arr)
274 return AVERROR(ENOMEM);
275 for (int r = min; r <= max; r++) {
276 uint8_t **arr = Ty->base_arr[r] = av_calloc(Ty->I, sizeof(uint8_t *));
277 if (!Ty->base_arr[r])
278 return AVERROR(ENOMEM);
279 for (int i = 0; i < Ty->I; i++) {
280 arr[i] = av_calloc(Ty->X + pre_pad_x, type_size);
281 if (!arr[i])
282 return AVERROR(ENOMEM);
283 if (mode == ERODE)
284 memset(arr[i], UINT8_MAX, pre_pad_x * type_size);
285 /* Shifting the X index such that negative indices correspond to
286 * the pre-padding.
287 */
288 arr[i] = &(arr[i][pre_pad_x * type_size]);
289 }
290 }
291
292 Ty->arr = &(Ty->base_arr[min - Ty->min_r]);
293
294 return 0;
295 }
296
free_lut(LUT * table)297 static void free_lut(LUT *table)
298 {
299 const int min = FFMAX(table->min_r, 0);
300 const int max = min + (table->max_r - table->min_r);
301
302 if (!table->base_arr)
303 return;
304
305 for (int r = min; r <= max; r++) {
306 if (!table->base_arr[r])
307 break;
308 for (int i = 0; i < table->I; i++) {
309 if (!table->base_arr[r][i])
310 break;
311 // The X index was also shifted, for padding purposes.
312 av_free(table->base_arr[r][i] - table->pre_pad_x * table->type_size);
313 }
314 av_freep(&table->base_arr[r]);
315 }
316 av_freep(&table->base_arr);
317 table->arr = NULL;
318 }
319
alloc_lut_if_necessary(LUT * Ty,IPlane * f,chord_set * SE,int y,int num,enum MorphModes mode)320 static int alloc_lut_if_necessary(LUT *Ty, IPlane *f, chord_set *SE,
321 int y, int num, enum MorphModes mode)
322 {
323 if (!Ty->arr || Ty->I != SE->Lnum ||
324 Ty->X != f->w ||
325 SE->minX < 0 && -SE->minX > Ty->pre_pad_x ||
326 Ty->min_r != SE->minY ||
327 Ty->max_r != SE->maxY + num - 1) {
328 int ret;
329
330 free_lut(Ty);
331
332 Ty->I = SE->Lnum;
333 Ty->X = f->w;
334 Ty->min_r = SE->minY;
335 Ty->max_r = SE->maxY + num - 1;
336 ret = alloc_lut(Ty, SE, f->type_size, mode);
337 if (ret < 0)
338 return ret;
339 }
340 return 0;
341 }
342
circular_swap(LUT * Ty)343 static void circular_swap(LUT *Ty)
344 {
345 /*
346 * Swap the pointers to r-indices in a circle. This is useful because
347 * Ty(r,i,x) = Ty-1(r+1,i,x) for r < ymax.
348 */
349 if (Ty->max_r - Ty->min_r > 0) {
350 uint8_t **Ty0 = Ty->arr[Ty->min_r];
351
352 for (int r = Ty->min_r; r < Ty->max_r; r++)
353 Ty->arr[r] = Ty->arr[r + 1];
354
355 Ty->arr[Ty->max_r] = Ty0;
356 }
357 }
358
compute_min_row(IPlane * f,LUT * Ty,chord_set * SE,int r,int y)359 static void compute_min_row(IPlane *f, LUT *Ty, chord_set *SE, int r, int y)
360 {
361 if (y + r >= 0 && y + r < f->h) {
362 memcpy(Ty->arr[r][0], f->img[y + r], Ty->X * Ty->type_size);
363 } else {
364 memset(Ty->arr[r][0], UINT8_MAX, Ty->X * Ty->type_size);
365 }
366
367 for (int i = 1; i < SE->Lnum; i++) {
368 int d = SE->R[i] - SE->R[i - 1];
369
370 f->min_out_place(Ty->arr[r][i] - Ty->pre_pad_x * f->type_size,
371 Ty->arr[r][i - 1] - Ty->pre_pad_x * f->type_size,
372 Ty->arr[r][i - 1] + (d - Ty->pre_pad_x) * f->type_size,
373 Ty->X + Ty->pre_pad_x - d);
374 memcpy(Ty->arr[r][i] + (Ty->X - d) * f->type_size,
375 Ty->arr[r][i - 1] + (Ty->X - d) * f->type_size,
376 d * f->type_size);
377 }
378 }
379
update_min_lut(IPlane * f,LUT * Ty,chord_set * SE,int y,int tid,int num)380 static void update_min_lut(IPlane *f, LUT *Ty, chord_set *SE, int y, int tid, int num)
381 {
382 for (int i = 0; i < num; i++)
383 circular_swap(Ty);
384
385 compute_min_row(f, Ty, SE, Ty->max_r - tid, y);
386 }
387
compute_min_lut(LUT * Ty,IPlane * f,chord_set * SE,int y,int num)388 static int compute_min_lut(LUT *Ty, IPlane *f, chord_set *SE, int y, int num)
389 {
390 int ret = alloc_lut_if_necessary(Ty, f, SE, y, num, ERODE);
391 if (ret < 0)
392 return ret;
393
394 for (int r = Ty->min_r; r <= Ty->max_r; r++)
395 compute_min_row(f, Ty, SE, r, y);
396
397 return 0;
398 }
399
compute_max_row(IPlane * f,LUT * Ty,chord_set * SE,int r,int y)400 static void compute_max_row(IPlane *f, LUT *Ty, chord_set *SE, int r, int y)
401 {
402 if (y + r >= 0 && y + r < f->h) {
403 memcpy(Ty->arr[r][0], f->img[y + r], Ty->X * Ty->type_size);
404 } else {
405 memset(Ty->arr[r][0], 0, Ty->X * Ty->type_size);
406 }
407
408 for (int i = 1; i < SE->Lnum; i++) {
409 int d = SE->R[i] - SE->R[i - 1];
410
411 f->max_out_place(Ty->arr[r][i] - Ty->pre_pad_x * f->type_size,
412 Ty->arr[r][i - 1] - Ty->pre_pad_x * f->type_size,
413 Ty->arr[r][i - 1] + (d - Ty->pre_pad_x) * f->type_size,
414 Ty->X + Ty->pre_pad_x - d);
415 memcpy(Ty->arr[r][i] + (Ty->X - d) * f->type_size,
416 Ty->arr[r][i - 1] + (Ty->X - d) * f->type_size,
417 d * f->type_size);
418 }
419 }
420
update_max_lut(IPlane * f,LUT * Ty,chord_set * SE,int y,int tid,int num)421 static void update_max_lut(IPlane *f, LUT *Ty, chord_set *SE, int y, int tid, int num)
422 {
423 for (int i = 0; i < num; i++)
424 circular_swap(Ty);
425
426 compute_max_row(f, Ty, SE, Ty->max_r - tid, y);
427 }
428
compute_max_lut(LUT * Ty,IPlane * f,chord_set * SE,int y,int num)429 static int compute_max_lut(LUT *Ty, IPlane *f, chord_set *SE, int y, int num)
430 {
431 int ret = alloc_lut_if_necessary(Ty, f, SE, y, num, DILATE);
432 if (ret < 0)
433 return ret;
434
435 for (int r = Ty->min_r; r <= Ty->max_r; r++)
436 compute_max_row(f, Ty, SE, r, y);
437
438 return 0;
439 }
440
line_dilate(IPlane * g,LUT * Ty,chord_set * SE,int y,int tid)441 static void line_dilate(IPlane *g, LUT *Ty, chord_set *SE, int y, int tid)
442 {
443 memset(g->img[y], 0, g->w * g->type_size);
444
445 for (int c = 0; c < SE->size; c++) {
446 g->max_in_place(g->img[y],
447 Ty->arr[SE->C[c].y + tid][SE->C[c].i] + SE->C[c].x * Ty->type_size,
448 av_clip(g->w - SE->C[c].x, 0, g->w));
449 }
450 }
451
line_erode(IPlane * g,LUT * Ty,chord_set * SE,int y,int tid)452 static void line_erode(IPlane *g, LUT *Ty, chord_set *SE, int y, int tid)
453 {
454 memset(g->img[y], UINT8_MAX, g->w * g->type_size);
455
456 for (int c = 0; c < SE->size; c++) {
457 g->min_in_place(g->img[y],
458 Ty->arr[SE->C[c].y + tid][SE->C[c].i] + SE->C[c].x * Ty->type_size,
459 av_clip(g->w - SE->C[c].x, 0, g->w));
460 }
461 }
462
dilate(IPlane * g,IPlane * f,chord_set * SE,LUT * Ty)463 static int dilate(IPlane *g, IPlane *f, chord_set *SE, LUT *Ty)
464 {
465 int ret = compute_max_lut(Ty, f, SE, 0, 1);
466 if (ret < 0)
467 return ret;
468
469 line_dilate(g, Ty, SE, 0, 0);
470 for (int y = 1; y < f->h; y++) {
471 update_max_lut(f, Ty, SE, y, 0, 1);
472 line_dilate(g, Ty, SE, y, 0);
473 }
474
475 return 0;
476 }
477
erode(IPlane * g,IPlane * f,chord_set * SE,LUT * Ty)478 static int erode(IPlane *g, IPlane *f, chord_set *SE, LUT *Ty)
479 {
480 int ret = compute_min_lut(Ty, f, SE, 0, 1);
481 if (ret < 0)
482 return ret;
483
484 line_erode(g, Ty, SE, 0, 0);
485 for (int y = 1; y < f->h; y++) {
486 update_min_lut(f, Ty, SE, y, 0, 1);
487 line_erode(g, Ty, SE, y, 0);
488 }
489
490 return 0;
491 }
492
difference(IPlane * g,IPlane * f)493 static void difference(IPlane *g, IPlane *f)
494 {
495 for (int y = 0; y < f->h; y++)
496 f->diff_in_place(g->img[y], f->img[y], f->w);
497 }
498
difference2(IPlane * g,IPlane * f)499 static void difference2(IPlane *g, IPlane *f)
500 {
501 for (int y = 0; y < f->h; y++)
502 f->diff_rin_place(g->img[y], f->img[y], f->w);
503 }
504
insert_chord_set(chord_set * chords,chord c)505 static int insert_chord_set(chord_set *chords, chord c)
506 {
507 // Checking if chord fits in dynamic array, resize if not.
508 if (chords->size == chords->cap) {
509 chords->C = av_realloc_f(chords->C, chords->cap * 2, sizeof(chord));
510 if (!chords->C)
511 return AVERROR(ENOMEM);
512 chords->cap *= 2;
513 }
514
515 // Add the chord to the dynamic array.
516 chords->C[chords->size].x = c.x;
517 chords->C[chords->size].y = c.y;
518 chords->C[chords->size++].l = c.l;
519
520 // Update minimum/maximum x/y offsets of the chord set.
521 chords->minX = FFMIN(chords->minX, c.x);
522 chords->maxX = FFMAX(chords->maxX, c.x);
523
524 chords->minY = FFMIN(chords->minY, c.y);
525 chords->maxY = FFMAX(chords->maxY, c.y);
526
527 return 0;
528 }
529
free_chord_set(chord_set * SE)530 static void free_chord_set(chord_set *SE)
531 {
532 av_freep(&SE->C);
533 SE->size = 0;
534 SE->cap = 0;
535
536 av_freep(&SE->R);
537 SE->Lnum = 0;
538 }
539
init_chordset(chord_set * chords)540 static int init_chordset(chord_set *chords)
541 {
542 chords->nb_elements = 0;
543 chords->size = 0;
544 chords->C = av_calloc(1, sizeof(chord));
545 if (!chords->C)
546 return AVERROR(ENOMEM);
547
548 chords->cap = 1;
549 chords->minX = INT16_MAX;
550 chords->maxX = INT16_MIN;
551 chords->minY = INT16_MAX;
552 chords->maxY = INT16_MIN;
553
554 return 0;
555 }
556
comp_chord_length(const void * p,const void * q)557 static int comp_chord_length(const void *p, const void *q)
558 {
559 chord a, b;
560 a = *((chord *)p);
561 b = *((chord *)q);
562
563 return (a.l > b.l) - (a.l < b.l);
564 }
565
comp_chord(const void * p,const void * q)566 static int comp_chord(const void *p, const void *q)
567 {
568 chord a, b;
569 a = *((chord *)p);
570 b = *((chord *)q);
571
572 return (a.y > b.y) - (a.y < b.y);
573 }
574
build_chord_set(IPlane * SE,chord_set * chords)575 static int build_chord_set(IPlane *SE, chord_set *chords)
576 {
577 const int mid = 1 << (SE->depth - 1);
578 int chord_length_index;
579 int chord_start, val, ret;
580 int centerX, centerY;
581 int r_cap = 1;
582 chord c;
583
584 ret = init_chordset(chords);
585 if (ret < 0)
586 return ret;
587 /*
588 * In erosion/dilation, the center of the IPlane has S.E. offset (0,0).
589 * Otherwise, the resulting IPlane would be shifted to the top-left.
590 */
591 centerX = (SE->w - 1) / 2;
592 centerY = (SE->h - 1) / 2;
593
594 /*
595 * Computing the set of chords C.
596 */
597 for (int y = 0; y < SE->h; y++) {
598 int x;
599
600 chord_start = -1;
601 for (x = 0; x < SE->w; x++) {
602 if (SE->type_size == 1) {
603 chords->nb_elements += (SE->img[y][x] >= mid);
604 //A chord is a run of non-zero pixels.
605 if (SE->img[y][x] >= mid && chord_start == -1) {
606 // Chord starts.
607 chord_start = x;
608 } else if (SE->img[y][x] < mid && chord_start != -1) {
609 // Chord ends before end of line.
610 c.x = chord_start - centerX;
611 c.y = y - centerY;
612 c.l = x - chord_start;
613 ret = insert_chord_set(chords, c);
614 if (ret < 0)
615 return AVERROR(ENOMEM);
616 chord_start = -1;
617 }
618 } else {
619 chords->nb_elements += (AV_RN16(&SE->img[y][x * 2]) >= mid);
620 //A chord is a run of non-zero pixels.
621 if (AV_RN16(&SE->img[y][x * 2]) >= mid && chord_start == -1) {
622 // Chord starts.
623 chord_start = x;
624 } else if (AV_RN16(&SE->img[y][x * 2]) < mid && chord_start != -1) {
625 // Chord ends before end of line.
626 c.x = chord_start - centerX;
627 c.y = y - centerY;
628 c.l = x - chord_start;
629 ret = insert_chord_set(chords, c);
630 if (ret < 0)
631 return AVERROR(ENOMEM);
632 chord_start = -1;
633 }
634 }
635 }
636 if (chord_start != -1) {
637 // Chord ends at end of line.
638 c.x = chord_start - centerX;
639 c.y = y - centerY;
640 c.l = x - chord_start;
641 ret = insert_chord_set(chords, c);
642 if (ret < 0)
643 return AVERROR(ENOMEM);
644 }
645 }
646
647 /*
648 * Computing the array of chord lengths R(i).
649 * This is needed because the lookup table will contain a row for each
650 * length index i.
651 */
652 qsort(chords->C, chords->size, sizeof(chord), comp_chord_length);
653 chords->R = av_calloc(1, sizeof(*chords->R));
654 if (!chords->R)
655 return AVERROR(ENOMEM);
656 chords->Lnum = 0;
657 val = 0;
658 r_cap = 1;
659
660 if (chords->size > 0) {
661 val = 1;
662 if (chords->Lnum >= r_cap) {
663 chords->R = av_realloc_f(chords->R, r_cap * 2, sizeof(*chords->R));
664 if (!chords->R)
665 return AVERROR(ENOMEM);
666 r_cap *= 2;
667 }
668 chords->R[chords->Lnum++] = 1;
669 val = 1;
670 }
671
672 for (int i = 0; i < chords->size; i++) {
673 if (val != chords->C[i].l) {
674 while (2 * val < chords->C[i].l && val != 0) {
675 if (chords->Lnum >= r_cap) {
676 chords->R = av_realloc_f(chords->R, r_cap * 2, sizeof(*chords->R));
677 if (!chords->R)
678 return AVERROR(ENOMEM);
679 r_cap *= 2;
680 }
681
682 chords->R[chords->Lnum++] = 2 * val;
683 val *= 2;
684 }
685 val = chords->C[i].l;
686
687 if (chords->Lnum >= r_cap) {
688 chords->R = av_realloc_f(chords->R, r_cap * 2, sizeof(*chords->R));
689 if (!chords->R)
690 return AVERROR(ENOMEM);
691 r_cap *= 2;
692 }
693 chords->R[chords->Lnum++] = val;
694 }
695 }
696
697 /*
698 * Setting the length indices of chords.
699 * These are needed so that the algorithm can, for each chord,
700 * access the lookup table at the correct length in constant time.
701 */
702 chord_length_index = 0;
703 for (int i = 0; i < chords->size; i++) {
704 while (chords->R[chord_length_index] < chords->C[i].l)
705 chord_length_index++;
706 chords->C[i].i = chord_length_index;
707 }
708
709 /*
710 * Chords are sorted on Y. This way, when a row of the lookup table or IPlane
711 * is cached, the next chord offset has a better chance of being on the
712 * same cache line.
713 */
714 qsort(chords->C, chords->size, sizeof(chord), comp_chord);
715
716 return 0;
717 }
718
free_iplane(IPlane * imp)719 static void free_iplane(IPlane *imp)
720 {
721 av_freep(&imp->img);
722 }
723
read_iplane(IPlane * imp,const uint8_t * dst,int dst_linesize,int w,int h,int R,int type_size,int depth)724 static int read_iplane(IPlane *imp, const uint8_t *dst, int dst_linesize,
725 int w, int h, int R, int type_size, int depth)
726 {
727 if (!imp->img)
728 imp->img = av_calloc(h, sizeof(*imp->img));
729 if (!imp->img)
730 return AVERROR(ENOMEM);
731
732 imp->w = w;
733 imp->h = h;
734 imp->range = R;
735 imp->depth = depth;
736 imp->type_size = type_size;
737 imp->max_out_place = type_size == 1 ? max_fun : max16_fun;
738 imp->min_out_place = type_size == 1 ? min_fun : min16_fun;
739 imp->diff_rin_place = type_size == 1 ? diff_fun : diff16_fun;
740 imp->max_in_place = type_size == 1 ? maxinplace_fun : maxinplace16_fun;
741 imp->min_in_place = type_size == 1 ? mininplace_fun : mininplace16_fun;
742 imp->diff_in_place = type_size == 1 ? diffinplace_fun : diffinplace16_fun;
743
744 for (int y = 0; y < h; y++)
745 imp->img[y] = (uint8_t *)dst + y * dst_linesize;
746
747 return 0;
748 }
749
config_input(AVFilterLink * inlink)750 static int config_input(AVFilterLink *inlink)
751 {
752 const AVPixFmtDescriptor *desc = av_pix_fmt_desc_get(inlink->format);
753 MorphoContext *s = inlink->dst->priv;
754
755 s->depth = desc->comp[0].depth;
756 s->type_size = (s->depth + 7) / 8;
757 s->nb_planes = desc->nb_components;
758 s->planewidth[1] = s->planewidth[2] = AV_CEIL_RSHIFT(inlink->w, desc->log2_chroma_w);
759 s->planewidth[0] = s->planewidth[3] = inlink->w;
760 s->planeheight[1] = s->planeheight[2] = AV_CEIL_RSHIFT(inlink->h, desc->log2_chroma_h);
761 s->planeheight[0] = s->planeheight[3] = inlink->h;
762
763 return 0;
764 }
765
config_input_structure(AVFilterLink * inlink)766 static int config_input_structure(AVFilterLink *inlink)
767 {
768 const AVPixFmtDescriptor *desc = av_pix_fmt_desc_get(inlink->format);
769 AVFilterContext *ctx = inlink->dst;
770 MorphoContext *s = inlink->dst->priv;
771
772 av_assert0(ctx->inputs[0]->format == ctx->inputs[1]->format);
773
774 s->splanewidth[1] = s->splanewidth[2] = AV_CEIL_RSHIFT(inlink->w, desc->log2_chroma_w);
775 s->splanewidth[0] = s->splanewidth[3] = inlink->w;
776 s->splaneheight[1] = s->splaneheight[2] = AV_CEIL_RSHIFT(inlink->h, desc->log2_chroma_h);
777 s->splaneheight[0] = s->splaneheight[3] = inlink->h;
778
779 return 0;
780 }
781
activate(AVFilterContext * ctx)782 static int activate(AVFilterContext *ctx)
783 {
784 MorphoContext *s = ctx->priv;
785 return ff_framesync_activate(&s->fs);
786 }
787
do_morpho(FFFrameSync * fs)788 static int do_morpho(FFFrameSync *fs)
789 {
790 AVFilterContext *ctx = fs->parent;
791 AVFilterLink *outlink = ctx->outputs[0];
792 MorphoContext *s = ctx->priv;
793 AVFrame *in = NULL, *structurepic = NULL;
794 AVFrame *out;
795 int ret;
796
797 ret = ff_framesync_dualinput_get(fs, &in, &structurepic);
798 if (ret < 0)
799 return ret;
800 if (!structurepic)
801 return ff_filter_frame(outlink, in);
802
803 out = ff_get_video_buffer(outlink, outlink->w, outlink->h);
804 if (!out) {
805 av_frame_free(&in);
806 return AVERROR(ENOMEM);
807 }
808 av_frame_copy_props(out, in);
809
810 for (int p = 0; p < s->nb_planes; p++) {
811 const uint8_t *src = in->data[p];
812 int src_linesize = in->linesize[p];
813 const uint8_t *ssrc = structurepic->data[p];
814 const int ssrc_linesize = structurepic->linesize[p];
815 uint8_t *dst = out->data[p];
816 int dst_linesize = out->linesize[p];
817 const int swidth = s->splanewidth[p];
818 const int sheight = s->splaneheight[p];
819 const int width = s->planewidth[p];
820 const int height = s->planeheight[p];
821 const int depth = s->depth;
822 int type_size = s->type_size;
823
824 if (ctx->is_disabled || !(s->planes & (1 << p))) {
825 copy:
826 av_image_copy_plane(out->data[p] + 0 * out->linesize[p],
827 out->linesize[p],
828 in->data[p] + 0 * in->linesize[p],
829 in->linesize[p],
830 width * ((s->depth + 7) / 8),
831 height);
832 continue;
833 }
834
835 if (!s->got_structure[p] || s->structures) {
836 free_chord_set(&s->SE[p]);
837
838 ret = read_iplane(&s->SEimg[p], ssrc, ssrc_linesize, swidth, sheight, 1, type_size, depth);
839 if (ret < 0)
840 goto fail;
841 ret = build_chord_set(&s->SEimg[p], &s->SE[p]);
842 if (ret < 0)
843 goto fail;
844 s->got_structure[p] = 1;
845 }
846
847 if (s->SE[p].minX == INT16_MAX ||
848 s->SE[p].minY == INT16_MAX ||
849 s->SE[p].maxX == INT16_MIN ||
850 s->SE[p].maxY == INT16_MIN)
851 goto copy;
852
853 ret = read_iplane(&s->f[p], src, src_linesize, width, height, 1, type_size, depth);
854 if (ret < 0)
855 goto fail;
856
857 ret = read_iplane(&s->g[p], dst, dst_linesize, s->f[p].w, s->f[p].h, s->f[p].range, type_size, depth);
858 if (ret < 0)
859 goto fail;
860
861 switch (s->mode) {
862 case ERODE:
863 ret = erode(&s->g[p], &s->f[p], &s->SE[p], &s->Ty[0][p]);
864 break;
865 case DILATE:
866 ret = dilate(&s->g[p], &s->f[p], &s->SE[p], &s->Ty[0][p]);
867 break;
868 case OPEN:
869 ret = read_iplane(&s->h[p], s->temp->data[p], s->temp->linesize[p], width, height, 1, type_size, depth);
870 if (ret < 0)
871 break;
872 ret = erode(&s->h[p], &s->f[p], &s->SE[p], &s->Ty[0][p]);
873 if (ret < 0)
874 break;
875 ret = dilate(&s->g[p], &s->h[p], &s->SE[p], &s->Ty[1][p]);
876 break;
877 case CLOSE:
878 ret = read_iplane(&s->h[p], s->temp->data[p], s->temp->linesize[p], width, height, 1, type_size, depth);
879 if (ret < 0)
880 break;
881 ret = dilate(&s->h[p], &s->f[p], &s->SE[p], &s->Ty[0][p]);
882 if (ret < 0)
883 break;
884 ret = erode(&s->g[p], &s->h[p], &s->SE[p], &s->Ty[1][p]);
885 break;
886 case GRADIENT:
887 ret = read_iplane(&s->h[p], s->temp->data[p], s->temp->linesize[p], width, height, 1, type_size, depth);
888 if (ret < 0)
889 break;
890 ret = dilate(&s->g[p], &s->f[p], &s->SE[p], &s->Ty[0][p]);
891 if (ret < 0)
892 break;
893 ret = erode(&s->h[p], &s->f[p], &s->SE[p], &s->Ty[1][p]);
894 if (ret < 0)
895 break;
896 difference(&s->g[p], &s->h[p]);
897 break;
898 case TOPHAT:
899 ret = read_iplane(&s->h[p], s->temp->data[p], s->temp->linesize[p], width, height, 1, type_size, depth);
900 if (ret < 0)
901 break;
902 ret = erode(&s->h[p], &s->f[p], &s->SE[p], &s->Ty[0][p]);
903 if (ret < 0)
904 break;
905 ret = dilate(&s->g[p], &s->h[p], &s->SE[p], &s->Ty[1][p]);
906 if (ret < 0)
907 break;
908 difference2(&s->g[p], &s->f[p]);
909 break;
910 case BLACKHAT:
911 ret = read_iplane(&s->h[p], s->temp->data[p], s->temp->linesize[p], width, height, 1, type_size, depth);
912 if (ret < 0)
913 break;
914 ret = dilate(&s->h[p], &s->f[p], &s->SE[p], &s->Ty[0][p]);
915 if (ret < 0)
916 break;
917 ret = erode(&s->g[p], &s->h[p], &s->SE[p], &s->Ty[1][p]);
918 if (ret < 0)
919 break;
920 difference(&s->g[p], &s->f[p]);
921 break;
922 default:
923 av_assert0(0);
924 }
925
926 if (ret < 0)
927 goto fail;
928 }
929
930 av_frame_free(&in);
931 out->pts = av_rescale_q(s->fs.pts, s->fs.time_base, outlink->time_base);
932 return ff_filter_frame(outlink, out);
933 fail:
934 av_frame_free(&out);
935 av_frame_free(&in);
936 return ret;
937 }
938
config_output(AVFilterLink * outlink)939 static int config_output(AVFilterLink *outlink)
940 {
941 AVFilterContext *ctx = outlink->src;
942 MorphoContext *s = ctx->priv;
943 AVFilterLink *mainlink = ctx->inputs[0];
944 int ret;
945
946 s->fs.on_event = do_morpho;
947 ret = ff_framesync_init_dualinput(&s->fs, ctx);
948 if (ret < 0)
949 return ret;
950 outlink->w = mainlink->w;
951 outlink->h = mainlink->h;
952 outlink->time_base = mainlink->time_base;
953 outlink->sample_aspect_ratio = mainlink->sample_aspect_ratio;
954 outlink->frame_rate = mainlink->frame_rate;
955
956 if ((ret = ff_framesync_configure(&s->fs)) < 0)
957 return ret;
958 outlink->time_base = s->fs.time_base;
959
960 s->temp = ff_get_video_buffer(outlink, outlink->w, outlink->h);
961 if (!s->temp)
962 return AVERROR(ENOMEM);
963
964 s->plane_f = av_calloc(outlink->w * outlink->h, sizeof(*s->plane_f));
965 s->plane_g = av_calloc(outlink->w * outlink->h, sizeof(*s->plane_g));
966 if (!s->plane_f || !s->plane_g)
967 return AVERROR(ENOMEM);
968
969 return 0;
970 }
971
uninit(AVFilterContext * ctx)972 static av_cold void uninit(AVFilterContext *ctx)
973 {
974 MorphoContext *s = ctx->priv;
975
976 for (int p = 0; p < 4; p++) {
977 free_iplane(&s->SEimg[p]);
978 free_iplane(&s->f[p]);
979 free_iplane(&s->g[p]);
980 free_iplane(&s->h[p]);
981 free_chord_set(&s->SE[p]);
982 free_lut(&s->Ty[0][p]);
983 free_lut(&s->Ty[1][p]);
984 }
985
986 ff_framesync_uninit(&s->fs);
987
988 av_frame_free(&s->temp);
989 av_freep(&s->plane_f);
990 av_freep(&s->plane_g);
991 }
992
993 static const AVFilterPad morpho_inputs[] = {
994 {
995 .name = "default",
996 .type = AVMEDIA_TYPE_VIDEO,
997 .config_props = config_input,
998 },
999 {
1000 .name = "structure",
1001 .type = AVMEDIA_TYPE_VIDEO,
1002 .config_props = config_input_structure,
1003 },
1004 };
1005
1006 static const AVFilterPad morpho_outputs[] = {
1007 {
1008 .name = "default",
1009 .type = AVMEDIA_TYPE_VIDEO,
1010 .config_props = config_output,
1011 },
1012 };
1013
1014 const AVFilter ff_vf_morpho = {
1015 .name = "morpho",
1016 .description = NULL_IF_CONFIG_SMALL("Apply Morphological filter."),
1017 .preinit = morpho_framesync_preinit,
1018 .priv_size = sizeof(MorphoContext),
1019 .priv_class = &morpho_class,
1020 .activate = activate,
1021 .uninit = uninit,
1022 FILTER_INPUTS(morpho_inputs),
1023 FILTER_OUTPUTS(morpho_outputs),
1024 FILTER_PIXFMTS_ARRAY(pix_fmts),
1025 .flags = AVFILTER_FLAG_SUPPORT_TIMELINE_INTERNAL,
1026 .process_command = ff_filter_process_command,
1027 };
1028