6 #define K2 7 /* from -.7 to +.7 inclusive, meaning .2 into each adjacent pixel */
22 * i0(x) is the modified Bessel function, Σ (x/2)^2L / (L!)²
23 * There are faster ways to calculate this, but we precompute
24 * into a table so let's keep it simple.
34 v += pow(x/2., 2*L)/pow(fac(L), 2);
39 kaiser(double x, double τ, double α)
43 return i0(α*sqrt(1-(x*x/(τ*τ))))/i0(α);
49 fprint(2, "usage: resample [-x xsize] [-y ysize] [imagefile]\n");
50 fprint(2, "\twhere size is an integer or a percentage in the form 25%%\n");
55 getint(char *s, int *percent)
59 *percent = (s[strlen(s)-1] == '%');
68 resamplex(uchar *in, int off, int d, int inx, uchar *out, int outx)
74 rat = (double)inx/(double)outx;
75 for(x=0; x<outx; x++){
77 /* don't resample if size unchanged */
78 out[off+x*d] = in[off+x*d];
83 for(k=-K2; k<=K2; k++){
90 v += in[off+i*d] * K[K2+k];
97 resampley(uchar **in, int off, int iny, uchar **out, int outy)
100 double Y, yy, v, rat;
102 rat = (double)iny/(double)outy;
103 for(y=0; y<outy; y++){
105 /* don't resample if size unchanged */
106 out[y][off] = in[y][off];
111 for(k=-K2; k<=K2; k++){
118 v += in[i][off] * K[K2+k];
134 resample(int xsize, int ysize, Memimage *m)
136 int i, j, d, bpl, nchan;
138 uchar **oscan, **nscan;
140 new = allocmemimage(Rect(0, 0, xsize, ysize), m->chan);
142 sysfatal("can't allocate new image: %r");
144 oscan = malloc(Dy(m->r)*sizeof(uchar*));
145 nscan = malloc(max(ysize, Dy(m->r))*sizeof(uchar*));
146 if(oscan == nil || nscan == nil)
147 sysfatal("can't allocate: %r");
149 /* unload original image into scan lines */
150 bpl = bytesperline(m->r, m->depth);
151 for(i=0; i<Dy(m->r); i++){
152 oscan[i] = malloc(bpl);
154 sysfatal("can't allocate: %r");
155 j = unloadmemimage(m, Rect(m->r.min.x, m->r.min.y+i, m->r.max.x, m->r.min.y+i+1), oscan[i], bpl);
157 sysfatal("unloadmemimage");
160 /* allocate scan lines for destination. we do y first, so need at least Dy(m->r) lines */
161 bpl = bytesperline(Rect(0, 0, xsize, Dy(m->r)), m->depth);
162 for(i=0; i<max(ysize, Dy(m->r)); i++){
163 nscan[i] = malloc(bpl);
165 sysfatal("can't allocate: %r");
169 nchan = d = m->depth/8;
170 if(m->chan == XRGB32)
172 for(i=0; i<Dy(m->r); i++){
173 for(j=0; j<nchan; j++)
174 resamplex(oscan[i], j, d, Dx(m->r), nscan[i], xsize);
177 nscan[i] = malloc(bpl);
179 sysfatal("can't allocate: %r");
183 for(i=0; i<xsize; i++)
184 for(j=0; j<nchan; j++)
185 resampley(oscan, d*i+j, Dy(m->r), nscan, ysize);
187 /* pack data into destination */
188 bpl = bytesperline(new->r, m->depth);
189 for(i=0; i<ysize; i++){
190 j = loadmemimage(new, Rect(0, i, xsize, i+1), nscan[i], bpl);
192 sysfatal("loadmemimage: %r");
198 main(int argc, char *argv[])
200 int i, fd, xsize, ysize, xpercent, ypercent;
202 Memimage *m, *new, *t1, *t2;
207 for(i=-K2; i<=K2; i++){
208 K[K2+i] = kaiser(i/10., K2/10., 4.);
209 // print("%g %g\n", i/10., K[K2+i]);
220 memset(&rparam, 0, sizeof rparam);
222 xpercent = ypercent = 0;
225 case 'a': /* compatibility; equivalent to just -x or -y */
226 if(xsize != 0 || ysize != 0)
228 xsize = getint(ARGF(), &xpercent);
237 xsize = getint(ARGF(), &xpercent);
244 ysize = getint(ARGF(), &ypercent);
252 if(xsize == 0 && ysize == 0)
261 fd = open(file, OREAD);
263 sysfatal("can't open %s: %r", file);
266 m = readmemimage(fd);
268 sysfatal("can't read %s: %r", file);
271 xsize = Dx(m->r)*xsize/100;
273 ysize = Dy(m->r)*ysize/100;
275 ysize = (xsize * Dy(m->r)) / Dx(m->r);
277 xsize = (ysize * Dx(m->r)) / Dy(m->r);
281 for(tchan = m->chan; tchan; tchan >>= 8)
282 if(TYPE(tchan) == CAlpha){
294 new = resample(xsize, ysize, m);
302 /* use library to convert to byte-per-chan form, then convert back */
303 t1 = allocmemimage(m->r, tchan);
305 sysfatal("can't allocate temporary image: %r");
306 memimagedraw(t1, t1->r, m, m->r.min, nil, ZP, S);
307 t2 = resample(xsize, ysize, t1);
309 new = allocmemimage(Rect(0, 0, xsize, ysize), m->chan);
311 sysfatal("can't allocate new image: %r");
312 /* should do error diffusion here */
313 memimagedraw(new, new->r, t2, t2->r.min, nil, ZP, S);
319 if(writememimage(1, new) < 0)
320 sysfatal("write error on output: %r");