The OpenD Programming Language

1 /++
2 	Image resizing support for [arsd.color.MemoryImage]. Handles up and down scaling.
3 	See [imageResize] for the main function, all others are lower level if you need
4 	more control.
5 
6 
7 	Note that this focuses more on quality than speed. You can tweak the `filterScale`
8 	argument to speed things up at the expense of quality though (lower number = faster).
9 
10 	I've found:
11 
12 	---
13 	auto size = calculateSizeKeepingAspectRatio(i.width, i.height, maxWidth, maxHeight);
14 	if(size.width != i.width || size.height != i.height) {
15 		i = imageResize(i, size.width, size.height, null, 1.0, 0.6);
16 	}
17 	---
18 
19 	Gives decent results balancing quality and speed. Compiling with ldc or gdc can also
20 	speed up your program.
21 
22 
23 
24 	Authors:
25 		Originally written in C by Rich Geldreich, ported to D by ketmar.
26 	License:
27 		Public Domain / Unlicense - http://unlicense.org/
28 +/
29 module arsd.imageresize;
30 
31 import arsd.color;
32 
33 // ////////////////////////////////////////////////////////////////////////// //
34 // Separable filtering image rescaler v2.21, Rich Geldreich - richgel99@gmail.com
35 //
36 // This is free and unencumbered software released into the public domain.
37 //
38 // Anyone is free to copy, modify, publish, use, compile, sell, or
39 // distribute this software, either in source code form or as a compiled
40 // binary, for any purpose, commercial or non-commercial, and by any
41 // means.
42 //
43 // In jurisdictions that recognize copyright laws, the author or authors
44 // of this software dedicate any and all copyright interest in the
45 // software to the public domain. We make this dedication for the benefit
46 // of the public at large and to the detriment of our heirs and
47 // successors. We intend this dedication to be an overt act of
48 // relinquishment in perpetuity of all present and future rights to this
49 // software under copyright law.
50 //
51 // THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
52 // EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF
53 // MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT.
54 // IN NO EVENT SHALL THE AUTHORS BE LIABLE FOR ANY CLAIM, DAMAGES OR
55 // OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE,
56 // ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR
57 // OTHER DEALINGS IN THE SOFTWARE.
58 //
59 // For more information, please refer to <http://unlicense.org/>
60 //
61 // Feb. 1996: Creation, losely based on a heavily bugfixed version of Schumacher's resampler in Graphics Gems 3.
62 // Oct. 2000: Ported to C++, tweaks.
63 // May 2001: Continous to discrete mapping, box filter tweaks.
64 // March 9, 2002: Kaiser filter grabbed from Jonathan Blow's GD magazine mipmap sample code.
65 // Sept. 8, 2002: Comments cleaned up a bit.
66 // Dec. 31, 2008: v2.2: Bit more cleanup, released as public domain.
67 // June 4, 2012: v2.21: Switched to unlicense.org, integrated GCC fixes supplied by Peter Nagy <petern@crytek.com>, Anteru at anteru.net, and clay@coge.net,
68 // added Codeblocks project (for testing with MinGW and GCC), VS2008 static code analysis pass.
69 // float or double
70 private:
71 
72 @system:
73 
74 //version = iresample_debug;
75 
76 
77 // ////////////////////////////////////////////////////////////////////////// //
78 public enum ImageResizeDefaultFilter = "lanczos4"; /// Default filter for image resampler.
79 public enum ImageResizeMaxDimension = 65536; /// Maximum image width/height for image resampler.
80 
81 
82 // ////////////////////////////////////////////////////////////////////////// //
83 /// Number of known image resizer filters.
84 public @property int imageResizeFilterCount () { pragma(inline, true); return NumFilters; }
85 
86 /// Get filter name. Will return `null` for invalid index.
87 public string imageResizeFilterName (long idx) { pragma(inline, true); return (idx >= 0 && idx < NumFilters ? gFilters.ptr[cast(uint)idx].name : null); }
88 
89 /// Find filter index by name. Will use default filter for invalid names.
90 public int imageResizeFindFilter (const(char)[] name, const(char)[] defaultFilter=ImageResizeDefaultFilter) {
91   int res = resamplerFindFilterInternal(name);
92   if (res >= 0) return res;
93   res = resamplerFindFilterInternal(defaultFilter);
94   if (res >= 0) return res;
95   res = resamplerFindFilterInternal("lanczos4");
96   assert(res >= 0);
97   return res;
98 }
99 
100 /++
101 	Calculates a new size that fits inside the maximums while keeping the original aspect ratio.
102 
103 	History:
104 		Added March 18, 2021 (dub v9.4)
105 +/
106 public Size calculateSizeKeepingAspectRatio(int currentWidth, int currentHeight, int maxWidth, int maxHeight) {
107 	if(currentWidth <= maxWidth && currentHeight <= maxHeight)
108 		return Size(currentWidth, currentHeight);
109 
110 	float shrinkage = 1.0;
111 
112 	if(currentWidth > maxWidth) {
113 		shrinkage = cast(float) maxWidth / currentWidth;
114 	}
115 	if(currentHeight > maxHeight) {
116 		auto shrinkage2 = cast(float) maxHeight / currentHeight;
117 		if(shrinkage2 < shrinkage)
118 			shrinkage = shrinkage2;
119 	}
120 
121 	return Size(cast(int) (currentWidth * shrinkage), cast(int) (currentHeight * shrinkage));
122 }
123 
124 // ////////////////////////////////////////////////////////////////////////// //
125 /// Resize image.
126 public TrueColorImage imageResize(int Components=4) (MemoryImage msrcimg, int dstwdt, int dsthgt, const(char)[] filter=null, float gamma=1.0f, float filterScale=1.0f) {
127   static assert(Components == 1 || Components == 3 || Components == 4, "invalid number of components in color");
128   return imageResize!Components(msrcimg, dstwdt, dsthgt, imageResizeFindFilter(filter), gamma, filterScale);
129 }
130 
131 /// ditto
132 public TrueColorImage imageResize(int Components=4) (MemoryImage msrcimg, int dstwdt, int dsthgt, int filter, float gamma=1.0f, float filterScale=1.0f) {
133   static assert(Components == 1 || Components == 3 || Components == 4, "invalid number of components in color");
134   if (msrcimg is null || msrcimg.width < 1 || msrcimg.height < 1 || msrcimg.width > ImageResizeMaxDimension || msrcimg.height > ImageResizeMaxDimension) {
135     throw new Exception("invalid source image");
136   }
137   if (dstwdt < 1 || dsthgt < 1 || dstwdt > ImageResizeMaxDimension || dsthgt > ImageResizeMaxDimension) throw new Exception("invalid destination image size");
138   auto resimg = new TrueColorImage(dstwdt, dsthgt);
139   scope(failure) .destroy(resimg);
140   if (auto tc = cast(TrueColorImage)msrcimg) {
141     imageResize!Components(
142       delegate (Color[] destrow, int y) { destrow[] = tc.imageData.colors[y*tc.width..(y+1)*tc.width]; },
143       delegate (int y, const(Color)[] row) { resimg.imageData.colors[y*resimg.width..(y+1)*resimg.width] = row[]; },
144       msrcimg.width, msrcimg.height, dstwdt, dsthgt, filter, gamma, filterScale
145     );
146   } else {
147     imageResize!Components(
148       delegate (Color[] destrow, int y) { foreach (immutable x, ref c; destrow) c = msrcimg.getPixel(cast(int)x, y); },
149       delegate (int y, const(Color)[] row) { resimg.imageData.colors[y*resimg.width..(y+1)*resimg.width] = row[]; },
150       msrcimg.width, msrcimg.height, dstwdt, dsthgt, filter, gamma, filterScale
151     );
152   }
153   return resimg;
154 }
155 
156 
157 private {
158   enum Linear2srgbTableSize = 4096;
159   enum InvLinear2srgbTableSize = cast(float)(1.0f/Linear2srgbTableSize);
160   float[256] srgb2linear = void;
161   ubyte[Linear2srgbTableSize] linear2srgb = void;
162   float lastGamma = float.nan;
163 }
164 
165 /// Resize image.
166 /// Partial gamma correction looks better on mips; set to 1.0 to disable gamma correction.
167 /// Filter scale: values < 1.0 cause aliasing, but create sharper looking mips (0.75f, for example).
168 public void imageResize(int Components=4) (
169   scope void delegate (Color[] destrow, int y) srcGetRow,
170   scope void delegate (int y, const(Color)[] row) dstPutRow,
171   int srcwdt, int srchgt, int dstwdt, int dsthgt,
172   int filter=-1, float gamma=1.0f, float filterScale=1.0f
173 ) {
174   static assert(Components == 1 || Components == 3 || Components == 4, "invalid number of components in color");
175   assert(srcGetRow !is null);
176   assert(dstPutRow !is null);
177 
178   if (srcwdt < 1 || srchgt < 1 || dstwdt < 1 || dsthgt < 1 ||
179       srcwdt > ImageResizeMaxDimension || srchgt > ImageResizeMaxDimension ||
180       dstwdt > ImageResizeMaxDimension || dsthgt > ImageResizeMaxDimension) throw new Exception("invalid image size");
181 
182   if (filter < 0 || filter >= NumFilters) {
183     filter = resamplerFindFilterInternal(ImageResizeDefaultFilter);
184     if (filter < 0) {
185       filter = resamplerFindFilterInternal("lanczos4");
186     }
187   }
188   assert(filter >= 0 && filter < NumFilters);
189 
190 
191   if (lastGamma != gamma) {
192     version(iresample_debug) { import core.stdc.stdio; stderr.fprintf("creating translation tables for gamma %f (previous gamma is %f)\n", gamma, lastGamma); }
193     foreach (immutable i, ref v; srgb2linear[]) {
194       import std.math : pow;
195       v = cast(float)pow(cast(int)i*1.0f/255.0f, gamma);
196     }
197     immutable float invSourceGamma = 1.0f/gamma;
198     foreach (immutable i, ref v; linear2srgb[]) {
199       import std.math : pow;
200       int k = cast(int)(255.0f*pow(cast(int)i*InvLinear2srgbTableSize, invSourceGamma)+0.5f);
201       if (k < 0) k = 0; else if (k > 255) k = 255;
202       v = cast(ubyte)k;
203     }
204     lastGamma = gamma;
205   }
206   version(iresample_debug) { import core.stdc.stdio; stderr.fprintf("filter is %d\n", filter); }
207 
208   ImageResampleWorker[Components] resamplers;
209   float[][Components] samples;
210   Color[] srcrow, dstrow;
211   scope(exit) {
212     foreach (ref rsm; resamplers[]) .destroy(rsm);
213     foreach (ref smr; samples[]) .destroy(smr);
214   }
215 
216   // now create a ImageResampleWorker instance for each component to process
217   // the first instance will create new contributor tables, which are shared by the resamplers
218   // used for the other components (a memory and slight cache efficiency optimization).
219   resamplers[0] = new ImageResampleWorker(srcwdt, srchgt, dstwdt, dsthgt, ImageResampleWorker.BoundaryClamp, 0.0f, 1.0f, filter, null, null, filterScale, filterScale);
220   samples[0].length = srcwdt;
221   srcrow.length = srcwdt;
222   dstrow.length = dstwdt;
223   foreach (immutable i; 1..Components) {
224     resamplers[i] = new ImageResampleWorker(srcwdt, srchgt, dstwdt, dsthgt, ImageResampleWorker.BoundaryClamp, 0.0f, 1.0f, filter, resamplers[0].getClistX(), resamplers[0].getClistY(), filterScale, filterScale);
225     samples[i].length = srcwdt;
226   }
227 
228   int dsty = 0;
229   foreach (immutable int srcy; 0..srchgt) {
230     // get row components
231     srcGetRow(srcrow, srcy);
232     {
233       auto scp = srcrow.ptr;
234       foreach (immutable x; 0..srcwdt) {
235         auto sc = *scp++;
236         samples.ptr[0].ptr[x] = srgb2linear.ptr[sc.r]; // first component
237         static if (Components > 1) samples.ptr[1].ptr[x] = srgb2linear.ptr[sc.g]; // second component
238         static if (Components > 2) samples.ptr[2].ptr[x] = srgb2linear.ptr[sc.b]; // thirs component
239         static if (Components == 4) samples.ptr[3].ptr[x] = sc.a*(1.0f/255.0f); // fourth component is alpha, and it is already linear
240       }
241     }
242 
243     foreach (immutable c; 0..Components) if (!resamplers.ptr[c].putLine(samples.ptr[c].ptr)) assert(0, "out of memory");
244 
245     for (;;) {
246       int compIdx = 0;
247       for (; compIdx < Components; ++compIdx) {
248         const(float)* outsmp = resamplers.ptr[compIdx].getLine();
249         if (outsmp is null) break;
250         auto dsc = dstrow.ptr;
251         // alpha?
252         static if (Components == 4) {
253           if (compIdx == 3) {
254             foreach (immutable x; 0..dstwdt) {
255               dsc.a = Color.clampToByte(cast(int)(255.0f*(*outsmp++)+0.5f));
256               ++dsc;
257             }
258             continue;
259           }
260         }
261         // color
262         auto dsb = (cast(ubyte*)dsc)+compIdx;
263         foreach (immutable x; 0..dstwdt) {
264           int j = cast(int)(Linear2srgbTableSize*(*outsmp++)+0.5f);
265           if (j < 0) j = 0; else if (j >= Linear2srgbTableSize) j = Linear2srgbTableSize-1;
266           *dsb = linear2srgb.ptr[j];
267           dsb += 4;
268         }
269       }
270       if (compIdx < Components) break;
271       // fill destination line
272       assert(dsty < dsthgt);
273       static if (Components != 4) {
274         auto dsc = dstrow.ptr;
275         foreach (immutable x; 0..dstwdt) {
276           static if (Components == 1) dsc.g = dsc.b = dsc.r;
277           dsc.a = 255;
278           ++dsc;
279         }
280       }
281       //version(iresample_debug) { import core.stdc.stdio; stderr.fprintf("writing dest row %d with %u components\n", dsty, Components); }
282       dstPutRow(dsty, dstrow);
283       ++dsty;
284     }
285   }
286 }
287 
288 
289 // ////////////////////////////////////////////////////////////////////////// //
290 public final class ImageResampleWorker {
291 nothrow @trusted @nogc:
292 public:
293   alias ResampleReal = float;
294   alias Sample = ResampleReal;
295 
296   static struct Contrib {
297     ResampleReal weight;
298     ushort pixel;
299   }
300 
301   static struct ContribList {
302     ushort n;
303     Contrib* p;
304   }
305 
306   alias BoundaryOp = int;
307   enum /*Boundary_Op*/ {
308     BoundaryWrap = 0,
309     BoundaryReflect = 1,
310     BoundaryClamp = 2,
311   }
312 
313   alias Status = int;
314   enum /*Status*/ {
315     StatusOkay = 0,
316     StatusOutOfMemory = 1,
317     StatusBadFilterName = 2,
318     StatusScanBufferFull = 3,
319   }
320 
321 private:
322   alias FilterFunc = ResampleReal function (ResampleReal) nothrow @trusted @nogc;
323 
324   int mIntermediateX;
325 
326   int mResampleSrcX;
327   int mResampleSrcY;
328   int mResampleDstX;
329   int mResampleDstY;
330 
331   BoundaryOp mBoundaryOp;
332 
333   Sample* mPdstBuf;
334   Sample* mPtmpBuf;
335 
336   ContribList* mPclistX;
337   ContribList* mPclistY;
338 
339   bool mClistXForced;
340   bool mClistYForced;
341 
342   bool mDelayXResample;
343 
344   int* mPsrcYCount;
345   ubyte* mPsrcYFlag;
346 
347   // The maximum number of scanlines that can be buffered at one time.
348   enum MaxScanBufSize = ImageResizeMaxDimension;
349 
350   static struct ScanBuf {
351     int[MaxScanBufSize] scanBufY;
352     Sample*[MaxScanBufSize] scanBufL;
353   }
354 
355   ScanBuf* mPscanBuf;
356 
357   int mCurSrcY;
358   int mCurDstY;
359 
360   Status mStatus;
361 
362   // The make_clist() method generates, for all destination samples,
363   // the list of all source samples with non-zero weighted contributions.
364   ContribList* makeClist(
365     int srcX, int dstX, BoundaryOp boundaryOp,
366     FilterFunc Pfilter,
367     ResampleReal filterSupport,
368     ResampleReal filterScale,
369     ResampleReal srcOfs)
370   {
371     import core.stdc.stdlib : calloc, free;
372     import std.math : floor, ceil;
373 
374     static struct ContribBounds {
375       // The center of the range in DISCRETE coordinates (pixel center = 0.0f).
376       ResampleReal center;
377       int left, right;
378     }
379 
380     ContribList* Pcontrib, PcontribRes;
381     Contrib* Pcpool;
382     Contrib* PcpoolNext;
383     ContribBounds* PcontribBounds;
384 
385     if ((Pcontrib = cast(ContribList*)calloc(dstX, ContribList.sizeof)) is null) return null;
386     scope(exit) if (Pcontrib !is null) free(Pcontrib);
387 
388     PcontribBounds = cast(ContribBounds*)calloc(dstX, ContribBounds.sizeof);
389     if (PcontribBounds is null) return null;
390     scope(exit) free(PcontribBounds);
391 
392     enum ResampleReal NUDGE = 0.5f;
393     immutable ResampleReal ooFilterScale = 1.0f/filterScale;
394     immutable ResampleReal xscale = dstX/cast(ResampleReal)srcX;
395 
396     if (xscale < 1.0f) {
397       int total = 0;
398       // Handle case when there are fewer destination samples than source samples (downsampling/minification).
399       // stretched half width of filter
400       immutable ResampleReal halfWidth = (filterSupport/xscale)*filterScale;
401       // Find the range of source sample(s) that will contribute to each destination sample.
402       foreach (immutable i; 0..dstX) {
403         // Convert from discrete to continuous coordinates, scale, then convert back to discrete.
404         ResampleReal center = (cast(ResampleReal)i+NUDGE)/xscale;
405         center -= NUDGE;
406         center += srcOfs;
407         immutable int left = castToInt(cast(ResampleReal)floor(center-halfWidth));
408         immutable int right = castToInt(cast(ResampleReal)ceil(center+halfWidth));
409         PcontribBounds[i].center = center;
410         PcontribBounds[i].left = left;
411         PcontribBounds[i].right = right;
412         total += (right-left+1);
413       }
414 
415       // Allocate memory for contributors.
416       if (total == 0 || ((Pcpool = cast(Contrib*)calloc(total, Contrib.sizeof)) is null)) return null;
417       //scope(failure) free(Pcpool);
418       //immutable int total = n;
419 
420       PcpoolNext = Pcpool;
421 
422       // Create the list of source samples which contribute to each destination sample.
423       foreach (immutable i; 0..dstX) {
424         int maxK = -1;
425         ResampleReal maxW = -1e+20f;
426 
427         ResampleReal center = PcontribBounds[i].center;
428         immutable int left = PcontribBounds[i].left;
429         immutable int right = PcontribBounds[i].right;
430 
431         Pcontrib[i].n = 0;
432         Pcontrib[i].p = PcpoolNext;
433         PcpoolNext += (right-left+1);
434         assert(PcpoolNext-Pcpool <= total);
435 
436         ResampleReal totalWeight0 = 0;
437         foreach (immutable j; left..right+1) totalWeight0 += Pfilter((center-cast(ResampleReal)j)*xscale*ooFilterScale);
438         immutable ResampleReal norm = cast(ResampleReal)(1.0f/totalWeight0);
439 
440         ResampleReal totalWeight1 = 0;
441         foreach (immutable j; left..right+1) {
442           immutable ResampleReal weight = Pfilter((center-cast(ResampleReal)j)*xscale*ooFilterScale)*norm;
443           if (weight == 0.0f) continue;
444           immutable int n = reflect(j, srcX, boundaryOp);
445           // Increment the number of source samples which contribute to the current destination sample.
446           immutable int k = Pcontrib[i].n++;
447           Pcontrib[i].p[k].pixel = cast(ushort)(n); // store src sample number
448           Pcontrib[i].p[k].weight = weight; // store src sample weight
449           totalWeight1 += weight; // total weight of all contributors
450           if (weight > maxW) {
451             maxW = weight;
452             maxK = k;
453           }
454         }
455         //assert(Pcontrib[i].n);
456         //assert(max_k != -1);
457         if (maxK == -1 || Pcontrib[i].n == 0) return null;
458         if (totalWeight1 != 1.0f) Pcontrib[i].p[maxK].weight += 1.0f-totalWeight1;
459       }
460     } else {
461       int total = 0;
462       // Handle case when there are more destination samples than source samples (upsampling).
463       immutable ResampleReal halfWidth = filterSupport*filterScale;
464       // Find the source sample(s) that contribute to each destination sample.
465       foreach (immutable i; 0..dstX) {
466         // Convert from discrete to continuous coordinates, scale, then convert back to discrete.
467         ResampleReal center = (cast(ResampleReal)i+NUDGE)/xscale;
468         center -= NUDGE;
469         center += srcOfs;
470         immutable int left = castToInt(cast(ResampleReal)floor(center-halfWidth));
471         immutable int right = castToInt(cast(ResampleReal)ceil(center+halfWidth));
472         PcontribBounds[i].center = center;
473         PcontribBounds[i].left = left;
474         PcontribBounds[i].right = right;
475         total += (right-left+1);
476       }
477 
478       // Allocate memory for contributors.
479       if (total == 0 || ((Pcpool = cast(Contrib*)calloc(total, Contrib.sizeof)) is null)) return null;
480       //scope(failure) free(Pcpool);
481 
482       PcpoolNext = Pcpool;
483 
484       // Create the list of source samples which contribute to each destination sample.
485       foreach (immutable i; 0..dstX) {
486         int maxK = -1;
487         ResampleReal maxW = -1e+20f;
488 
489         ResampleReal center = PcontribBounds[i].center;
490         immutable int left = PcontribBounds[i].left;
491         immutable int right = PcontribBounds[i].right;
492 
493         Pcontrib[i].n = 0;
494         Pcontrib[i].p = PcpoolNext;
495         PcpoolNext += (right-left+1);
496         assert(PcpoolNext-Pcpool <= total);
497 
498         ResampleReal totalWeight0 = 0;
499         foreach (immutable j; left..right+1) totalWeight0 += Pfilter((center-cast(ResampleReal)j)*ooFilterScale);
500         immutable ResampleReal norm = cast(ResampleReal)(1.0f/totalWeight0);
501 
502         ResampleReal totalWeight1 = 0;
503         foreach (immutable j; left..right+1) {
504           immutable ResampleReal weight = Pfilter((center-cast(ResampleReal)j)*ooFilterScale)*norm;
505           if (weight == 0.0f) continue;
506           immutable int n = reflect(j, srcX, boundaryOp);
507           // Increment the number of source samples which contribute to the current destination sample.
508           immutable int k = Pcontrib[i].n++;
509           Pcontrib[i].p[k].pixel = cast(ushort)(n); // store src sample number
510           Pcontrib[i].p[k].weight = weight; // store src sample weight
511           totalWeight1 += weight; // total weight of all contributors
512           if (weight > maxW) {
513             maxW = weight;
514             maxK = k;
515           }
516         }
517         //assert(Pcontrib[i].n);
518         //assert(max_k != -1);
519         if (maxK == -1 || Pcontrib[i].n == 0) return null;
520         if (totalWeight1 != 1.0f) Pcontrib[i].p[maxK].weight += 1.0f-totalWeight1;
521       }
522     }
523     // don't free return value
524     PcontribRes = Pcontrib;
525     Pcontrib = null;
526     return PcontribRes;
527   }
528 
529   static int countOps (const(ContribList)* Pclist, int k) {
530     int t = 0;
531     foreach (immutable i; 0..k) t += Pclist[i].n;
532     return t;
533   }
534 
535   private ResampleReal mLo;
536   private ResampleReal mHi;
537 
538   ResampleReal clampSample (ResampleReal f) const {
539     pragma(inline, true);
540     if (f < mLo) f = mLo; else if (f > mHi) f = mHi;
541     return f;
542   }
543 
544 public:
545   // src_x/src_y - Input dimensions
546   // dst_x/dst_y - Output dimensions
547   // boundary_op - How to sample pixels near the image boundaries
548   // sample_low/sample_high - Clamp output samples to specified range, or disable clamping if sample_low >= sample_high
549   // Pclist_x/Pclist_y - Optional pointers to contributor lists from another instance of a ImageResampleWorker
550   // src_x_ofs/src_y_ofs - Offset input image by specified amount (fractional values okay)
551   this(
552     int srcX, int srcY,
553     int dstX, int dstY,
554     BoundaryOp boundaryOp=BoundaryClamp,
555     ResampleReal sampleLow=0.0f, ResampleReal sampleHigh=0.0f,
556     int PfilterIndex=-1,
557     ContribList* PclistX=null,
558     ContribList* PclistY=null,
559     ResampleReal filterXScale=1.0f,
560     ResampleReal filterYScale=1.0f,
561     ResampleReal srcXOfs=0.0f,
562     ResampleReal srcYOfs=0.0f)
563   {
564     import core.stdc.stdlib : calloc, malloc;
565 
566     int i, j;
567     ResampleReal support;
568     FilterFunc func;
569 
570     assert(srcX > 0);
571     assert(srcY > 0);
572     assert(dstX > 0);
573     assert(dstY > 0);
574 
575     mLo = sampleLow;
576     mHi = sampleHigh;
577 
578     mDelayXResample = false;
579     mIntermediateX = 0;
580     mPdstBuf = null;
581     mPtmpBuf = null;
582     mClistXForced = false;
583     mPclistX = null;
584     mClistYForced = false;
585     mPclistY = null;
586     mPsrcYCount = null;
587     mPsrcYFlag = null;
588     mPscanBuf = null;
589     mStatus = StatusOkay;
590 
591     mResampleSrcX = srcX;
592     mResampleSrcY = srcY;
593     mResampleDstX = dstX;
594     mResampleDstY = dstY;
595 
596     mBoundaryOp = boundaryOp;
597 
598     if ((mPdstBuf = cast(Sample*)malloc(mResampleDstX*Sample.sizeof)) is null) {
599       mStatus = StatusOutOfMemory;
600       return;
601     }
602 
603     if (PfilterIndex < 0 || PfilterIndex >= NumFilters) {
604       PfilterIndex = resamplerFindFilterInternal(ImageResizeDefaultFilter);
605       if (PfilterIndex < 0 || PfilterIndex >= NumFilters) {
606         mStatus = StatusBadFilterName;
607         return;
608       }
609     }
610 
611     func = gFilters[PfilterIndex].func;
612     support = gFilters[PfilterIndex].support;
613 
614     // Create contributor lists, unless the user supplied custom lists.
615     if (PclistX is null) {
616       mPclistX = makeClist(mResampleSrcX, mResampleDstX, mBoundaryOp, func, support, filterXScale, srcXOfs);
617       if (mPclistX is null) {
618         mStatus = StatusOutOfMemory;
619         return;
620       }
621     } else {
622       mPclistX = PclistX;
623       mClistXForced = true;
624     }
625 
626     if (PclistY is null) {
627       mPclistY = makeClist(mResampleSrcY, mResampleDstY, mBoundaryOp, func, support, filterYScale, srcYOfs);
628       if (mPclistY is null) {
629         mStatus = StatusOutOfMemory;
630         return;
631       }
632     } else {
633       mPclistY = PclistY;
634       mClistYForced = true;
635     }
636 
637     if ((mPsrcYCount = cast(int*)calloc(mResampleSrcY, int.sizeof)) is null) {
638       mStatus = StatusOutOfMemory;
639       return;
640     }
641 
642     if ((mPsrcYFlag = cast(ubyte*)calloc(mResampleSrcY, ubyte.sizeof)) is null) {
643       mStatus = StatusOutOfMemory;
644       return;
645     }
646 
647     // Count how many times each source line contributes to a destination line.
648     for (i = 0; i < mResampleDstY; ++i) {
649       for (j = 0; j < mPclistY[i].n; ++j) {
650         ++mPsrcYCount[resamplerRangeCheck(mPclistY[i].p[j].pixel, mResampleSrcY)];
651       }
652     }
653 
654     if ((mPscanBuf = cast(ScanBuf*)malloc(ScanBuf.sizeof)) is null) {
655       mStatus = StatusOutOfMemory;
656       return;
657     }
658 
659     for (i = 0; i < MaxScanBufSize; ++i) {
660       mPscanBuf.scanBufY.ptr[i] = -1;
661       mPscanBuf.scanBufL.ptr[i] = null;
662     }
663 
664     mCurSrcY = mCurDstY = 0;
665     {
666       // Determine which axis to resample first by comparing the number of multiplies required
667       // for each possibility.
668       int xOps = countOps(mPclistX, mResampleDstX);
669       int yOps = countOps(mPclistY, mResampleDstY);
670 
671       // Hack 10/2000: Weight Y axis ops a little more than X axis ops.
672       // (Y axis ops use more cache resources.)
673       int xyOps = xOps*mResampleSrcY+(4*yOps*mResampleDstX)/3;
674       int yxOps = (4*yOps*mResampleSrcX)/3+xOps*mResampleDstY;
675 
676       // Now check which resample order is better. In case of a tie, choose the order
677       // which buffers the least amount of data.
678       if (xyOps > yxOps || (xyOps == yxOps && mResampleSrcX < mResampleDstX)) {
679         mDelayXResample = true;
680         mIntermediateX = mResampleSrcX;
681       } else {
682         mDelayXResample = false;
683         mIntermediateX = mResampleDstX;
684       }
685     }
686 
687     if (mDelayXResample) {
688       if ((mPtmpBuf = cast(Sample*)malloc(mIntermediateX*Sample.sizeof)) is null) {
689         mStatus = StatusOutOfMemory;
690         return;
691       }
692     }
693   }
694 
695   ~this () {
696      import core.stdc.stdlib : free;
697 
698      if (mPdstBuf !is null) {
699        free(mPdstBuf);
700        mPdstBuf = null;
701      }
702 
703      if (mPtmpBuf !is null) {
704        free(mPtmpBuf);
705        mPtmpBuf = null;
706      }
707 
708      // Don't deallocate a contibutor list if the user passed us one of their own.
709      if (mPclistX !is null && !mClistXForced) {
710        free(mPclistX.p);
711        free(mPclistX);
712        mPclistX = null;
713      }
714      if (mPclistY !is null && !mClistYForced) {
715        free(mPclistY.p);
716        free(mPclistY);
717        mPclistY = null;
718      }
719 
720      if (mPsrcYCount !is null) {
721        free(mPsrcYCount);
722        mPsrcYCount = null;
723      }
724 
725      if (mPsrcYFlag !is null) {
726        free(mPsrcYFlag);
727        mPsrcYFlag = null;
728      }
729 
730      if (mPscanBuf !is null) {
731        foreach (immutable i; 0..MaxScanBufSize) if (mPscanBuf.scanBufL.ptr[i] !is null) free(mPscanBuf.scanBufL.ptr[i]);
732        free(mPscanBuf);
733        mPscanBuf = null;
734      }
735   }
736 
737   // Reinits resampler so it can handle another frame.
738   void restart () {
739     import core.stdc.stdlib : free;
740     if (StatusOkay != mStatus) return;
741     mCurSrcY = mCurDstY = 0;
742     foreach (immutable i; 0..mResampleSrcY) {
743       mPsrcYCount[i] = 0;
744       mPsrcYFlag[i] = false;
745     }
746     foreach (immutable i; 0..mResampleDstY) {
747       foreach (immutable j; 0..mPclistY[i].n) {
748         ++mPsrcYCount[resamplerRangeCheck(mPclistY[i].p[j].pixel, mResampleSrcY)];
749       }
750     }
751     foreach (immutable i; 0..MaxScanBufSize) {
752       mPscanBuf.scanBufY.ptr[i] = -1;
753       free(mPscanBuf.scanBufL.ptr[i]);
754       mPscanBuf.scanBufL.ptr[i] = null;
755     }
756   }
757 
758   // false on out of memory.
759   bool putLine (const(Sample)* Psrc) {
760     int i;
761 
762     if (mCurSrcY >= mResampleSrcY) return false;
763 
764     // Does this source line contribute to any destination line? if not, exit now.
765     if (!mPsrcYCount[resamplerRangeCheck(mCurSrcY, mResampleSrcY)]) {
766       ++mCurSrcY;
767       return true;
768     }
769 
770     // Find an empty slot in the scanline buffer. (FIXME: Perf. is terrible here with extreme scaling ratios.)
771     for (i = 0; i < MaxScanBufSize; ++i) if (mPscanBuf.scanBufY.ptr[i] == -1) break;
772 
773     // If the buffer is full, exit with an error.
774     if (i == MaxScanBufSize) {
775       mStatus = StatusScanBufferFull;
776       return false;
777     }
778 
779     mPsrcYFlag[resamplerRangeCheck(mCurSrcY, mResampleSrcY)] = true;
780     mPscanBuf.scanBufY.ptr[i] = mCurSrcY;
781 
782     // Does this slot have any memory allocated to it?
783     if (!mPscanBuf.scanBufL.ptr[i]) {
784       import core.stdc.stdlib : malloc;
785       if ((mPscanBuf.scanBufL.ptr[i] = cast(Sample*)malloc(mIntermediateX*Sample.sizeof)) is null) {
786         mStatus = StatusOutOfMemory;
787         return false;
788       }
789     }
790 
791     // Resampling on the X axis first?
792     if (mDelayXResample) {
793       import core.stdc.string : memcpy;
794       assert(mIntermediateX == mResampleSrcX);
795       // Y-X resampling order
796       memcpy(mPscanBuf.scanBufL.ptr[i], Psrc, mIntermediateX*Sample.sizeof);
797     } else {
798       assert(mIntermediateX == mResampleDstX);
799       // X-Y resampling order
800       resampleX(mPscanBuf.scanBufL.ptr[i], Psrc);
801     }
802 
803     ++mCurSrcY;
804 
805     return true;
806   }
807 
808   // null if no scanlines are currently available (give the resampler more scanlines!)
809   const(Sample)* getLine () {
810     // if all the destination lines have been generated, then always return null
811     if (mCurDstY == mResampleDstY) return null;
812     // check to see if all the required contributors are present, if not, return null
813     foreach (immutable i; 0..mPclistY[mCurDstY].n) {
814       if (!mPsrcYFlag[resamplerRangeCheck(mPclistY[mCurDstY].p[i].pixel, mResampleSrcY)]) return null;
815     }
816     resampleY(mPdstBuf);
817     ++mCurDstY;
818     return mPdstBuf;
819   }
820 
821   @property Status status () const { pragma(inline, true); return mStatus; }
822 
823   // returned contributor lists can be shared with another ImageResampleWorker
824   void getClists (ContribList** ptrClistX, ContribList** ptrClistY) {
825     if (ptrClistX !is null) *ptrClistX = mPclistX;
826     if (ptrClistY !is null) *ptrClistY = mPclistY;
827   }
828 
829   @property ContribList* getClistX () { pragma(inline, true); return mPclistX; }
830   @property ContribList* getClistY () { pragma(inline, true); return mPclistY; }
831 
832   // filter accessors
833   static @property auto filters () {
834     static struct FilterRange {
835     pure nothrow @trusted @nogc:
836       int idx;
837       @property bool empty () const { pragma(inline, true); return (idx >= NumFilters); }
838       @property string front () const { pragma(inline, true); return (idx < NumFilters ? gFilters[idx].name : null); }
839       void popFront () { if (idx < NumFilters) ++idx; }
840       int length () const { return cast(int)NumFilters; }
841       alias opDollar = length;
842     }
843     return FilterRange();
844   }
845 
846 private:
847   /* Ensure that the contributing source sample is
848   * within bounds. If not, reflect, clamp, or wrap.
849   */
850   int reflect (in int j, in int srcX, in BoundaryOp boundaryOp) {
851     int n;
852     if (j < 0) {
853       if (boundaryOp == BoundaryReflect) {
854         n = -j;
855         if (n >= srcX) n = srcX-1;
856       } else if (boundaryOp == BoundaryWrap) {
857         n = posmod(j, srcX);
858       } else {
859         n = 0;
860       }
861     } else if (j >= srcX) {
862       if (boundaryOp == BoundaryReflect) {
863         n = (srcX-j)+(srcX-1);
864         if (n < 0) n = 0;
865       } else if (boundaryOp == BoundaryWrap) {
866         n = posmod(j, srcX);
867       } else {
868         n = srcX-1;
869       }
870     } else {
871       n = j;
872     }
873     return n;
874   }
875 
876   void resampleX (Sample* Pdst, const(Sample)* Psrc) {
877     assert(Pdst);
878     assert(Psrc);
879 
880     Sample total;
881     ContribList *Pclist = mPclistX;
882     Contrib *p;
883 
884     for (int i = mResampleDstX; i > 0; --i, ++Pclist) {
885       int j = void;
886       for (j = Pclist.n, p = Pclist.p, total = 0; j > 0; --j, ++p) total += Psrc[p.pixel]*p.weight;
887       *Pdst++ = total;
888     }
889   }
890 
891   void scaleYMov (Sample* Ptmp, const(Sample)* Psrc, ResampleReal weight, int dstX) {
892     // Not += because temp buf wasn't cleared.
893     for (int i = dstX; i > 0; --i) *Ptmp++ = *Psrc++*weight;
894   }
895 
896   void scaleYAdd (Sample* Ptmp, const(Sample)* Psrc, ResampleReal weight, int dstX) {
897     for (int i = dstX; i > 0; --i) (*Ptmp++) += *Psrc++*weight;
898   }
899 
900   void clamp (Sample* Pdst, int n) {
901     while (n > 0) {
902       *Pdst = clampSample(*Pdst);
903       ++Pdst;
904       --n;
905     }
906   }
907 
908   void resampleY (Sample* Pdst) {
909     Sample* Psrc;
910     ContribList* Pclist = &mPclistY[mCurDstY];
911 
912     Sample* Ptmp = mDelayXResample ? mPtmpBuf : Pdst;
913     assert(Ptmp);
914 
915     // process each contributor
916     foreach (immutable i; 0..Pclist.n) {
917       // locate the contributor's location in the scan buffer -- the contributor must always be found!
918       int j = void;
919       for (j = 0; j < MaxScanBufSize; ++j) if (mPscanBuf.scanBufY.ptr[j] == Pclist.p[i].pixel) break;
920       assert(j < MaxScanBufSize);
921       Psrc = mPscanBuf.scanBufL.ptr[j];
922       if (!i) {
923         scaleYMov(Ptmp, Psrc, Pclist.p[i].weight, mIntermediateX);
924       } else {
925         scaleYAdd(Ptmp, Psrc, Pclist.p[i].weight, mIntermediateX);
926       }
927 
928       /* If this source line doesn't contribute to any
929        * more destination lines then mark the scanline buffer slot
930        * which holds this source line as free.
931        * (The max. number of slots used depends on the Y
932        * axis sampling factor and the scaled filter width.)
933        */
934 
935       if (--mPsrcYCount[resamplerRangeCheck(Pclist.p[i].pixel, mResampleSrcY)] == 0) {
936         mPsrcYFlag[resamplerRangeCheck(Pclist.p[i].pixel, mResampleSrcY)] = false;
937         mPscanBuf.scanBufY.ptr[j] = -1;
938       }
939     }
940 
941     // now generate the destination line
942     if (mDelayXResample) {
943       // X was resampling delayed until after Y resampling
944       assert(Pdst != Ptmp);
945       resampleX(Pdst, Ptmp);
946     } else {
947       assert(Pdst == Ptmp);
948     }
949 
950     if (mLo < mHi) clamp(Pdst, mResampleDstX);
951   }
952 }
953 
954 
955 // ////////////////////////////////////////////////////////////////////////// //
956 private nothrow @trusted @nogc:
957 int resamplerRangeCheck (int v, int h) {
958   version(assert) {
959     //import std.conv : to;
960     //assert(v >= 0 && v < h, "invalid v ("~to!string(v)~"), should be in [0.."~to!string(h)~")");
961     assert(v >= 0 && v < h); // alas, @nogc
962     return v;
963   } else {
964     pragma(inline, true);
965     return v;
966   }
967 }
968 
969 enum M_PI = 3.14159265358979323846;
970 
971 // Float to int cast with truncation.
972 int castToInt (ImageResampleWorker.ResampleReal i) { pragma(inline, true); return cast(int)i; }
973 
974 // (x mod y) with special handling for negative x values.
975 int posmod (int x, int y) {
976   pragma(inline, true);
977   if (x >= 0) {
978     return (x%y);
979   } else {
980     int m = (-x)%y;
981     if (m != 0) m = y-m;
982     return m;
983   }
984 }
985 
986 // To add your own filter, insert the new function below and update the filter table.
987 // There is no need to make the filter function particularly fast, because it's
988 // only called during initializing to create the X and Y axis contributor tables.
989 
990 /* pulse/Fourier window */
991 enum BoxFilterSupport = 0.5f;
992 ImageResampleWorker.ResampleReal boxFilter (ImageResampleWorker.ResampleReal t) {
993   // make_clist() calls the filter function with t inverted (pos = left, neg = right)
994   if (t >= -0.5f && t < 0.5f) return 1.0f; else return 0.0f;
995 }
996 
997 /* box (*) box, bilinear/triangle */
998 enum TentFilterSupport = 1.0f;
999 ImageResampleWorker.ResampleReal tentFilter (ImageResampleWorker.ResampleReal t) {
1000   if (t < 0.0f) t = -t;
1001   if (t < 1.0f) return 1.0f-t; else return 0.0f;
1002 }
1003 
1004 /* box (*) box (*) box */
1005 enum BellSupport = 1.5f;
1006 ImageResampleWorker.ResampleReal bellFilter (ImageResampleWorker.ResampleReal t) {
1007   if (t < 0.0f) t = -t;
1008   if (t < 0.5f) return (0.75f-(t*t));
1009   if (t < 1.5f) { t = (t-1.5f); return (0.5f*(t*t)); }
1010   return (0.0f);
1011 }
1012 
1013 /* box (*) box (*) box (*) box */
1014 enum BSplineSupport = 2.0f;
1015 ImageResampleWorker.ResampleReal BSplineFilter (ImageResampleWorker.ResampleReal t) {
1016   if (t < 0.0f) t = -t;
1017   if (t < 1.0f) { immutable ImageResampleWorker.ResampleReal tt = t*t; return ((0.5f*tt*t)-tt+(2.0f/3.0f)); }
1018   if (t < 2.0f) { t = 2.0f-t; return ((1.0f/6.0f)*(t*t*t)); }
1019   return 0.0f;
1020 }
1021 
1022 // Dodgson, N., "Quadratic Interpolation for Image Resampling"
1023 enum QuadraticSupport = 1.5f;
1024 ImageResampleWorker.ResampleReal quadratic (ImageResampleWorker.ResampleReal t, in ImageResampleWorker.ResampleReal R) {
1025   pragma(inline, true);
1026   if (t < 0.0f) t = -t;
1027   if (t < QuadraticSupport) {
1028     immutable ImageResampleWorker.ResampleReal tt = t*t;
1029     if (t <= 0.5f) return (-2.0f*R)*tt+0.5f*(R+1.0f);
1030     return (R*tt)+(-2.0f*R-0.5f)*t+(3.0f/4.0f)*(R+1.0f);
1031   }
1032   return 0.0f;
1033 }
1034 
1035 ImageResampleWorker.ResampleReal quadraticInterpFilter (ImageResampleWorker.ResampleReal t) {
1036   return quadratic(t, 1.0f);
1037 }
1038 
1039 ImageResampleWorker.ResampleReal quadraticApproxFilter (ImageResampleWorker.ResampleReal t) {
1040   return quadratic(t, 0.5f);
1041 }
1042 
1043 ImageResampleWorker.ResampleReal quadraticMixFilter (ImageResampleWorker.ResampleReal t) {
1044   return quadratic(t, 0.8f);
1045 }
1046 
1047 // Mitchell, D. and A. Netravali, "Reconstruction Filters in Computer Graphics."
1048 // Computer Graphics, Vol. 22, No. 4, pp. 221-228.
1049 // (B, C)
1050 // (1/3, 1/3)  - Defaults recommended by Mitchell and Netravali
1051 // (1, 0)    - Equivalent to the Cubic B-Spline
1052 // (0, 0.5)   - Equivalent to the Catmull-Rom Spline
1053 // (0, C)   - The family of Cardinal Cubic Splines
1054 // (B, 0)   - Duff's tensioned B-Splines.
1055 ImageResampleWorker.ResampleReal mitchell (ImageResampleWorker.ResampleReal t, in ImageResampleWorker.ResampleReal B, in ImageResampleWorker.ResampleReal C) {
1056   ImageResampleWorker.ResampleReal tt = t*t;
1057   if (t < 0.0f) t = -t;
1058   if (t < 1.0f) {
1059     t = (((12.0f-9.0f*B-6.0f*C)*(t*tt))+
1060          ((-18.0f+12.0f*B+6.0f*C)*tt)+
1061          (6.0f-2.0f*B));
1062     return (t/6.0f);
1063   }
1064   if (t < 2.0f) {
1065     t = (((-1.0f*B-6.0f*C)*(t*tt))+
1066          ((6.0f*B+30.0f*C)*tt)+
1067          ((-12.0f*B-48.0f*C)*t)+
1068          (8.0f*B+24.0f*C));
1069     return (t/6.0f);
1070   }
1071   return 0.0f;
1072 }
1073 
1074 enum MitchellSupport = 2.0f;
1075 ImageResampleWorker.ResampleReal mitchellFilter (ImageResampleWorker.ResampleReal t) {
1076   return mitchell(t, 1.0f/3.0f, 1.0f/3.0f);
1077 }
1078 
1079 enum CatmullRomSupport = 2.0f;
1080 ImageResampleWorker.ResampleReal catmullRomFilter (ImageResampleWorker.ResampleReal t) {
1081   return mitchell(t, 0.0f, 0.5f);
1082 }
1083 
1084 double sinc (double x) {
1085   pragma(inline, true);
1086   import std.math : sin;
1087   x *= M_PI;
1088   if (x < 0.01f && x > -0.01f) return 1.0f+x*x*(-1.0f/6.0f+x*x*1.0f/120.0f);
1089   return sin(x)/x;
1090 }
1091 
1092 ImageResampleWorker.ResampleReal clean (double t) {
1093   pragma(inline, true);
1094   import std.math : abs;
1095   enum EPSILON = cast(ImageResampleWorker.ResampleReal)0.0000125f;
1096   if (abs(t) < EPSILON) return 0.0f;
1097   return cast(ImageResampleWorker.ResampleReal)t;
1098 }
1099 
1100 //static double blackman_window(double x)
1101 //{
1102 //  return 0.42f+0.50f*cos(M_PI*x)+0.08f*cos(2.0f*M_PI*x);
1103 //}
1104 
1105 double blackmanExactWindow (double x) {
1106   pragma(inline, true);
1107   import std.math : cos;
1108   return 0.42659071f+0.49656062f*cos(M_PI*x)+0.07684867f*cos(2.0f*M_PI*x);
1109 }
1110 
1111 enum BlackmanSupport = 3.0f;
1112 ImageResampleWorker.ResampleReal blackmanFilter (ImageResampleWorker.ResampleReal t) {
1113   if (t < 0.0f) t = -t;
1114   if (t < 3.0f) {
1115     //return clean(sinc(t)*blackman_window(t/3.0f));
1116     return clean(sinc(t)*blackmanExactWindow(t/3.0f));
1117   }
1118   return (0.0f);
1119 }
1120 
1121 // with blackman window
1122 enum GaussianSupport = 1.25f;
1123 ImageResampleWorker.ResampleReal gaussianFilter (ImageResampleWorker.ResampleReal t) {
1124   import std.math : exp, sqrt;
1125   if (t < 0) t = -t;
1126   if (t < GaussianSupport) return clean(exp(-2.0f*t*t)*sqrt(2.0f/M_PI)*blackmanExactWindow(t/GaussianSupport));
1127   return 0.0f;
1128 }
1129 
1130 // Windowed sinc -- see "Jimm Blinn's Corner: Dirty Pixels" pg. 26.
1131 enum Lanczos3Support = 3.0f;
1132 ImageResampleWorker.ResampleReal lanczos3Filter (ImageResampleWorker.ResampleReal t) {
1133   if (t < 0.0f) t = -t;
1134   if (t < 3.0f) return clean(sinc(t)*sinc(t/3.0f));
1135   return (0.0f);
1136 }
1137 
1138 enum Lanczos4Support = 4.0f;
1139 ImageResampleWorker.ResampleReal lanczos4Filter (ImageResampleWorker.ResampleReal t) {
1140   if (t < 0.0f) t = -t;
1141   if (t < 4.0f) return clean(sinc(t)*sinc(t/4.0f));
1142   return (0.0f);
1143 }
1144 
1145 enum Lanczos6Support = 6.0f;
1146 ImageResampleWorker.ResampleReal lanczos6Filter (ImageResampleWorker.ResampleReal t) {
1147   if (t < 0.0f) t = -t;
1148   if (t < 6.0f) return clean(sinc(t)*sinc(t/6.0f));
1149   return (0.0f);
1150 }
1151 
1152 enum Lanczos12Support = 12.0f;
1153 ImageResampleWorker.ResampleReal lanczos12Filter (ImageResampleWorker.ResampleReal t) {
1154   if (t < 0.0f) t = -t;
1155   if (t < 12.0f) return clean(sinc(t)*sinc(t/12.0f));
1156   return (0.0f);
1157 }
1158 
1159 double bessel0 (double x) {
1160   enum EpsilonRatio = cast(double)1E-16;
1161   double xh = 0.5*x;
1162   double sum = 1.0;
1163   double pow = 1.0;
1164   int k = 0;
1165   double ds = 1.0;
1166   // FIXME: Shouldn't this stop after X iterations for max. safety?
1167   while (ds > sum*EpsilonRatio) {
1168     ++k;
1169     pow = pow*(xh/k);
1170     ds = pow*pow;
1171     sum = sum+ds;
1172   }
1173   return sum;
1174 }
1175 
1176 enum KaiserAlpha = cast(ImageResampleWorker.ResampleReal)4.0;
1177 double kaiser (double alpha, double halfWidth, double x) {
1178   pragma(inline, true);
1179   import std.math : sqrt;
1180   immutable double ratio = (x/halfWidth);
1181   return bessel0(alpha*sqrt(1-ratio*ratio))/bessel0(alpha);
1182 }
1183 
1184 enum KaiserSupport = 3;
1185 static ImageResampleWorker.ResampleReal kaiserFilter (ImageResampleWorker.ResampleReal t) {
1186   if (t < 0.0f) t = -t;
1187   if (t < KaiserSupport) {
1188     import std.math : exp, log;
1189     // db atten
1190     immutable ImageResampleWorker.ResampleReal att = 40.0f;
1191     immutable ImageResampleWorker.ResampleReal alpha = cast(ImageResampleWorker.ResampleReal)(exp(log(cast(double)0.58417*(att-20.96))*0.4)+0.07886*(att-20.96));
1192     //const ImageResampleWorker.Resample_Real alpha = KAISER_ALPHA;
1193     return cast(ImageResampleWorker.ResampleReal)clean(sinc(t)*kaiser(alpha, KaiserSupport, t));
1194   }
1195   return 0.0f;
1196 }
1197 
1198 // filters[] is a list of all the available filter functions.
1199 struct FilterInfo {
1200   string name;
1201   ImageResampleWorker.FilterFunc func;
1202   ImageResampleWorker.ResampleReal support;
1203 }
1204 
1205 static immutable FilterInfo[16] gFilters = [
1206    FilterInfo("box",              &boxFilter,             BoxFilterSupport),
1207    FilterInfo("tent",             &tentFilter,            TentFilterSupport),
1208    FilterInfo("bell",             &bellFilter,            BellSupport),
1209    FilterInfo("bspline",          &BSplineFilter,         BSplineSupport),
1210    FilterInfo("mitchell",         &mitchellFilter,        MitchellSupport),
1211    FilterInfo("lanczos3",         &lanczos3Filter,        Lanczos3Support),
1212    FilterInfo("blackman",         &blackmanFilter,        BlackmanSupport),
1213    FilterInfo("lanczos4",         &lanczos4Filter,        Lanczos4Support),
1214    FilterInfo("lanczos6",         &lanczos6Filter,        Lanczos6Support),
1215    FilterInfo("lanczos12",        &lanczos12Filter,       Lanczos12Support),
1216    FilterInfo("kaiser",           &kaiserFilter,          KaiserSupport),
1217    FilterInfo("gaussian",         &gaussianFilter,        GaussianSupport),
1218    FilterInfo("catmullrom",       &catmullRomFilter,      CatmullRomSupport),
1219    FilterInfo("quadratic_interp", &quadraticInterpFilter, QuadraticSupport),
1220    FilterInfo("quadratic_approx", &quadraticApproxFilter, QuadraticSupport),
1221    FilterInfo("quadratic_mix",    &quadraticMixFilter,    QuadraticSupport),
1222 ];
1223 
1224 enum NumFilters = cast(int)gFilters.length;
1225 
1226 
1227 bool rsmStringEqu (const(char)[] s0, const(char)[] s1) {
1228   for (;;) {
1229     if (s0.length && (s0.ptr[0] <= ' ' || s0.ptr[0] == '_')) { s0 = s0[1..$]; continue; }
1230     if (s1.length && (s1.ptr[0] <= ' ' || s1.ptr[0] == '_')) { s1 = s1[1..$]; continue; }
1231     if (s0.length == 0) {
1232       while (s1.length && (s1.ptr[0] <= ' ' || s1.ptr[0] == '_')) s1 = s1[1..$];
1233       return (s1.length == 0);
1234     }
1235     if (s1.length == 0) {
1236       while (s0.length && (s0.ptr[0] <= ' ' || s0.ptr[0] == '_')) s0 = s0[1..$];
1237       return (s0.length == 0);
1238     }
1239     assert(s0.length && s1.length);
1240     char c0 = s0.ptr[0];
1241     char c1 = s1.ptr[0];
1242     if (c0 >= 'A' && c0 <= 'Z') c0 += 32; // poor man's tolower
1243     if (c1 >= 'A' && c1 <= 'Z') c1 += 32; // poor man's tolower
1244     if (c0 != c1) return false;
1245     s0 = s0[1..$];
1246     s1 = s1[1..$];
1247   }
1248 }
1249 
1250 
1251 int resamplerFindFilterInternal (const(char)[] name) {
1252   if (name.length) {
1253     foreach (immutable idx, const ref fi; gFilters[]) if (rsmStringEqu(name, fi.name)) return cast(int)idx;
1254   }
1255   return -1;
1256 }