Rename GP_Context -> GP_Pixmap
[gfxprim.git] / libs / filters / GP_Median.c
bloba0327dc855c716e9a3c20d4f06916a9e2654f5f8
1 /*****************************************************************************
2 * This file is part of gfxprim library. *
3 * *
4 * Gfxprim is free software; you can redistribute it and/or *
5 * modify it under the terms of the GNU Lesser General Public *
6 * License as published by the Free Software Foundation; either *
7 * version 2.1 of the License, or (at your option) any later version. *
8 * *
9 * Gfxprim is distributed in the hope that it will be useful, *
10 * but WITHOUT ANY WARRANTY; without even the implied warranty of *
11 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU *
12 * Lesser General Public License for more details. *
13 * *
14 * You should have received a copy of the GNU Lesser General Public *
15 * License along with gfxprim; if not, write to the Free Software *
16 * Foundation, Inc., 51 Franklin Street, Fifth Floor, *
17 * Boston, MA 02110-1301 USA *
18 * *
19 * Copyright (C) 2009-2013 Cyril Hrubis <metan@ucw.cz> *
20 * *
21 *****************************************************************************/
23 #include <errno.h>
25 #include "core/GP_Pixmap.h"
26 #include "core/GP_GetPutPixel.h"
27 #include "core/GP_TempAlloc.h"
28 #include "core/GP_Clamp.h"
30 #include "core/GP_Debug.h"
32 #include "GP_Median.h"
34 #include <string.h>
36 struct hist8 {
37 unsigned int coarse[16];
38 unsigned int fine[16][16];
41 struct hist8u {
42 unsigned int coarse[16];
43 unsigned int fine[16][16];
44 unsigned int lx[16];
47 static inline void hist8_inc(struct hist8 *h, unsigned int x, unsigned int val)
49 h[x].coarse[val>>4]++;
50 h[x].fine[val>>4][val&0x0f]++;
53 static inline void hist8_dec(struct hist8 *h, unsigned int x, unsigned int val)
55 h[x].coarse[val>>4]--;
56 h[x].fine[val>>4][val&0x0f]--;
59 static inline void hist8_add(struct hist8u *out, struct hist8 *in, unsigned int x)
61 int i;
63 for (i = 0; i < 16; i += 4) {
64 out->coarse[i + 0] += in[x].coarse[i + 0];
65 out->coarse[i + 1] += in[x].coarse[i + 1];
66 out->coarse[i + 2] += in[x].coarse[i + 2];
67 out->coarse[i + 3] += in[x].coarse[i + 3];
71 static inline void hist8_add_fine(struct hist8u *out, struct hist8 *in, unsigned int x)
73 int i, j;
75 for (i = 0; i < 16; i += 4) {
76 out->coarse[i + 0] += in[x].coarse[i + 0];
77 out->coarse[i + 1] += in[x].coarse[i + 1];
78 out->coarse[i + 2] += in[x].coarse[i + 2];
79 out->coarse[i + 3] += in[x].coarse[i + 3];
81 for (j = 0; j < 16; j++) {
82 out->fine[i + 0][j] += in[x].fine[i + 0][j];
83 out->fine[i + 1][j] += in[x].fine[i + 1][j];
84 out->fine[i + 2][j] += in[x].fine[i + 2][j];
85 out->fine[i + 3][j] += in[x].fine[i + 3][j];
90 static inline void hist8_sub(struct hist8u *out, struct hist8 *in, unsigned int x)
92 int i;
94 for (i = 0; i < 16; i += 4) {
95 out->coarse[i + 0] -= in[x].coarse[i + 0];
96 out->coarse[i + 1] -= in[x].coarse[i + 1];
97 out->coarse[i + 2] -= in[x].coarse[i + 2];
98 out->coarse[i + 3] -= in[x].coarse[i + 3];
103 * Updates only one specified fine part of the histogram.
105 * The structure hist8u remebers when was particular fine row updated so we either
106 * generate it from scatch or update depending on the number of needed operations.
108 static inline void hist8_update(struct hist8u *h, unsigned int i,
109 struct hist8 *row, unsigned int x, unsigned int xmed)
111 unsigned int j, k;
112 unsigned int lx = h->lx[i];
113 unsigned int dx = x - lx;
115 if (dx > 2*xmed) {
116 /* if last update was long ago clear it and load again */
117 for (j = 0; j < 16; j+=4) {
118 h->fine[i][j + 0] = 0;
119 h->fine[i][j + 1] = 0;
120 h->fine[i][j + 2] = 0;
121 h->fine[i][j + 3] = 0;
124 for (j = 0; j < 16; j++)
125 for (k = 0; k <= 2*xmed; k++)
126 h->fine[i][j + 0] += row[x + k].fine[i][j + 0];
127 } else {
128 /* update only missing bits */
129 for (j = 0; j < 16; j++) {
130 for (k = 0; k < dx; k++) {
131 h->fine[i][j] -= row[lx + k].fine[i][j];
132 h->fine[i][j] += row[lx + k + 2*xmed + 1].fine[i][j];
138 h->lx[i] = x;
142 static inline unsigned int hist8_median(struct hist8u *h, struct hist8 *row,
143 unsigned int x, int xmed, unsigned int trigger)
145 unsigned int i, j;
146 unsigned int acc = 0;
148 for (i = 0; i < 16; i++) {
149 acc += h->coarse[i];
151 if (acc >= trigger) {
152 acc -= h->coarse[i];
154 /* update fine on position i */
155 hist8_update(h, i, row, x, xmed);
157 for (j = 0; j < 16; j++) {
158 acc += h->fine[i][j];
160 if (acc >= trigger)
161 return (i<<4) | j;
166 GP_BUG("Trigger not reached");
167 return 0;
170 static int GP_FilterMedian_Raw(const GP_Pixmap *src,
171 GP_Coord x_src, GP_Coord y_src,
172 GP_Size w_src, GP_Size h_src,
173 GP_Pixmap *dst,
174 GP_Coord x_dst, GP_Coord y_dst,
175 int xmed, int ymed,
176 GP_ProgressCallback *callback)
178 int i, x, y;
179 unsigned int trigger = ((2*xmed+1)*(2*ymed+1))/2;
181 if (src->pixel_type != GP_PIXEL_RGB888) {
182 errno = ENOSYS;
183 return -1;
186 GP_DEBUG(1, "Median filter size %ux%u xmed=%u ymed=%u",
187 w_src, h_src, 2 * xmed + 1, 2 * ymed + 1);
189 /* The buffer is w + 2*xmed + 1 size because we read the last value but we don't use it */
190 unsigned int size = (w_src + 2 * xmed + 1);
192 /* Create and initalize arrays for row of histograms */
193 GP_TempAllocCreate(temp, 3 * sizeof(struct hist8) * size + 3 * sizeof(struct hist8u));
195 struct hist8 *R = GP_TempAllocGet(temp, sizeof(struct hist8) * size);
196 struct hist8 *G = GP_TempAllocGet(temp, sizeof(struct hist8) * size);
197 struct hist8 *B = GP_TempAllocGet(temp, sizeof(struct hist8) * size);
199 memset(R, 0, sizeof(*R) * size);
200 memset(G, 0, sizeof(*G) * size);
201 memset(B, 0, sizeof(*B) * size);
203 struct hist8u *XR = GP_TempAllocGet(temp, sizeof(struct hist8u));
204 struct hist8u *XG = GP_TempAllocGet(temp, sizeof(struct hist8u));
205 struct hist8u *XB = GP_TempAllocGet(temp, sizeof(struct hist8u));
207 /* Prefill row of histograms */
208 for (x = 0; x < (int)w_src + 2*xmed; x++) {
209 int xi = GP_CLAMP(x_src + x - xmed, 0, (int)src->w - 1);
211 for (y = -ymed; y <= ymed; y++) {
212 int yi = GP_CLAMP(y_src + y, 0, (int)src->h - 1);
214 GP_Pixel pix = GP_GetPixel_Raw_24BPP(src, xi, yi);
216 hist8_inc(R, x, GP_Pixel_GET_R_RGB888(pix));
217 hist8_inc(G, x, GP_Pixel_GET_G_RGB888(pix));
218 hist8_inc(B, x, GP_Pixel_GET_B_RGB888(pix));
222 /* Apply the median filter */
223 for (y = 0; y < (int)h_src; y++) {
224 memset(XR, 0, sizeof(*XR));
225 memset(XG, 0, sizeof(*XG));
226 memset(XB, 0, sizeof(*XB));
228 /* Compute first histogram */
229 for (i = 0; i <= 2*xmed; i++) {
230 hist8_add_fine(XR, R, i);
231 hist8_add_fine(XG, G, i);
232 hist8_add_fine(XB, B, i);
235 /* Generate row */
236 for (x = 0; x < (int)w_src; x++) {
237 int r = hist8_median(XR, R, x, xmed, trigger);
238 int g = hist8_median(XG, G, x, xmed, trigger);
239 int b = hist8_median(XB, B, x, xmed, trigger);
241 GP_PutPixel_Raw_24BPP(dst, x_dst + x, y_dst + y,
242 GP_Pixel_CREATE_RGB888(r, g, b));
244 /* Recompute histograms */
245 hist8_sub(XR, R, x);
246 hist8_sub(XG, G, x);
247 hist8_sub(XB, B, x);
249 hist8_add(XR, R, (x + 2 * xmed + 1));
250 hist8_add(XG, G, (x + 2 * xmed + 1));
251 hist8_add(XB, B, (x + 2 * xmed + 1));
254 /* Recompute histograms, remove y - ymed pixel add y + ymed + 1 */
255 for (x = 0; x < (int)w_src + 2*xmed; x++) {
256 int xi = GP_CLAMP(x_src + x - xmed, 0, (int)src->w - 1);
257 int yi = GP_CLAMP(y_src + y - ymed, 0, (int)src->h - 1);
259 GP_Pixel pix = GP_GetPixel_Raw_24BPP(src, xi, yi);
261 hist8_dec(R, x, GP_Pixel_GET_R_RGB888(pix));
262 hist8_dec(G, x, GP_Pixel_GET_G_RGB888(pix));
263 hist8_dec(B, x, GP_Pixel_GET_B_RGB888(pix));
265 yi = GP_MIN(y_src + y + ymed + 1, (int)src->h - 1);
267 pix = GP_GetPixel_Raw_24BPP(src, xi, yi);
269 hist8_inc(R, x, GP_Pixel_GET_R_RGB888(pix));
270 hist8_inc(G, x, GP_Pixel_GET_G_RGB888(pix));
271 hist8_inc(B, x, GP_Pixel_GET_B_RGB888(pix));
274 if (GP_ProgressCallbackReport(callback, y, h_src, w_src)) {
275 GP_TempAllocFree(temp);
276 return 1;
280 GP_TempAllocFree(temp);
281 GP_ProgressCallbackDone(callback);
283 return 0;
286 int GP_FilterMedianEx(const GP_Pixmap *src,
287 GP_Coord x_src, GP_Coord y_src,
288 GP_Size w_src, GP_Size h_src,
289 GP_Pixmap *dst,
290 GP_Coord x_dst, GP_Coord y_dst,
291 int xmed, int ymed,
292 GP_ProgressCallback *callback)
294 GP_CHECK(src->pixel_type == dst->pixel_type);
296 /* Check that destination is large enough */
297 GP_CHECK(x_dst + (GP_Coord)w_src <= (GP_Coord)dst->w);
298 GP_CHECK(y_dst + (GP_Coord)h_src <= (GP_Coord)dst->h);
300 GP_CHECK(xmed >= 0 && ymed >= 0);
302 return GP_FilterMedian_Raw(src, x_src, y_src, w_src, h_src,
303 dst, x_dst, y_dst, xmed, ymed, callback);
306 GP_Pixmap *GP_FilterMedianExAlloc(const GP_Pixmap *src,
307 GP_Coord x_src, GP_Coord y_src,
308 GP_Size w_src, GP_Size h_src,
309 int xmed, int ymed,
310 GP_ProgressCallback *callback)
312 int ret;
314 GP_CHECK(xmed >= 0 && ymed >= 0);
316 GP_Pixmap *dst = GP_PixmapAlloc(w_src, h_src, src->pixel_type);
318 if (dst == NULL)
319 return NULL;
321 ret = GP_FilterMedian_Raw(src, x_src, y_src, w_src, h_src,
322 dst, 0, 0, xmed, ymed, callback);
324 if (ret) {
325 GP_PixmapFree(dst);
326 return NULL;
329 return dst;