1 /* This file is an image processing operation for GEGL 2 * 3 * GEGL is free software; you can redistribute it and/or 4 * modify it under the terms of the GNU Lesser General Public 5 * License as published by the Free Software Foundation; either 6 * version 3 of the License, or (at your option) any later version. 7 * 8 * GEGL is distributed in the hope that it will be useful, 9 * but WITHOUT ANY WARRANTY; without even the implied warranty of 10 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU 11 * Lesser General Public License for more details. 12 * 13 * You should have received a copy of the GNU Lesser General Public 14 * License along with GEGL; if not, see <https://www.gnu.org/licenses/>. 15 * 16 * Copyright (C) 1997 Lauri Alanko <la@iki.fi> 17 * Copyright 2011 Robert Sasu (sasu.robert@gmail.com) 18 */ 19 20 #include "config.h" 21 #include <glib/gi18n-lib.h> 22 23 #ifdef GEGL_PROPERTIES 24 25 property_double (a1, _("(1,1)"), 0.0) 26 property_double (a2, _("(1,2)"), 0.0) 27 property_double (a3, _("(1,3)"), 0.0) 28 property_double (a4, _("(1,4)"), 0.0) 29 property_double (a5, _("(1,5)"), 0.0) 30 property_double (b1, _("(2,1)"), 0.0) 31 property_double (b2, _("(2,2)"), 0.0) 32 property_double (b3, _("(2,3)"), 0.0) 33 property_double (b4, _("(2,4)"), 0.0) 34 property_double (b5, _("(2,5)"), 0.0) 35 property_double (c1, _("(3,1)"), 0.0) 36 property_double (c2, _("(3,2)"), 0.0) 37 property_double (c3, _("(3,3)"), 1.0) 38 property_double (c4, _("(3,4)"), 0.0) 39 property_double (c5, _("(3,5)"), 0.0) 40 property_double (d1, _("(4,1)"), 0.0) 41 property_double (d2, _("(4,2)"), 0.0) 42 property_double (d3, _("(4,3)"), 0.0) 43 property_double (d4, _("(4,4)"), 0.0) 44 property_double (d5, _("(4,5)"), 0.0) 45 property_double (e1, _("(5,1)"), 0.0) 46 property_double (e2, _("(5,2)"), 0.0) 47 property_double (e3, _("(5,3)"), 0.0) 48 property_double (e4, _("(5,4)"), 0.0) 49 property_double (e5, _("(5,5)"), 0.0) 50 51 property_double (divisor, _("Divisor"), 1.0) 52 ui_range (-1000.0, 1000.0) 53 ui_meta ("sensitive", "! normalize") 54 55 property_double (offset, _("Offset"), 0.0) 56 value_range (-1.0, 1.0) 57 ui_meta ("sensitive", "! normalize") 58 59 property_boolean (red, _("Red channel"), TRUE) 60 property_boolean (green, _("Green channel"), TRUE) 61 property_boolean (blue, _("Blue channel"), TRUE) 62 property_boolean (alpha, _("Alpha channel"), TRUE) 63 64 property_boolean (normalize, _("Normalize"), TRUE) 65 property_boolean (alpha_weight, _("Alpha-weighting"), TRUE) 66 67 property_enum (border, _("Border"), 68 GeglAbyssPolicy, gegl_abyss_policy, 69 GEGL_ABYSS_CLAMP) 70 71 #else 72 73 #define GEGL_OP_AREA_FILTER 74 #define GEGL_OP_NAME convolution_matrix 75 #define GEGL_OP_C_SOURCE convolution-matrix.c 76 77 #include "gegl-op.h" 78 #include <stdio.h> 79 80 #define MAX_MATRIX_SIZE 5 81 82 static gboolean 83 enough_with_3x3 (GeglProperties *o) 84 { 85 if (o->a1 == 0.0 && 86 o->a2 == 0.0 && 87 o->a3 == 0.0 && 88 o->a4 == 0.0 && 89 o->a5 == 0.0 && 90 o->b1 == 0.0 && 91 o->b5 == 0.0 && 92 o->c1 == 0.0 && 93 o->c5 == 0.0 && 94 o->d1 == 0.0 && 95 o->d5 == 0.0 && 96 o->e1 == 0.0 && 97 o->e2 == 0.0 && 98 o->e3 == 0.0 && 99 o->e4 == 0.0 && 100 o->e5 == 0.0) 101 { 102 return TRUE; 103 } 104 return FALSE; 105 } 106 107 static void 108 prepare (GeglOperation *operation) 109 { 110 const Babl *space = gegl_operation_get_source_space (operation, "input"); 111 112 GeglOperationAreaFilter *op_area = GEGL_OPERATION_AREA_FILTER (operation); 113 114 GeglProperties *o = GEGL_PROPERTIES (operation); 115 if (enough_with_3x3 (o)) 116 op_area->left = op_area->right = op_area->top = op_area->bottom = 1; /* 3 */ 117 else 118 op_area->left = op_area->right = op_area->top = op_area->bottom = 2; /* 5 */ 119 120 gegl_operation_set_format (operation, "output", 121 babl_format_with_space ("RGBA float", space)); 122 } 123 124 static void 125 make_matrix (GeglProperties *o, 126 gfloat matrix[MAX_MATRIX_SIZE][MAX_MATRIX_SIZE], 127 gint matrix_size) 128 { 129 if (matrix_size == 3) 130 { 131 matrix[0][0] = o->b2; 132 matrix[0][1] = o->b3; 133 matrix[0][2] = o->b4; 134 135 matrix[1][0] = o->c2; 136 matrix[1][1] = o->c3; 137 matrix[1][2] = o->c4; 138 139 matrix[2][0] = o->d2; 140 matrix[2][1] = o->d3; 141 matrix[2][2] = o->d4; 142 } 143 else 144 { 145 matrix[0][0] = o->a1; 146 matrix[0][1] = o->a2; 147 matrix[0][2] = o->a3; 148 matrix[0][3] = o->a4; 149 matrix[0][4] = o->a5; 150 151 matrix[1][0] = o->b1; 152 matrix[1][1] = o->b2; 153 matrix[1][2] = o->b3; 154 matrix[1][3] = o->b4; 155 matrix[1][4] = o->b5; 156 157 matrix[2][0] = o->c1; 158 matrix[2][1] = o->c2; 159 matrix[2][2] = o->c3; 160 matrix[2][3] = o->c4; 161 matrix[2][4] = o->c5; 162 163 matrix[3][0] = o->d1; 164 matrix[3][1] = o->d2; 165 matrix[3][2] = o->d3; 166 matrix[3][3] = o->d4; 167 matrix[3][4] = o->d5; 168 169 matrix[4][0] = o->e1; 170 matrix[4][1] = o->e2; 171 matrix[4][2] = o->e3; 172 matrix[4][3] = o->e4; 173 matrix[4][4] = o->e5; 174 } 175 } 176 177 static gboolean 178 normalize_div_off (gfloat matrix[MAX_MATRIX_SIZE][MAX_MATRIX_SIZE], 179 gint matrix_size, 180 gfloat *divisor, 181 gfloat *offset) 182 { 183 gint x, y; 184 gboolean valid = FALSE; 185 gfloat sum = 0.0; 186 187 for (y = 0; y < matrix_size; y++) 188 for (x = 0; x < matrix_size; x++) 189 { 190 sum += matrix[x][y]; 191 if (matrix[x][y] != 0.0) 192 valid = TRUE; 193 } 194 195 if (sum > 0) 196 { 197 *offset = 0.0; 198 *divisor = sum; 199 } 200 else if (sum < 0) 201 { 202 *offset = 1.0; 203 *divisor = -sum; 204 } 205 else 206 { 207 *offset = 0.5; 208 *divisor = 1; 209 } 210 211 return valid; 212 } 213 214 static gint inline 215 matrix_center_offset (const GeglRectangle *extended, 216 gint matrix_size) 217 { 218 return (extended->width + 1) * (matrix_size / 2); 219 } 220 221 static void inline 222 convolve_pixel_componentwise (GeglProperties *o, 223 gfloat *src_buf, 224 gfloat *dst_buf, 225 const GeglRectangle *result, 226 const GeglRectangle *extended, 227 gfloat matrix[MAX_MATRIX_SIZE][MAX_MATRIX_SIZE], 228 gint matrix_size, 229 gint d_offset, 230 gint ss_offset, 231 gint xx, 232 gint yy, 233 gfloat matrixsum, 234 gfloat inv_divisor, 235 gfloat offset) 236 { 237 gint i; 238 gint s_stride = (extended->width - matrix_size) * 4; 239 for (i = 0; i < 4; i++) 240 { 241 gfloat sum = 0.0; 242 gint s_offset = ss_offset; 243 244 if ((i == 0 && o->red) || 245 (i == 1 && o->green) || 246 (i == 2 && o->blue) || 247 (i == 3 && o->alpha)) 248 { 249 gint x, y; 250 for (y = 0; y < matrix_size; y++) 251 { 252 for (x = 0; x < matrix_size; x++) 253 { 254 sum += matrix[x][y] * src_buf[s_offset + i]; 255 s_offset += 4; 256 } 257 s_offset += s_stride; 258 } 259 sum = sum * inv_divisor; 260 sum += offset; 261 dst_buf[d_offset + i] = sum; 262 } 263 else 264 { 265 s_offset += 4 * matrix_center_offset (extended, matrix_size); 266 dst_buf[d_offset + i] = src_buf[s_offset + i]; 267 } 268 } 269 } 270 271 static void inline 272 convolve_pixel (GeglProperties *o, 273 gfloat *src_buf, 274 gfloat *dst_buf, 275 const GeglRectangle *result, 276 const GeglRectangle *extended, 277 gfloat matrix[MAX_MATRIX_SIZE][MAX_MATRIX_SIZE], 278 gint matrix_size, 279 gint d_offset, 280 gint ss_offset, 281 gint xx, 282 gint yy, 283 gfloat matrixsum, 284 gfloat inv_divisor, 285 gfloat offset) 286 { 287 gint i; 288 gint s_stride = (extended->width - matrix_size) * 4; 289 for (i = 0; i < 4; i++) 290 { 291 gfloat sum = 0.0; 292 gint s_offset = ss_offset; 293 gint x, y; 294 for (y = 0; y < matrix_size; y++) 295 { 296 for (x = 0; x < matrix_size; x++) 297 { 298 sum += matrix[x][y] * src_buf[s_offset + i]; 299 s_offset += 4; 300 } 301 s_offset += s_stride; 302 } 303 sum = sum * inv_divisor; 304 sum += offset; 305 dst_buf[d_offset + i] = sum; 306 } 307 } 308 309 static void inline 310 convolve_pixel_alpha_weight (GeglProperties *o, 311 gfloat *src_buf, 312 gfloat *dst_buf, 313 const GeglRectangle *result, 314 const GeglRectangle *extended, 315 gfloat matrix[MAX_MATRIX_SIZE][MAX_MATRIX_SIZE], 316 gint matrix_size, 317 gint d_offset, 318 gint ss_offset, 319 gint xx, 320 gint yy, 321 gfloat matrixsum, 322 gfloat inv_divisor, 323 gfloat offset) 324 { 325 gint i; 326 gint s_stride = (extended->width - matrix_size) * 4; 327 328 for (i = 0; i < 3; i++) 329 { 330 gfloat sum = 0.0; 331 gint s_offset = ss_offset; 332 gint x, y; 333 for (y = 0; y < matrix_size; y++) 334 { 335 for (x = 0; x < matrix_size; x++) 336 { 337 sum += matrix[x][y] * src_buf[s_offset + i] * src_buf[s_offset + 3]; 338 s_offset += 4; 339 } 340 s_offset += s_stride; 341 } 342 sum = sum * inv_divisor + offset; 343 dst_buf[d_offset + i] = sum; 344 } 345 { 346 gfloat sum = 0.0; 347 gfloat alphasum = 0.0; 348 gint s_offset = ss_offset; 349 gint x, y; 350 351 for (y = 0; y < matrix_size; y++) 352 { 353 for (x = 0; x < matrix_size; x++) 354 { 355 float val = matrix[x][y] * src_buf[s_offset + i]; 356 sum += val; 357 alphasum += fabsf (val); 358 s_offset += 4; 359 } 360 s_offset += s_stride; 361 } 362 363 if (alphasum > 0.0) 364 { 365 sum = sum * inv_divisor; 366 sum = sum * matrixsum / alphasum + offset; 367 } 368 else 369 sum = offset; 370 dst_buf[d_offset + i] = sum; 371 } 372 } 373 374 static void inline 375 convolve_pixel_alpha_weight_componentwise (GeglProperties *o, 376 gfloat *src_buf, 377 gfloat *dst_buf, 378 const GeglRectangle *result, 379 const GeglRectangle *extended, 380 gfloat matrix[MAX_MATRIX_SIZE][MAX_MATRIX_SIZE], 381 gint matrix_size, 382 gint d_offset, 383 gint ss_offset, 384 gint xx, 385 gint yy, 386 gfloat matrixsum, 387 gfloat inv_divisor, 388 gfloat offset) 389 { 390 gint i; 391 gint s_stride = (extended->width - matrix_size) * 4; 392 393 for (i = 0; i < 3; i++) 394 { 395 gfloat sum = 0.0; 396 gint s_offset = ss_offset; 397 398 if ((i == 0 && o->red) || 399 (i == 1 && o->green) || 400 (i == 2 && o->blue)) 401 { 402 gint x, y; 403 for (y = 0; y < matrix_size; y++) 404 { 405 for (x = 0; x < matrix_size; x++) 406 { 407 sum += matrix[x][y] * src_buf[s_offset + i] * src_buf[s_offset + 3]; 408 s_offset += 4; 409 } 410 s_offset += s_stride; 411 } 412 sum = sum * inv_divisor + offset; 413 } 414 else 415 { 416 s_offset += 4 * matrix_center_offset (extended, matrix_size); 417 sum = src_buf[s_offset + i]; 418 } 419 dst_buf[d_offset + i] = sum; 420 } 421 { 422 gfloat sum = 0.0; 423 gfloat alphasum = 0.0; 424 gint s_offset = ss_offset; 425 426 if (o->alpha) 427 { 428 gint x, y; 429 430 for (y = 0; y < matrix_size; y++) 431 { 432 for (x = 0; x < matrix_size; x++) 433 { 434 float val = matrix[x][y] * src_buf[s_offset + i]; 435 sum += val; 436 alphasum += fabsf (val); 437 s_offset += 4; 438 } 439 s_offset += s_stride; 440 } 441 442 443 if (alphasum != 0) 444 { 445 sum = sum * inv_divisor; 446 sum = sum * matrixsum / alphasum + offset; 447 } 448 else 449 sum = offset; 450 } 451 else 452 { 453 s_offset += 4 * matrix_center_offset (extended, matrix_size); 454 sum = src_buf[s_offset + i]; 455 } 456 dst_buf[d_offset + i] = sum; 457 } 458 } 459 460 static gboolean 461 process (GeglOperation *operation, 462 GeglBuffer *input, 463 GeglBuffer *output, 464 const GeglRectangle *result, 465 gint level) 466 { 467 GeglProperties *o = GEGL_PROPERTIES (operation); 468 GeglOperationAreaFilter *op_area = GEGL_OPERATION_AREA_FILTER (operation); 469 const Babl *format = gegl_operation_get_format (operation, "output"); 470 GeglRectangle rect; 471 gfloat *src_buf; 472 gfloat *dst_buf; 473 gfloat matrix[MAX_MATRIX_SIZE][MAX_MATRIX_SIZE]={{0,}}; 474 gfloat matrixsum = 0.0; 475 gfloat divisor = o->divisor; 476 gfloat offset = o->offset; 477 gfloat inv_divisor; 478 gint x, y; 479 gint matrix_size = MAX_MATRIX_SIZE; 480 481 if (enough_with_3x3 (o)) 482 matrix_size = 3; 483 484 make_matrix (o, matrix, matrix_size); 485 486 if (o->normalize) 487 normalize_div_off (matrix, matrix_size, &divisor, &offset); 488 inv_divisor = 1.0 / divisor; 489 490 for (x = 0; x < matrix_size; x++) 491 for (y = 0; y < matrix_size; y++) 492 matrixsum += fabsf (matrix[x][y]); 493 494 rect.x = result->x - op_area->left; 495 rect.width = result->width + op_area->left + op_area->right; 496 rect.y = result->y - op_area->top; 497 rect.height = result->height + op_area->top + op_area->bottom; 498 499 src_buf = g_new (gfloat, rect.width * rect.height * 4); 500 dst_buf = g_new (gfloat, result->width * result->height * 4); 501 502 gegl_buffer_get (input, &rect, 1.0, format, src_buf, 503 GEGL_AUTO_ROWSTRIDE, o->border); 504 505 if (o->divisor != 0) 506 { 507 gint ss_offset = (result->y - matrix_size/2 - rect.y) * rect.width * 4 + 508 (result->x - matrix_size/2 - rect.x) * 4; 509 int d_offset = 0; 510 if (o->alpha_weight) 511 { 512 if (o->red == FALSE || o->green == FALSE || o->blue == FALSE || o->alpha == FALSE) 513 { 514 gint x, y; 515 for (y = result->y; y < result->height + result->y; y++) 516 { 517 for (x = result->x; x < result->width + result->x; x++) 518 { 519 convolve_pixel_alpha_weight_componentwise (o, src_buf, dst_buf, result, &rect, 520 matrix, matrix_size, d_offset, ss_offset, x, y, matrixsum, inv_divisor, offset); 521 d_offset += 4; 522 ss_offset += 4; 523 } 524 ss_offset += (rect.width - result->width) * 4; 525 } 526 } 527 else 528 { 529 gint x, y; 530 for (y = result->y; y < result->height + result->y; y++) 531 { 532 for (x = result->x; x < result->width + result->x; x++) 533 { 534 convolve_pixel_alpha_weight (o, src_buf, dst_buf, result, &rect, 535 matrix, matrix_size, d_offset, ss_offset, x, y, matrixsum, inv_divisor, offset); 536 d_offset += 4; 537 ss_offset += 4; 538 } 539 ss_offset += (rect.width - result->width) * 4; 540 } 541 } 542 } 543 else 544 { 545 if (o->red == FALSE || o->green == FALSE || o->blue == FALSE || o->alpha == FALSE) 546 { 547 gint x, y; 548 for (y = result->y; y < result->height + result->y; y++) 549 { 550 for (x = result->x; x < result->width + result->x; x++) 551 { 552 convolve_pixel_componentwise (o, src_buf, dst_buf, result, &rect, 553 matrix, matrix_size, d_offset, ss_offset, x, y, matrixsum, inv_divisor, offset); 554 d_offset += 4; 555 ss_offset += 4; 556 } 557 ss_offset += (rect.width - result->width) * 4; 558 } 559 } 560 else 561 { 562 gint x, y; 563 for (y = result->y; y < result->height + result->y; y++) 564 { 565 for (x = result->x; x < result->width + result->x; x++) 566 { 567 convolve_pixel (o, src_buf, dst_buf, result, &rect, 568 matrix, matrix_size, d_offset, ss_offset, x, y, matrixsum, inv_divisor, offset); 569 d_offset += 4; 570 ss_offset += 4; 571 } 572 ss_offset += (rect.width - result->width) * 4; 573 } 574 } 575 } 576 gegl_buffer_set (output, result, 0, format, 577 dst_buf, GEGL_AUTO_ROWSTRIDE); 578 } 579 else 580 { 581 gegl_buffer_set (output, &rect, 0, format, 582 src_buf, GEGL_AUTO_ROWSTRIDE); 583 } 584 585 g_free (src_buf); 586 g_free (dst_buf); 587 588 return TRUE; 589 } 590 591 static GeglRectangle 592 get_bounding_box (GeglOperation *operation) 593 { 594 GeglRectangle result = { 0, }; 595 GeglRectangle *in_rect; 596 597 in_rect = gegl_operation_source_get_bounding_box (operation, "input"); 598 if (!in_rect) 599 return result; 600 601 return *in_rect; 602 } 603 604 static GeglAbyssPolicy 605 get_abyss_policy (GeglOperation *operation, 606 const gchar *input_pad) 607 { 608 GeglProperties *o = GEGL_PROPERTIES (operation); 609 610 return o->border; 611 } 612 613 static void 614 gegl_op_class_init (GeglOpClass *klass) 615 { 616 GeglOperationClass *operation_class; 617 GeglOperationFilterClass *filter_class; 618 GeglOperationAreaFilterClass *area_class; 619 620 operation_class = GEGL_OPERATION_CLASS (klass); 621 filter_class = GEGL_OPERATION_FILTER_CLASS (klass); 622 area_class = GEGL_OPERATION_AREA_FILTER_CLASS (klass); 623 624 area_class->get_abyss_policy = get_abyss_policy; 625 filter_class->process = process; 626 operation_class->prepare = prepare; 627 operation_class->get_bounding_box = get_bounding_box; 628 629 gegl_operation_class_set_keys (operation_class, 630 "categories", "generic", 631 "name", "gegl:convolution-matrix", 632 "reference-hash", "22d2d47a2da3d3e7cd402ea9fa1a3a25", 633 "reference-hashB", "4eddc0aaa970a59ee8a813627874cdf3", 634 "title", _("Convolution Matrix"), 635 "description", _("Apply a generic 5x5 convolution matrix"), 636 NULL); 637 } 638 639 #endif 640