Basic Image AlgorithmS Library 2.8.0

Rescale.cpp

00001 /* This file is part of the BIAS library (Basic ImageAlgorithmS).
00002 
00003    Copyright (C) 2003-2009    (see file CONTACT for details)
00004    Multimediale Systeme der Informationsverarbeitung
00005    Institut fuer Informatik
00006    Christian-Albrechts-Universitaet Kiel
00007 
00008    BIAS is free software; you can redistribute it and/or modify
00009    it under the terms of the GNU Lesser General Public License as published by
00010    the Free Software Foundation; either version 2.1 of the License, or
00011    (at your option) any later version.
00012 
00013    BIAS is distributed in the hope that it will be useful,
00014    but WITHOUT ANY WARRANTY; without even the implied warranty of
00015    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
00016    GNU Lesser General Public License for more details.
00017 
00018    You should have received a copy of the GNU Lesser General Public License
00019    along with BIAS; if not, write to the Free Software
00020    Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA*/
00021 #include "Rescale.hh"
00022 
00023 #ifdef BIAS_DEBUG
00024 #include <Base/Image/ImageIO.hh>
00025 #endif
00026 
00027 #include <Base/Common/BIASpragma.hh>
00028 #include <Base/Common/CompareFloatingPoint.hh>
00029 #include <Filter/Mean.hh>
00030 #include <Filter/Gauss.hh>
00031 #include <Filter/Binomial.hh>
00032 #include <Filter/GaussThreshold.hh>
00033 
00034 using namespace BIAS;
00035 using namespace std;
00036 
00037 namespace BIAS {
00038 
00039 //////////////////////////////////////////////////////////////////////////
00040 //                 implementation
00041 //////////////////////////////////////////////////////////////////////////
00042 
00043 template <class InputImageType, class OutputImageType>
00044 Rescale<InputImageType, OutputImageType>::
00045 Rescale()
00046   : FilterBase<InputImageType, OutputImageType>()
00047 {
00048   _RescaleMeanSigmaFactor=0.9;
00049   _RescaleFactor= 2.0;
00050   Mean<InputImageType, OutputImageType> *pMean =
00051     new Mean<InputImageType, OutputImageType>;
00052   _lowpass = pMean;
00053   // set offset to 0.5 if we have an even mean filter:
00054   _LastPositionOffset = 0.5 * double((pMean->GetSize()+1)%2);
00055 }
00056 
00057 template <class InputImageType, class OutputImageType>
00058 Rescale<InputImageType, OutputImageType>::
00059 Rescale(const Rescale<InputImageType, OutputImageType>& other)
00060 {
00061    _lowpass = NULL;
00062    (*this) = other;
00063 }
00064 
00065 template <class InputStorageType, class OutputStorageType>
00066 Rescale<InputStorageType, OutputStorageType>&
00067 Rescale<InputStorageType, OutputStorageType>::
00068 operator=(const Rescale<InputStorageType, OutputStorageType>& other) {
00069   _RescaleFactor =  other._RescaleFactor;
00070   _RescaleMeanSigmaFactor = other._RescaleMeanSigmaFactor;
00071   _LastPositionOffset = other._LastPositionOffset;
00072   if (_lowpass) delete _lowpass;
00073   if (other._lowpass)
00074     _lowpass = other._lowpass->Clone();
00075   else
00076     _lowpass = NULL;
00077   return (*this);
00078 }
00079 
00080 
00081 template <class InputImageType, class OutputImageType>
00082 Rescale<InputImageType, OutputImageType>::~Rescale()
00083 {
00084   if (_lowpass) delete _lowpass;
00085 }
00086 
00087 template <class InputImageType, class OutputImageType>
00088 int Rescale<InputImageType, OutputImageType>::
00089 Filter(const Image<InputImageType>& src, Image<OutputImageType>& dst)
00090 {
00091   BIASDOUT(D_FILTERBASE_CALLSTACK, "Rescale::Filter\n");
00092   int res=0;
00093   if (_RescaleFactor<1.0){
00094     res=Upsample(src, dst);
00095   } else if (_RescaleFactor>1.0){
00096     res=Downsample(src, dst);
00097   } else {
00098     dst=src;
00099     res=0;
00100     _LastPositionOffset=0.0;
00101   }
00102   return res;
00103 }
00104 
00105 //////////////////////////////////////////////////////////////////////////
00106 //               downsampling functions
00107 //////////////////////////////////////////////////////////////////////////
00108 
00109 template <class InputImageType, class OutputImageType>
00110 int Rescale<InputImageType, OutputImageType>::
00111 Downsample(const Image<InputImageType>& src, Image<OutputImageType>& dst)
00112 {
00113   BIASDOUT(D_FILTERBASE_CALLSTACK, "Rescale::Downsample(src, dst)\n");
00114   double pow2 = log(_RescaleFactor)/log(2.0);
00115 
00116   if (dynamic_cast<Mean<InputImageType, OutputImageType>*>(_lowpass)){
00117     unsigned newwidth, newheight;
00118     if (fabs(_RescaleFactor-2.0)<=DBL_EPSILON){ // factor is 2
00119       newwidth=src.GetWidth()>>1;
00120       newheight=src.GetHeight()>>1;
00121       if (dst.GetWidth()!=newwidth || dst.GetHeight()!=newheight
00122           || dst.GetChannelCount()!=src.GetChannelCount()){
00123         if (!dst.IsEmpty()) dst.Release();
00124         dst.Init(newwidth, newheight, src.GetChannelCount(),
00125                  src.GetStorageType(), src.IsInterleaved() );
00126       }
00127       return DownsampleBy2(src, dst);
00128     } else if (fabs(_RescaleFactor-4.0)<=DBL_EPSILON ){ // factor is 4
00129       newwidth=src.GetWidth()>>2;
00130       newheight=src.GetHeight()>>2;
00131       if (dst.GetWidth()!=newwidth || dst.GetHeight()!=newheight
00132           || dst.GetChannelCount()!=src.GetChannelCount()){
00133         if (!dst.IsEmpty()) dst.Release();
00134         dst.Init(newwidth, newheight, src.GetChannelCount(),
00135                  src.GetStorageType(), src.IsInterleaved());
00136       }
00137       return DownsampleBy4(src, dst);
00138     } else if ( pow2 - floor(pow2) <DBL_EPSILON ) { // factor is a power of 2
00139       return DownsampleBPoT(src,dst);
00140     }
00141   }
00142   return  Downsample(src, dst, _RescaleFactor);
00143 }
00144 
00145 
00146 template <class InputImageType, class OutputImageType>
00147 int Rescale<InputImageType, OutputImageType>::
00148 Downsample(const Image<InputImageType>& src, Image<OutputImageType>& dst,
00149            const double factor)
00150 {
00151   const double old_width = (double)src.GetWidth();
00152   const double old_height = (double)src.GetHeight();
00153   const unsigned newwidth = (unsigned)floor(old_width/factor);
00154   const unsigned newheight = (unsigned)floor(old_height/factor);
00155 
00156   Image<OutputImageType> meanim;
00157   int res = 0;
00158   if ((res=_ApplyMeanFilter(src, meanim, factor))!=0){
00159     BIASERR("error applying mean filter");
00160     return res;
00161   }
00162   int tlx, tly, brx, bry;
00163   meanim.GetROI()->GetCorners(tlx, tly, brx, bry);
00164 
00165   tlx=(int)ceil((double)tlx/factor);
00166   tly=(int)ceil((double)tly/factor);
00167   brx=(int)floor((double)brx/factor);
00168   bry=(int)floor((double)bry/factor);
00169 
00170   dst.Init(newwidth, newheight, src.GetChannelCount(), src.GetStorageType(),
00171            src.IsInterleaved());
00172 
00173   dst.GetROI()->SetCorners(tlx, tly, brx, bry);
00174   return _FillInterpolated(meanim, factor, dst);
00175 }
00176 
00177 
00178 
00179 template <class InputImageType, class OutputImageType>
00180 int Rescale<InputImageType, OutputImageType>::
00181 Downsample(const Image<InputImageType>& src, Image<OutputImageType>& dst,
00182            const unsigned newwidth, const unsigned newheight)
00183 { 
00184 
00185   BIASDOUT(D_FILTERBASE_CALLSTACK, "Rescale::Downsample(src, dst, width, height)\n");
00186 
00187   if (src.GetWidth() == newwidth && src.GetHeight() == newheight) {
00188       dst = src;
00189       return 0;
00190   }
00191 
00192 
00193   const double hfactor = (double)src.GetWidth() / (double)newwidth;
00194   const double vfactor = (double)src.GetHeight() / (double)newheight;
00195   const double factor = (vfactor>hfactor)?vfactor:hfactor;
00196 
00197   Image<OutputImageType> meanim;
00198   if (_ApplyMeanFilter(src, meanim, factor)!=0){
00199     BIASERR("error applying mean filter");
00200   }
00201 
00202   Mean<InputImageType, OutputImageType> *pMean =
00203     dynamic_cast<Mean<InputImageType, OutputImageType> *>(_lowpass);
00204   if (pMean) {
00205     // set offset to 0.5 if we have an even mean filter:
00206     _LastPositionOffset = 0.5 * double((pMean->GetSize()+1)%2);
00207   } else {
00208     _LastPositionOffset = 0.0;
00209   }
00210 
00211   int tlx, tly, brx, bry;
00212   meanim.GetROI()->GetCorners(tlx, tly, brx, bry);
00213 
00214   tlx=(int)ceil((double)tlx/hfactor);
00215   tly=(int)ceil((double)tly/vfactor);
00216   brx=(int)floor((double)brx/hfactor);
00217   bry=(int)floor((double)bry/vfactor);
00218 
00219   if (dst.GetWidth()!=newwidth || dst.GetHeight()!=newheight){
00220     if (!dst.IsEmpty()) dst.Release();
00221     dst.Init(newwidth, newheight, src.GetChannelCount(), src.GetStorageType(),
00222              src.IsInterleaved());
00223   }
00224   dst.GetROI()->SetCorners(tlx, tly, brx, bry);
00225 
00226   int res=0;
00227   if (src.GetChannelCount()==1 && src.GetColorModel()==ImageBase::CM_Grey){
00228     res=_FillInterpolatedGrey(meanim, dst);
00229   } else {
00230     res=_FillInterpolatedColor(meanim, dst);
00231   }
00232 
00233   return res;
00234 }
00235 
00236 /*
00237 template <class InputImageType, class OutputImageType>
00238 int Rescale<InputImageType, OutputImageType>::
00239 DownsampleROI(const Image<InputImageType>& src, Image<OutputImageType>& dst,
00240               const unsigned newwidth, const unsigned newheight)
00241 {
00242   BIASDOUT(D_FILTERBASE_CALLSTACK, "DownsampleROI(width, height)\n");
00243   unsigned int src_width = src.GetWidth();
00244   unsigned int src_height = src.GetHeight();
00245 
00246 
00247   // The ROI boundaries
00248   unsigned int ul_x;
00249   unsigned int ul_y;
00250   unsigned int lr_x;
00251   unsigned int lr_y;
00252   src.GetROI()->GetCorners(ul_x, ul_y, lr_x, lr_y);
00253 
00254 #ifdef BIAS_DEBUG
00255   if(ul_x==0 && ul_y==0 && lr_x==src_width && lr_y==src_height){
00256     BIASERR("No ROI defined in source image. Nothing to downsample.");
00257     return -1;
00258   }
00259 #endif
00260 
00261   unsigned int roi_width = lr_x - ul_x + 1;
00262 
00263   // So what is the factor anyway? (Based on newwidth)
00264   float factor = 1.0 / ((float) newwidth / (float) roi_width);
00265 
00266 #ifdef BIAS_DEBUG
00267   unsigned int roi_height = lr_y - ul_y + 1;
00268   // Scaling factor in x and y direction equal?
00269   if (  roi_width * newheight != roi_height * newwidth) {
00270     BIASERR("Downsample to new width AND height for ROI width and height must "
00271             <<"be of equal factor in both directions. ");
00272     // ... or adapt the mean filter!
00273     return -1;
00274   }
00275 
00276   if (factor<1.0){
00277     BIASERR("New width and height would suggest an upsampling. Nothing done "
00278             <<factor);
00279     return -1;
00280   }
00281 #endif
00282 
00283   unsigned int halfsize = (factor <= 2.0) ? 1 :
00284     (unsigned int) ceil((factor-1.0)/2.0);
00285 
00286   // Check the ROI boundaries
00287   if (ul_x < halfsize ||
00288       ul_y < halfsize ||
00289       lr_x > src_width - 1 - halfsize ||
00290       lr_y > src_height - 1 - halfsize
00291       )
00292     {
00293       BIASERR("ROI is to near to the edge of the source image.");
00294       return -1;
00295     }
00296 
00297   // Create an image containing the ROI plus edge
00298   Image<InputImageType> roiwithedge, srcroi=src;
00299   srcroi.GetROI()->SetCorners(ul_x - halfsize, ul_y - halfsize, lr_x + halfsize,
00300                 lr_y + halfsize);
00301   srcroi.GetCopyOfROI2(roiwithedge);
00302 
00303   Image<OutputImageType> meanim;
00304   if (_ApplyMeanFilter(roiwithedge, meanim, factor)!=0){
00305     BIASERR("error applying mean filter");
00306   }
00307 
00308   meanim.GetROI()->SetCorners(halfsize, halfsize,
00309     meanim.GetWidth() - 1 - halfsize,
00310     meanim.GetHeight() - 1 - halfsize);
00311 
00312   Image<OutputImageType> filtered_image;
00313   meanim.GetCopyOfROI(filtered_image);
00314 
00315   return Upsample(filtered_image, dst, newwidth, newheight);
00316 }
00317 */
00318 ////////////////////////////////////////////////////////////////////////////
00319 
00320 template <class InputImageType, class OutputImageType>
00321 int Rescale<InputImageType, OutputImageType>::
00322 DownsampleBy2(const Image<InputImageType>& src, Image<OutputImageType>& dst)
00323 {
00324   BIASDOUT(D_FILTERBASE_CALLSTACK, "Rescale::DownsampleBy2(src, dst)\n");
00325   int res=0;
00326   unsigned nwidth=src.GetWidth()>>1;
00327   unsigned nheight=src.GetHeight()>>1;
00328   if (dst.GetWidth()!=nwidth || dst.GetHeight()!=nheight){
00329     if (!dst.IsEmpty()) dst.Release();
00330     dst.Init(nwidth, nheight, src.GetChannelCount(), src.GetStorageType(),
00331              src.IsInterleaved());
00332   }
00333   if (src.GetColorModel()==ImageBase::CM_Grey){
00334     res = DownsampleBy2Grey(src, dst);
00335   } else {
00336     if (src.IsInterleaved() && src.GetChannelCount()==3){
00337       res = DownsampleBy2Color(src, dst);
00338     } else {
00339       // this is slow
00340       int newwidth=src.GetWidth()>>1;
00341       int newheight=src.GetHeight()>>1;
00342       res = Downsample(src, dst, newwidth, newheight);
00343     }
00344   }
00345   // see DownsampleBy2 for details
00346   _LastPositionOffset=0.5;
00347   return res;
00348 }
00349 
00350 template <class InputImageType, class OutputImageType>
00351 int Rescale<InputImageType, OutputImageType>::
00352 DownsampleBy2Grey(const Image<InputImageType>& Source,
00353                   Image<OutputImageType>& Destination)
00354 {
00355   BIASDOUT(D_FILTERBASE_CALLSTACK, "Rescale::DownsampleBy2Grey\n");
00356   unsigned minx, miny, maxx, maxy, dx, dy, mdx, mdy;
00357   Source.GetROI()->GetCorners(minx, miny, maxx, maxy);
00358 
00359   // make sure the image is not shifted
00360   if (minx%2!=0) minx++;
00361   if (maxx%2!=0) maxx--;
00362   if (miny%2!=0) miny++;
00363   if (maxy%2!=0) maxy--;
00364 
00365   dx=minx>>1;
00366   dy=miny>>1;
00367   mdx=maxx>>1;
00368   mdy=(maxy>>1)-1;
00369 
00370   Destination.GetROI()->SetCorners(dx, dy, mdx, mdy+1);
00371 
00372   register double sum;
00373   register const InputImageType *srcu, *srcl;
00374   register OutputImageType *dst, *dstend, *dstlend;
00375   register int srcwidth=Source.GetWidth(), dstwidth=Destination.GetWidth(),
00376     dststep;
00377   dst = Destination.GetImageData();
00378 
00379   srcl = srcu = Source.GetImageData()+minx+miny*srcwidth;
00380   srcl += srcwidth;
00381   dst = Destination.GetImageData() + dy*dstwidth;
00382   dstlend = dst + mdx;
00383   dstend = Destination.GetImageData() + mdx + mdy*dstwidth;
00384   dst+=dx;
00385 
00386   srcwidth+=(srcwidth-maxx+minx);
00387   dststep=dstwidth-mdx+dx;
00388 
00389   if (srcwidth%2!=0){
00390     srcwidth+=1;
00391   }
00392   while (dst < dstend){
00393     while (dst < dstlend){
00394       sum = 0.0;
00395       sum += (double)*srcu++;
00396       sum += (double)*srcu++;
00397       sum += (double)*srcl++;
00398       sum += (double)*srcl++;
00399       *dst++ = (OutputImageType)(sum / 4.0);
00400     }
00401     srcu += srcwidth;
00402     srcl += srcwidth;
00403     dst += dststep;
00404     dstlend += dstwidth;
00405   }
00406   // see DownsampleBy2 for details
00407   _LastPositionOffset=0.5;
00408   return 0;
00409 }
00410 
00411 
00412 template <>
00413 int Rescale<unsigned char, unsigned char>::
00414 DownsampleBy2Grey(const Image<unsigned char>& Source,
00415                   Image<unsigned char>& Destination) {
00416   BIASDOUT(D_FILTERBASE_CALLSTACK, "Rescale::DownsampleBy2Grey\n");
00417   register unsigned  sum;
00418   unsigned minx, miny, maxx, maxy, dx, dy, mdx, mdy;
00419   Source.GetROI()->GetCorners(minx, miny, maxx, maxy);
00420 
00421   // make sure the image is not shifted
00422   if (minx%2!=0) minx++;
00423   if (maxx%2!=0) maxx--;
00424   if (miny%2!=0) miny++;
00425   if (maxy%2!=0) maxy--;
00426 
00427   dx=minx>>1;
00428   dy=miny>>1;
00429   mdx=maxx>>1;
00430   mdy=(maxy>>1)-1;
00431 
00432   Destination.GetROI()->SetCorners(dx, dy, mdx, mdy+1);
00433 
00434   register const unsigned char *srcu, *srcl;
00435   register unsigned char *dst, *dstend, *dstlend;
00436   register int srcwidth=Source.GetWidth(), dstwidth=Destination.GetWidth(),
00437     dststep;
00438 
00439   srcl = srcu = Source.GetImageData()+minx+miny*srcwidth;
00440   srcl += srcwidth;
00441   dst = Destination.GetImageData() + dy*dstwidth;
00442   dstlend = dst + mdx;
00443   dstend = Destination.GetImageData() + mdx + mdy*dstwidth;
00444   dst+=dx;
00445 
00446   srcwidth+=(srcwidth-maxx+minx);
00447   dststep=dstwidth-mdx+dx;
00448 
00449   if (srcwidth%2!=0){
00450     srcwidth+=1;
00451   }
00452   while (dst < dstend){
00453     while (dst < dstlend){
00454       sum = 0;
00455       sum += *srcu++;
00456       sum += *srcu++;
00457       sum += *srcl++;
00458       sum += *srcl++;
00459       *dst++ = sum >> 2;
00460     }
00461     srcu += srcwidth;
00462     srcl += srcwidth;
00463     dst += dststep;
00464     dstlend += dstwidth;
00465   }
00466   // see DownsampleBy2 for details
00467   _LastPositionOffset=0.5;
00468 
00469   return 0;
00470 }
00471 
00472 
00473 #ifdef BUILD_IMAGE_USHORT
00474 template <>
00475 int Rescale<unsigned short, unsigned short>::
00476 DownsampleBy2Grey(const Image<unsigned short>& Source,
00477                       Image<unsigned short>& Destination)
00478 {
00479   BIASDOUT(D_FILTERBASE_CALLSTACK, "Rescale::DownsampleBy2Grey\n");
00480   register unsigned  sum;
00481   unsigned minx, miny, maxx, maxy, dx, dy, mdx, mdy;
00482   Source.GetROI()->GetCorners(minx, miny, maxx, maxy);
00483 
00484   // make sure the image is not shifted
00485   if (minx%2!=0) minx++;
00486   if (maxx%2!=0) maxx--;
00487   if (miny%2!=0) miny++;
00488   if (maxy%2!=0) maxy--;
00489 
00490   dx=minx>>1;
00491   dy=miny>>1;
00492   mdx=maxx>>1;
00493   mdy=(maxy>>1)-1;
00494 
00495   Destination.GetROI()->SetCorners(dx, dy, mdx, mdy+1);
00496 
00497   register const unsigned short *srcu, *srcl;
00498   register unsigned short *dst, *dstend, *dstlend;
00499   register int srcwidth=Source.GetWidth(), dstwidth=Destination.GetWidth(),
00500     dststep;
00501 
00502   srcl = srcu = Source.GetImageData()+minx+miny*srcwidth;
00503   srcl += srcwidth;
00504   dst = Destination.GetImageData() + dy*dstwidth;
00505   dstlend = dst + mdx;
00506   dstend = Destination.GetImageData() + mdx + mdy*dstwidth;
00507   dst+=dx;
00508 
00509   srcwidth+=(srcwidth-maxx+minx);
00510   dststep=dstwidth-mdx+dx;
00511 
00512   if (srcwidth%2!=0){
00513     srcwidth+=1;
00514   }
00515   while (dst < dstend){
00516     while (dst < dstlend){
00517       sum = 0;
00518       sum += *srcu++;
00519       sum += *srcu++;
00520       sum += *srcl++;
00521       sum += *srcl++;
00522       *dst++ = sum >> 2;
00523     }
00524     srcu += srcwidth;
00525     srcl += srcwidth;
00526     dst += dststep;
00527     dstlend += dstwidth;
00528   }
00529   // see DownsampleBy2 for details
00530   _LastPositionOffset=0.5;
00531 
00532   return 0;
00533 }
00534 #endif
00535 
00536 template <class InputImageType, class OutputImageType>
00537 int Rescale<InputImageType, OutputImageType>::
00538 DownsampleBy2Color(const Image<InputImageType>& Source,
00539                    Image<OutputImageType>& Destination)
00540 {
00541   BIASDOUT(D_FILTERBASE_CALLSTACK, "Rescale::DownsampleBy2Color\n");
00542   BIASASSERT(Source.IsPlanar()==Destination.IsPlanar());
00543 #ifdef BIAS_DEBUG
00544   if (this->DebugLevelIsSet(D_RESCALE_WRITE_IMG)) {
00545     ImageBase im=Source;
00546     //ImageIO::Save("filtered-source.mip", im);
00547     ImageIO::Save("filtered-source.mip", im);
00548   }
00549   unsigned tlx, tly, brx, bry;
00550   Source.GetROI()->GetCorners(tlx, tly, brx, bry);
00551   if (tlx!=0 || tly!=0 || brx!=Source.GetWidth() || bry!=Source.GetHeight()){
00552     BIASERR("ROI not yet supported");
00553     //BIASASSERT(false);
00554     BIASABORT;
00555   }
00556 #endif
00557 
00558 
00559 
00560   register double sum1, sum2, sum3;
00561   register const InputImageType *srcu, *srcl;
00562   register OutputImageType   *dst, *dstend, *dstlend;
00563   register int srcwidth=Source.GetWidth()*3, dstwidth=Destination.GetWidth()*3;
00564   dst = Destination.GetImageData();
00565   srcl = srcu = Source.GetImageData();
00566   srcl += srcwidth;
00567   dstend = dst + Destination.GetPixelCount() * 3;
00568   dstlend = dst + dstwidth * 3;
00569 
00570   if (srcwidth%2!=0){
00571     srcwidth+=3;
00572   }
00573   while (dst < dstend){
00574     while (dst < dstlend){
00575       sum1 = sum2 = sum3 = 0.0;
00576       sum1 += (double)*srcu++;
00577       sum2 += (double)*srcu++;
00578       sum3 += (double)*srcu++;
00579       sum1 += (double)*srcu++;
00580       sum2 += (double)*srcu++;
00581       sum3 += (double)*srcu++;
00582       sum1 += (double)*srcl++;
00583       sum2 += (double)*srcl++;
00584       sum3 += (double)*srcl++;
00585       sum1 += (double)*srcl++;
00586       sum2 += (double)*srcl++;
00587       sum3 += (double)*srcl++;
00588       *dst++ = (OutputImageType)(sum1 / 4.0);
00589       *dst++ = (OutputImageType)(sum2 / 4.0);
00590       *dst++ = (OutputImageType)(sum3 / 4.0);
00591     }
00592     srcu += srcwidth;
00593     srcl += srcwidth;
00594     dstlend += dstwidth;
00595   }
00596   // see DownsampleBy2 for details
00597   _LastPositionOffset=0.5;
00598   return 0;
00599 }
00600 
00601 template <>
00602 int Rescale<unsigned char, unsigned char>::
00603 DownsampleBy2Color(const Image<unsigned char>& Source,
00604                    Image<unsigned char>& Destination)
00605 {
00606   BIASDOUT(D_FILTERBASE_CALLSTACK, "Rescale::DownsampleBy2Color\n");
00607   BIASASSERT(Source.IsPlanar()==Destination.IsPlanar());
00608 #ifdef BIAS_DEBUG
00609   if (this->DebugLevelIsSet(D_RESCALE_WRITE_IMG)) {
00610     ImageBase im=Source;
00611     //ImageIO::Save("filtered-source.mip", im);
00612     ImageIO::Save("filtered-source.mip", im);
00613   }
00614   unsigned tlx, tly, brx, bry;
00615   Source.GetROI()->GetCorners(tlx, tly, brx, bry);
00616   if (tlx!=0 || tly!=0 || brx!=Source.GetWidth() || bry!=Source.GetHeight()){
00617     BIASERR("ROI not yet supported");
00618     //BIASASSERT(false);
00619     BIASABORT;
00620   }
00621 #endif
00622   register unsigned int sum1, sum2, sum3;
00623   register const unsigned char *srcu, *srcl;
00624   register unsigned char *dst, *dstend, *dstlend;
00625   register int srcwidth=Source.GetWidth() * 3,
00626     dstwidth=Destination.GetWidth() *3;
00627   dst = Destination.GetImageData();
00628   srcl = srcu = Source.GetImageData();
00629   srcl += srcwidth;
00630   dstend = dst + Destination.GetPixelCount() * 3;
00631   dstlend = dst + dstwidth;
00632 
00633   if (srcwidth%2!=0){
00634     srcwidth+=3;
00635   }
00636   while (dst < dstend){
00637     while (dst < dstlend){
00638       sum1 = sum2 = sum3 = 0;
00639       sum1 += *srcu++;
00640       sum2 += *srcu++;
00641       sum3 += *srcu++;
00642       sum1 += *srcu++;
00643       sum2 += *srcu++;
00644       sum3 += *srcu++;
00645       sum1 += *srcl++;
00646       sum2 += *srcl++;
00647       sum3 += *srcl++;
00648       sum1 += *srcl++;
00649       sum2 += *srcl++;
00650       sum3 += *srcl++;
00651       *dst++ = sum1 >> 2;
00652       *dst++ = sum2 >> 2;
00653       *dst++ = sum3 >> 2;
00654     }
00655     srcu += srcwidth;
00656     srcl += srcwidth;
00657     dstlend += dstwidth;
00658   }
00659 
00660 //   for (unsigned int y=0;y<Destination.GetHeight(); y++){
00661 //     for (unsigned int x=0;x<Destination.GetWidth(); x++){
00662 //       for (unsigned int c=0;c<3;c++){
00663 //  sum1 = srcu[y *
00664 //       }
00665 //     }
00666 //   }
00667 
00668   // see DownsampleBy2 for details
00669   _LastPositionOffset=0.5;
00670 
00671   return 0;
00672 }
00673 
00674 #ifdef BUILD_IMAGE_USHORT
00675 template <>
00676 int Rescale<unsigned short, unsigned short>::
00677 DownsampleBy2Color(const Image<unsigned short>& Source,
00678                    Image<unsigned short>& Destination)
00679 {
00680   BIASDOUT(D_FILTERBASE_CALLSTACK, "Rescale::DownsampleBy2Color\n");
00681 #ifdef BIAS_DEBUG
00682   if (this->DebugLevelIsSet(D_RESCALE_WRITE_IMG)) {
00683     ImageBase im=Source;
00684     //ImageIO::Save("filtered-source.mip", im);
00685     ImageIO::Save("filtered-source.mip", im);
00686   }
00687   unsigned tlx, tly, brx, bry;
00688   Source.GetROI()->GetCorners(tlx, tly, brx, bry);
00689   if (tlx!=0 || tly!=0 || brx!=Source.GetWidth() || bry!=Source.GetHeight()){
00690     BIASERR("ROI not yet supported");
00691     BIASASSERT(false);
00692   }
00693 #endif
00694   register unsigned int sum1, sum2, sum3;
00695   register const unsigned short *srcu, *srcl;
00696   register unsigned short *dst, *dstend, *dstlend;
00697   register int srcwidth=Source.GetWidth() * 3,
00698     dstwidth=Destination.GetWidth() *3;
00699   dst = Destination.GetImageData();
00700   srcl = srcu = Source.GetImageData();
00701   srcl += srcwidth;
00702   dstend = dst + Destination.GetPixelCount() * 3;
00703   dstlend = dst + dstwidth;
00704 
00705   if (srcwidth%2!=0){
00706     srcwidth+=3;
00707   }
00708   while (dst < dstend){
00709     while (dst < dstlend){
00710       sum1 = sum2 = sum3 = 0;
00711       sum1 += *srcu++;
00712       sum2 += *srcu++;
00713       sum3 += *srcu++;
00714       sum1 += *srcu++;
00715       sum2 += *srcu++;
00716       sum3 += *srcu++;
00717       sum1 += *srcl++;
00718       sum2 += *srcl++;
00719       sum3 += *srcl++;
00720       sum1 += *srcl++;
00721       sum2 += *srcl++;
00722       sum3 += *srcl++;
00723       *dst++ = sum1 >> 2;
00724       *dst++ = sum2 >> 2;
00725       *dst++ = sum3 >> 2;
00726     }
00727     srcu += srcwidth;
00728     srcl += srcwidth;
00729     dstlend += dstwidth;
00730   }
00731   // see DownsampleBy2 for details
00732   _LastPositionOffset=0.5;
00733   return 0;
00734 }
00735 #endif
00736 
00737 /////////////////////////////////////////////////////////////////////////////
00738 
00739 template <class InputImageType, class OutputImageType>
00740 int Rescale<InputImageType, OutputImageType>::
00741 DownsampleBy4(const Image<InputImageType>& src, Image<OutputImageType>& dst)
00742 {
00743   BIASDOUT(D_FILTERBASE_CALLSTACK, "Rescale::DownsampleBy4\n");
00744   int res=0;
00745   unsigned nwidth=(int)src.GetWidth()>>2;
00746   unsigned nheight=(int)src.GetHeight()>>2;
00747   if (dst.GetWidth()!=nwidth || dst.GetHeight()!=nheight){
00748     if (!dst.IsEmpty()) dst.Release();
00749     dst.Init(nwidth, nheight, src.GetChannelCount(), src.GetStorageType(),
00750              src.IsInterleaved());
00751   }
00752   if (src.GetColorModel()==ImageBase::CM_Grey){
00753     res = DownsampleBy4Grey(src, dst);
00754   } else {
00755     if (src.IsInterleaved() && src.GetChannelCount()==3){
00756       double oldFactor = GetFactor();
00757       SetFactor(4.0);
00758       res = DownsampleBPoT(src, dst);
00759       SetFactor(oldFactor);
00760     } else {
00761       // this is slow
00762       int newwidth=src.GetWidth()>>2;
00763       int newheight=src.GetHeight()>>2;
00764       res = Downsample(src, dst, newwidth, newheight);
00765     }
00766   }
00767   // see DownsampleBy2 for details
00768   _LastPositionOffset=1.5;
00769   return res;
00770 }
00771 
00772 template <class InputImageType, class OutputImageType>
00773 int Rescale<InputImageType, OutputImageType>::
00774 DownsampleBy4Grey(const Image<InputImageType>& Source,
00775                   Image<OutputImageType>& Destination)
00776 {
00777   BIASDOUT(D_FILTERBASE_CALLSTACK, "Rescale::DownsampleBy4Grey\n");
00778   BIASASSERT(Source.GetChannelCount()==1 &&
00779              Source.GetColorModel()==ImageBase::CM_Grey);
00780   BIASASSERT(Destination.GetChannelCount()==1 &&
00781              Destination.GetColorModel()==ImageBase::CM_Grey);
00782   BIASASSERT(Source.GetWidth()>>2==Destination.GetWidth());
00783   BIASASSERT(Source.GetHeight()>>2==Destination.GetHeight());
00784 
00785   /// \todo warning kerneltype sum needed here
00786 
00787   double  sum;
00788   unsigned minx, miny, maxx, maxy, dx, dy, mdx, mdy;
00789   Source.GetROI()->GetCorners(minx, miny, maxx, maxy);
00790 
00791   // make sure the image is not shifted
00792   if (minx%4!=0) minx+=(4-(minx%4));
00793   if (maxx%4!=0) maxx-=(maxx%4);
00794   if (miny%4!=0) miny+=(4-(miny%4));
00795   if (maxy%4!=0) maxy-=(maxy%4);
00796 
00797   dx=minx>>2;
00798   dy=miny>>2;
00799   mdx=maxx>>2;
00800   mdy=maxy>>2;
00801 
00802   Destination.GetROI()->SetCorners(dx, dy, mdx, mdy);
00803 
00804   register const InputImageType *src1, *src2, *src3, *src4;
00805   register OutputImageType *dst, *dstend, *dstlend;
00806   register int srcwidth=(int)Source.GetWidth();
00807   register const int dstwidth=(int)Destination.GetWidth();
00808   register int dststep;
00809 
00810   src1 = Source.GetImageData()+minx+miny*srcwidth;
00811   src2 = src1 + srcwidth;
00812   src3 = src2 + srcwidth;
00813   src4 = src3 + srcwidth;
00814 
00815   dst = Destination.GetImageData() + dy*dstwidth;
00816   dstlend = dst + mdx;
00817   dst+=dx;
00818   dstend = Destination.GetImageData() + mdx + (mdy-1)*dstwidth;
00819 
00820   srcwidth+=(3*srcwidth-maxx+minx);
00821   dststep=dstwidth-mdx+dx;
00822 
00823   if (srcwidth%2!=0){
00824     srcwidth+=1;
00825   }
00826 
00827   /// \todo warning kerneltype div needed here
00828 
00829   const double div=(1.0/16.0);
00830   while (dst < dstend){
00831     while (dst < dstlend){
00832       sum = 0.0;
00833       sum += *src1++;
00834       sum += *src1++;
00835       sum += *src1++;
00836       sum += *src1++;
00837 
00838       sum += *src2++;
00839       sum += *src2++;
00840       sum += *src2++;
00841       sum += *src2++;
00842 
00843       sum += *src3++;
00844       sum += *src3++;
00845       sum += *src3++;
00846       sum += *src3++;
00847 
00848       sum += *src4++;
00849       sum += *src4++;
00850       sum += *src4++;
00851       sum += *src4++;
00852 
00853       *dst++ = (OutputImageType)(sum * div);
00854     }
00855     src1 += srcwidth;
00856     src2 += srcwidth;
00857     src3 += srcwidth;
00858     src4 += srcwidth;
00859     dst += dststep;
00860     dstlend += dstwidth;
00861   }
00862   // see DownsampleBy2 for details
00863   _LastPositionOffset=1.5;
00864 
00865   return 0;
00866 }
00867 
00868 
00869 template <>
00870 int Rescale<unsigned char, unsigned char>::
00871 DownsampleBy4Grey(const Image<unsigned char>& Source,
00872                   Image<unsigned char>& Destination)
00873 {
00874   BIASDOUT(D_FILTERBASE_CALLSTACK, "Rescale::DownsampleBy4Grey\n");
00875   BIASASSERT(Source.GetChannelCount()==1 &&
00876              Source.GetColorModel()==ImageBase::CM_Grey);
00877   BIASASSERT(Destination.GetChannelCount()==1 &&
00878              Destination.GetColorModel()==ImageBase::CM_Grey);
00879   BIASASSERT(Source.GetWidth()>>2==Destination.GetWidth());
00880   BIASASSERT(Source.GetHeight()>>2==Destination.GetHeight());
00881   register unsigned  sum;
00882   unsigned minx, miny, maxx, maxy, dx, dy, mdx, mdy;
00883   Source.GetROI()->GetCorners(minx, miny, maxx, maxy);
00884 
00885   // make sure the image is not shifted
00886   if (minx%4!=0) minx+=(4-(minx%4));
00887   if (maxx%4!=0) maxx-=(maxx%4);
00888   if (miny%4!=0) miny+=(4-(miny%4));
00889   if (maxy%4!=0) maxy-=(maxy%4);
00890 
00891   dx=minx>>2;
00892   dy=miny>>2;
00893   mdx=maxx>>2;
00894   mdy=maxy>>2;
00895 
00896   Destination.GetROI()->SetCorners(dx, dy, mdx, mdy);
00897 
00898   register const unsigned char *src1, *src2, *src3, *src4;
00899   register unsigned char *dst, *dstend, *dstlend;
00900   register int srcwidth=(int)Source.GetWidth();
00901   register const int dstwidth=(int)Destination.GetWidth();
00902   register int dststep;
00903 
00904   src1 = Source.GetImageData()+minx+miny*srcwidth;
00905   src2 = src1 + srcwidth;
00906   src3 = src2 + srcwidth;
00907   src4 = src3 + srcwidth;
00908 
00909   dst = Destination.GetImageData() + dy*dstwidth;
00910   dstlend = dst + mdx;
00911   dst+=dx;
00912   dstend = Destination.GetImageData() + mdx + (mdy-1)*dstwidth;
00913 
00914   srcwidth+=(3*srcwidth-maxx+minx);
00915   dststep=dstwidth-mdx+dx;
00916 
00917   if (srcwidth%2!=0){
00918     srcwidth+=1;
00919   }
00920   while (dst < dstend){
00921     while (dst < dstlend){
00922       sum = 0;
00923       sum += *src1++;
00924       sum += *src1++;
00925       sum += *src1++;
00926       sum += *src1++;
00927 
00928       sum += *src2++;
00929       sum += *src2++;
00930       sum += *src2++;
00931       sum += *src2++;
00932 
00933       sum += *src3++;
00934       sum += *src3++;
00935       sum += *src3++;
00936       sum += *src3++;
00937 
00938       sum += *src4++;
00939       sum += *src4++;
00940       sum += *src4++;
00941       sum += *src4++;
00942 
00943       *dst++ = sum >> 4;
00944     }
00945     src1 += srcwidth;
00946     src2 += srcwidth;
00947     src3 += srcwidth;
00948     src4 += srcwidth;
00949     dst += dststep;
00950     dstlend += dstwidth;
00951   }
00952   // see DownsampleBy2 for details
00953   _LastPositionOffset=1.5;
00954 
00955   return 0;
00956 }
00957 
00958 #ifdef BUILD_IMAGE_USHORT
00959 template <>
00960 int Rescale<unsigned short, unsigned short>::
00961 DownsampleBy4Grey(const Image<unsigned short>& Source,
00962                       Image<unsigned short>& Destination)
00963 {
00964   BIASDOUT(D_FILTERBASE_CALLSTACK, "Rescale::DownsampleBy4Grey\n");
00965   BIASASSERT(Source.GetChannelCount()==1 &&
00966              Source.GetColorModel()==ImageBase::CM_Grey);
00967   BIASASSERT(Destination.GetChannelCount()==1 &&
00968              Destination.GetColorModel()==ImageBase::CM_Grey);
00969   BIASASSERT(Source.GetWidth()>>2==Destination.GetWidth());
00970   BIASASSERT(Source.GetHeight()>>2==Destination.GetHeight());
00971   register unsigned  sum;
00972   unsigned minx, miny, maxx, maxy, dx, dy, mdx, mdy;
00973   Source.GetROI()->GetCorners(minx, miny, maxx, maxy);
00974 
00975   // make sure the image is not shifted
00976   if (minx%4!=0) minx+=(4-(minx%4));
00977   if (maxx%4!=0) maxx-=(maxx%4);
00978   if (miny%4!=0) miny+=(4-(miny%4));
00979   if (maxy%4!=0) maxy-=(maxy%4);
00980 
00981   dx=minx>>2;
00982   dy=miny>>2;
00983   mdx=maxx>>2;
00984   mdy=maxy>>2;
00985 
00986   Destination.GetROI()->SetCorners(dx, dy, mdx, mdy);
00987 
00988   register const unsigned short *src1, *src2, *src3, *src4;
00989   register unsigned short *dst, *dstend, *dstlend;
00990   register int srcwidth=(int)Source.GetWidth();
00991   register const int dstwidth=(int)Destination.GetWidth();
00992   register int dststep;
00993 
00994   src1 = Source.GetImageData()+minx+miny*srcwidth;
00995   src2 = src1 + srcwidth;
00996   src3 = src2 + srcwidth;
00997   src4 = src3 + srcwidth;
00998 
00999   dst = Destination.GetImageData() + dy*dstwidth;
01000   dstlend = dst + mdx;
01001   dst+=dx;
01002   dstend = Destination.GetImageData() + mdx + (mdy-1)*dstwidth;
01003 
01004   srcwidth+=(3*srcwidth-maxx+minx);
01005   dststep=dstwidth-mdx+dx;
01006 
01007   if (srcwidth%2!=0){
01008     srcwidth+=1;
01009   }
01010   while (dst < dstend){
01011     while (dst < dstlend){
01012       sum = 0;
01013       sum += *src1++;
01014       sum += *src1++;
01015       sum += *src1++;
01016       sum += *src1++;
01017 
01018       sum += *src2++;
01019       sum += *src2++;
01020       sum += *src2++;
01021       sum += *src2++;
01022 
01023       sum += *src3++;
01024       sum += *src3++;
01025       sum += *src3++;
01026       sum += *src3++;
01027 
01028       sum += *src4++;
01029       sum += *src4++;
01030       sum += *src4++;
01031       sum += *src4++;
01032 
01033       *dst++ = sum >> 4;
01034     }
01035     src1 += srcwidth;
01036     src2 += srcwidth;
01037     src3 += srcwidth;
01038     src4 += srcwidth;
01039     dst += dststep;
01040     dstlend += dstwidth;
01041   }
01042   // see DownsampleBy2 for details
01043   _LastPositionOffset=1.5;
01044 
01045   return 0;
01046 }
01047 #endif
01048 
01049 
01050 
01051 
01052 
01053 template <class InputImageType, class OutputImageType>
01054 int Rescale<InputImageType, OutputImageType>::
01055 DownsampleBPoT(const Image<InputImageType>& Source,
01056                Image<OutputImageType>& Destination)
01057 {
01058   BIASDOUT(D_FILTERBASE_CALLSTACK, "Rescale::DownsampleBPoT\n");
01059   int res = 0;
01060   unsigned int channels = Source.GetChannelCount();
01061   int PowerOfTwo = (int)rint(log(_RescaleFactor)/log(2.0));
01062   if (PowerOfTwo == 0){
01063 
01064     // No downsample required, just copy
01065     Destination = Source;
01066     _LastPositionOffset=0.0;
01067 
01068   } else { // downsample
01069 
01070     unsigned int newwidth, newheight;
01071     unsigned int oldwidth,oldheight;
01072 
01073     oldwidth = Source.GetWidth();
01074     oldheight = Source.GetHeight();
01075 
01076     newwidth = (oldwidth >> PowerOfTwo);
01077     newheight = (oldheight >> PowerOfTwo);
01078 
01079     if (Destination.GetWidth() != newwidth
01080         || Destination.GetHeight() != newheight
01081         ||  Destination.GetChannelCount() != channels ){
01082       if ( ! Destination.IsEmpty())
01083         Destination.Release();
01084       Destination.Init(newwidth,newheight,channels,
01085                        Source.GetStorageType(),
01086                        Source.IsInterleaved());
01087     }
01088     OutputImageType *sink;
01089     const InputImageType *src;
01090     bool IdenticalImages = false;
01091 
01092     sink = Destination.GetImageData();
01093     src = Source.GetImageData();
01094 
01095     if ((void*)src==(void*)sink){
01096       IdenticalImages = true;
01097       sink = new OutputImageType[newwidth * newheight];
01098     }
01099 
01100     unsigned min_oldx, min_oldy, max_oldx, max_oldy;
01101     const double div = pow(2.0,PowerOfTwo+PowerOfTwo);
01102     const int sze = pow(2.0,PowerOfTwo);
01103     double sum;
01104     if (Source.IsInterleaved()){
01105 
01106       for (unsigned y=0; y<newheight; y++){
01107         min_oldy = y<<PowerOfTwo;
01108         max_oldy = min_oldy + sze;
01109         for (unsigned x=0; x<newwidth; x++){
01110           min_oldx = x<<PowerOfTwo;
01111           max_oldx = min_oldx + sze;
01112           for (unsigned c=0; c<channels; c++){
01113             sum=0.0;
01114             // now sum up the values in the old image
01115             for (unsigned oy=min_oldy; oy<max_oldy; oy++){
01116               for (unsigned ox=min_oldx; ox<max_oldx; ox++){
01117                 sum += (double) src[oy*oldwidth*channels + ox*channels+c];
01118               }
01119             }
01120             // now set the new value
01121             sink[y*newwidth*channels+x*channels+c] =
01122               Cast<OutputImageType>(sum / div);
01123           }
01124         }
01125       }
01126 
01127     } else {
01128 
01129       BIASERR("Rescale::DownsampleBPoT unfinished for planar images");
01130       BIASABORT;
01131 
01132     }
01133 
01134 
01135     if (IdenticalImages) {
01136       Destination.Release();
01137       Destination.Init(newwidth, newheight, channels);
01138       Destination.ReleaseImageDataPointer();
01139 
01140       Destination.RedirectImageDataPointer(sink);
01141     }
01142     _LastPositionOffset= pow(2.0, PowerOfTwo-1.0)-0.5;
01143   }
01144 
01145   return res;
01146 }
01147 
01148 
01149 
01150 
01151 //////////////////////////////////////////////////////////////////////////
01152 //               upsampling functions
01153 //////////////////////////////////////////////////////////////////////////
01154 
01155 template <class InputImageType, class OutputImageType>
01156 int Rescale<InputImageType, OutputImageType>::
01157 Upsample(const Image<InputImageType>& src, Image<OutputImageType>& dst)
01158 {
01159   BIASDOUT(D_FILTERBASE_CALLSTACK, "Rescale::Upsample(src, dst)\n");
01160 
01161   return Upsample(src, dst,
01162                   (unsigned int)floor((double)src.GetWidth()/
01163                                       _RescaleFactor),
01164                   (unsigned int)floor((double)src.GetHeight()/
01165                                       _RescaleFactor));
01166 }
01167 
01168 template <class InputImageType, class OutputImageType>
01169 int Rescale<InputImageType, OutputImageType>::
01170 Upsample(const Image<InputImageType>& src, Image<OutputImageType>& dst,
01171          const unsigned newwidth, const unsigned newheight)
01172 {
01173   BIASDOUT(D_FILTERBASE_CALLSTACK, "Rescale::Upsample(src, dst, w, h)\n");
01174 #ifdef BIAS_DEBUG
01175   unsigned tlx, tly, brx, bry;
01176   src.GetROI()->GetCorners(tlx, tly, brx, bry);
01177 
01178   if (tlx!=0 || tly!=0 || brx!=src.GetWidth() || bry!=src.GetHeight()){
01179     BIASERR("Upsampling with ROI not yet supported, ignoring.");
01180     //do not abort here!! aborts suck!
01181     //BIASASSERT(false);
01182     //BIASABORT;
01183   }
01184 #endif
01185   BIASASSERT(src.GetWidth()<=newwidth);
01186   BIASASSERT(src.GetHeight()<=newheight);
01187 
01188   int res=0;
01189   float oldx, oldy;
01190   float scalex, scaley;
01191   unsigned int channelcount=src.GetChannelCount();
01192   bool identicalimages=false;
01193   OutputImageType *sink;
01194 
01195   if ((void *)src.GetImageData() == (void *)dst.GetImageData())
01196     identicalimages=true;
01197 
01198   if (!identicalimages){
01199     if (!src.SamePixelAndChannelCount(dst))
01200       if (!dst.IsEmpty()) dst.Release();
01201 
01202     if (dst.IsEmpty())
01203       dst.Init(newwidth, newheight, channelcount,  src.GetStorageType(),
01204                src.IsInterleaved());
01205 
01206     sink = dst.GetImageData();
01207   } else {
01208     sink = new OutputImageType[ newwidth *  newheight * channelcount];
01209   }
01210 
01211   scalex = (float)(src.GetWidth() - 1)/ (float) (newwidth);
01212   scaley = (float)(src.GetHeight() - 1) / (float) (newheight);
01213   if (src.IsInterleaved()) {
01214     for (unsigned int y = 0; y<newheight; y++)
01215       for (unsigned int x = 0; x<newwidth; x++)
01216         for (unsigned int c = 0; c<channelcount; c++) {
01217           oldx = (float)x * scalex;
01218           oldy = (float)y * scaley;
01219 //           oldx = (float)x / factor;
01220 //           oldy = (float)y / factor;
01221           sink[(y * newwidth + x)* channelcount + c] =
01222       // (OutputImageType)rint(src.BilinearInterpolation(oldx, oldy, c));
01223             (OutputImageType)src.BilinearInterpolation(oldx, oldy, c);
01224         }
01225   } else { // planar
01226     BIASERR("Upscale only implemented on interleaved images");
01227     res = -1;
01228   }
01229   if (identicalimages && res==0){
01230     dst.Release();
01231     dst.Init(newwidth, newheight, channelcount);
01232 
01233     dst.ReleaseImageDataPointer();
01234     dst.RedirectImageDataPointer(sink);
01235   }
01236 
01237   return res;
01238 }
01239 
01240 template <class InputImageType, class OutputImageType>
01241 int Rescale<InputImageType, OutputImageType>::
01242 UpsampleGrey(const Image<InputImageType>& src, Image<OutputImageType>& dst,
01243              const unsigned newwidth, const unsigned newheight)
01244 {
01245   BIASDOUT(D_FILTERBASE_CALLSTACK, "Rescale::UpsampleGrey(src, dst, w, h)\n");
01246 #ifdef BIAS_DEBUG
01247   unsigned tlx, tly, brx, bry;
01248   src.GetROI()->GetCorners(tlx, tly, brx, bry);
01249   if (tlx!=0 || tly!=0 || brx!=src.GetWidth() || bry!=src.GetHeight()){
01250     BIASERR("ROI not yet supported for UpsampleGrey, ignoring");
01251     //BIASASSERT(false);
01252     //BIASABORT;
01253   }
01254 #endif
01255 #ifdef BIAS_DEBUG
01256   if (src.GetChannelCount()!=1 || src.GetColorModel()!=ImageBase::CM_Grey){
01257     BIASERR("UpsampleGrey only for grey images");
01258     return -1;
01259     //BIASASSERT(false);
01260     //BIASABORT;
01261   }
01262 #endif
01263   BIASASSERT(src.GetWidth()<newwidth);
01264   BIASASSERT(src.GetHeight()<newheight);
01265   BIASASSERT(src.GetChannelCount()==1);
01266   BIASASSERT(src.GetColorModel()==ImageBase::CM_Grey);
01267   int res=0;
01268   float oldx, oldy;
01269   float scalex, scaley;
01270   bool identicalimages=false;
01271   OutputImageType *sink;
01272 
01273   if ((void *)src.GetImageData() == (void *)dst.GetImageData())
01274     identicalimages=true;
01275 
01276   if (!identicalimages){
01277     if (!src.SamePixelAndChannelCount(dst))
01278       if (!dst.IsEmpty()) dst.Release();
01279 
01280     if (dst.IsEmpty())
01281       dst.Init(newwidth, newheight, 1);
01282 
01283     sink = dst.GetImageData();
01284   } else {
01285     sink = new OutputImageType[ newwidth *  newheight ];
01286   }
01287 
01288   scalex = (float)(src.GetWidth() - 1)/ (float) (newwidth);
01289   scaley = (float)(src.GetHeight() - 1) / (float) (newheight);
01290   if (src.IsInterleaved()) {
01291     for (unsigned int y = 0; y<newheight; y++)
01292       for (unsigned int x = 0; x<newwidth; x++){
01293         oldx = (float)x * scalex;
01294         oldy = (float)y * scaley;
01295         sink[(y * newwidth + x)] =
01296           src.FastBilinearInterpolationGrey(oldx, oldy);
01297       }
01298   } else { // planar
01299     BIASERR("UpscaleGrey only implemented on interleaved images");
01300     res = -1;
01301   }
01302   if (identicalimages && res==0){
01303     dst.Release();
01304     dst.Init(newwidth, newheight, 1);
01305 
01306     dst.ReleaseImageDataPointer();
01307     dst.RedirectImageDataPointer(sink);
01308   }
01309 
01310   return res;
01311 }
01312 
01313 
01314 // bicubic weights for factor 2 upsampling
01315 #define WY0 -0.03125
01316 #define WY1  0.28125
01317 
01318 #define WX0  -0.125
01319 #define WX1  1.125
01320 
01321 ////////////////////////////////////////////////////////////////////////////
01322 
01323 template <class InputImageType, class OutputImageType>
01324 int Rescale<InputImageType, OutputImageType>::
01325 UpsampleBy2Grey(const Image<InputImageType>& Source,
01326                 Image<OutputImageType>& Destination,
01327                 InterpolationType interpolation_method,
01328                 const double& valmin,
01329                 const double& valmax) {
01330   BIASDOUT(D_FILTERBASE_CALLSTACK, "Rescale::UpsampleBy2Grey\n");
01331   register OutputImageType *dstu=NULL, *dstl=NULL;
01332   register const InputImageType *src=NULL,  *srclend=NULL;
01333   const int srcwidth = Source.GetWidth(),
01334     dstwidth = Destination.GetWidth();
01335 #ifdef BIAS_DEBUG
01336   if (srcwidth*2!=dstwidth) {
01337     BIASERR("wrong image width :"<<srcwidth<<"*2 != "<<dstwidth);
01338     BIASABORT;
01339   }
01340 #endif
01341   dstu = dstl = Destination.GetImageData();
01342   src = Source.GetImageData();
01343   dstl += dstwidth;
01344 
01345   srclend = src + srcwidth;
01346 
01347   switch (interpolation_method) {
01348   case NearestNeighbour: {
01349     register const InputImageType* srcend = src + Source.GetPixelCount();
01350     while (src < srcend){
01351       while (src < srclend){
01352         *dstu++ = *dstl++ = *src;
01353         *dstu++ = *dstl++ = *src++;
01354       }
01355       dstu += dstwidth;
01356       dstl += dstwidth;
01357       srclend += srcwidth;
01358     }
01359     return 0;
01360   }
01361   case Bilinear: {
01362 #ifdef BIAS_DEBUG
01363     if (srcwidth<2) {
01364       BIASERR("need image of more than one line !");
01365       BIASABORT;
01366     }
01367 #endif
01368 
01369     register const InputImageType* srcend =
01370       src + Source.GetPixelCount()- srcwidth;
01371     // main part
01372     while (src < srcend){
01373       const InputImageType* c_srclend = srclend-1;
01374       while (src < c_srclend){
01375         *dstl++ =  OutputImageType(0.5 * (src[0] + src[srcwidth]));
01376         *dstu++ = *src;
01377         *dstu++ = OutputImageType(0.5  * (src[0] + src[1]));
01378         *dstl++ = OutputImageType(0.25 * (src[0] + src[1]+
01379                                           src[srcwidth] + src[srcwidth+1]));
01380         src++;
01381       }
01382       // last column
01383       *dstl++ = OutputImageType(0.5 * (src[0] + src[srcwidth]));
01384       *dstu++ = *src;
01385       *dstu++ = *src;
01386       *dstl = dstl[-1]; dstl++;
01387       src++;
01388 
01389       // dst: skip already completed line
01390       dstu += dstwidth;
01391       dstl += dstwidth;
01392       srclend += srcwidth;
01393     }
01394     srclend-=1;
01395     // last line only copied
01396     while (src < srclend){
01397       *dstu++ = *dstl++ = *src++;
01398       *dstu++ = *dstl++ = OutputImageType(0.5*(src[-1]+src[0]));
01399     }
01400     // last pixel (nearest neighbour)
01401     *dstu++ = *dstl++ = *src;
01402     *dstu++ = *dstl++ = *src;
01403     return 0;
01404   }
01405   case Bicubic: {
01406 #ifdef BIAS_DEBUG
01407     if (srcwidth<4) {
01408       BIASERR("need image of more than three lines !");
01409       BIASABORT;
01410     }
01411 #endif
01412 
01413     // first line bilinear
01414     const InputImageType* c_firstsrclend = srclend-1;
01415     while (src < c_firstsrclend){
01416       *dstl++ = OutputImageType(0.5 * (*src + src[srcwidth]));
01417       *dstu++ = *src;
01418       *dstu++ = OutputImageType(0.5 * (src[0] + src[1]));
01419       *dstl++ = OutputImageType(0.25 * (src[0] + src[1] +
01420                                         src[srcwidth] + src[srcwidth+1]));
01421       src++;
01422     }
01423     // last column
01424     *dstl++ = OutputImageType(0.5 * (src[0] + src[srcwidth]));
01425     *dstu++ = *src;
01426     *dstu++ = *src;
01427     *dstl   = dstl[-1]; dstl++;
01428     src++;
01429     dstu += dstwidth;
01430     dstl += dstwidth;
01431     srclend += srcwidth;
01432 
01433     register const InputImageType* srcend =
01434       Source.GetImageData() + Source.GetPixelCount()- (srcwidth*2+2);
01435 
01436     const int width2 = 2*srcwidth;
01437     const int mwidth = -srcwidth;
01438     // main part
01439     while (src < srcend){
01440       // first column bilinear
01441       *dstl++ = OutputImageType(0.5 * (*src + src[srcwidth]));
01442       *dstu++ = *src;
01443       *dstu++ = OutputImageType(0.5  * (src[0] + src[1]));
01444       *dstl++ = OutputImageType(0.25 * (src[0] + src[1] +
01445                                         src[srcwidth] + src[srcwidth+1]));
01446       src++;
01447 
01448       const InputImageType* c_srclend = srclend - 2;
01449       while (src < c_srclend){
01450         // lower left
01451         register double val = 2.0*((WY0*(src[mwidth    ]+src[width2    ])) +
01452                                    (WY1*(src[ 0        ]+src[srcwidth  ])));
01453         if (val>=valmax) *dstl++ = OutputImageType(valmax);
01454         else if (val<=valmin)  *dstl++ = OutputImageType(valmin);
01455         else  *dstl++ = OutputImageType(val);
01456 
01457         // upper left is original pixel
01458         *dstu++ = OutputImageType(*src);
01459 
01460         // upper right
01461         val = (WX0*(src[-1] + src[2])+WX1*(src[0]+src[1]))*0.5;
01462 
01463         if (val>=valmax) *dstu++ = OutputImageType(valmax);
01464         else if (val<=valmin)  *dstu++ = OutputImageType(valmin);
01465         else  *dstu++ = OutputImageType(val);
01466 
01467         // lower right
01468         val = ((WX0*(src[mwidth - 1] + src[mwidth + 2]) +
01469                 WX1*(src[mwidth    ] + src[mwidth + 1])) +
01470                (WX0*(src[width2 - 1] + src[width2 + 2]) +
01471                 WX1*(src[width2    ] + src[width2 + 1])))*WY0 +
01472           ((WX0*(src[-1        ] + src[ 2        ]) +
01473             WX1*(src[ 0        ] + src[ 1        ])) +
01474            (WX0*(src[srcwidth  - 1] + src[srcwidth  + 2]) +
01475             WX1*(src[srcwidth     ] + src[srcwidth  + 1])))*WY1;
01476 
01477         if (val>=valmax) *dstl++ = OutputImageType(valmax);
01478         else if (val<=valmin)  *dstl++ = OutputImageType(valmin);
01479         else  *dstl++ = OutputImageType(val);
01480 
01481         src++;
01482       }
01483 
01484       // last columns
01485       // bilinear
01486       *dstl++ = OutputImageType(0.5 * (*src + src[srcwidth]));
01487       *dstu++ = *src;
01488       *dstu++ = OutputImageType(0.5 * (src[0] + src[1]));
01489       *dstl++ = OutputImageType(0.25 * (src[0] + src[1] +
01490                                         src[srcwidth] + src[srcwidth+1]));
01491       src++;
01492       // copy
01493       *dstl++ = OutputImageType(0.5 * (src[0] + src[srcwidth]));
01494       *dstu++ = *src;
01495       *dstu++ = *src;
01496       *dstl   = dstl[-1]; dstl++;
01497       src++;
01498 
01499       dstu += dstwidth;
01500       dstl += dstwidth;
01501       srclend += srcwidth;
01502     }
01503 
01504     // bilinear
01505     register const InputImageType* srcendlin =
01506       Source.GetImageData() + Source.GetPixelCount()- srcwidth;
01507     // main part
01508     while (src < srcendlin){
01509       const InputImageType* c_srclend = srclend - 1;
01510       while (src < c_srclend){
01511         *dstl++ = OutputImageType(0.5 * (*src + src[srcwidth]));
01512         *dstu++ = *src++;
01513         *dstu++ = OutputImageType(0.5 * (src[-1] + src[0]));
01514         *dstl++ = OutputImageType(0.25 * (src[-1] + src[0] +
01515                                           src[srcwidth-1] + src[srcwidth]));
01516       }
01517       *dstl++ = OutputImageType(0.5 * (src[0] + src[srcwidth]));
01518       *dstu++ = *src;
01519       *dstu++ = *src;
01520       *dstl   = dstl[-1]; dstl++;
01521       src++;
01522 
01523       dstu += dstwidth;
01524       dstl += dstwidth;
01525       srclend += srcwidth;
01526     }
01527 
01528     srclend-=1;
01529     // last line only copied
01530     while (src < srclend){
01531       *dstu++ = *dstl++ = *src++;
01532       *dstu++ = *dstl++ = OutputImageType(0.5*(src[-1]+src[0]));
01533     }
01534     // last pixel (nearest neighbour)
01535     *dstu++ = *dstl++ = *src;
01536     *dstu++ = *dstl++ = *src;
01537     return 0;
01538   }
01539   default:
01540     BIASERR("interpolation not implemented in rescale");
01541     BIASABORT;
01542   }
01543   return 0;
01544 }
01545 
01546 
01547 template <class InputImageType, class OutputImageType>
01548 int Rescale<InputImageType, OutputImageType>::
01549 UpsampleBy2RGBInterleaved(const Image<InputImageType>& Source,
01550                           Image<OutputImageType>& Destination,
01551                           InterpolationType interpolation_method,
01552                           const double& valmin,
01553                           const double& valmax) {
01554   BIASDOUT(D_FILTERBASE_CALLSTACK, "Rescale::UpsampleBy2RGBInterleaved\n");
01555   register OutputImageType *dstu=NULL, *dstl=NULL;
01556   register const InputImageType *src=NULL,  *srclend=NULL;
01557   const int srcwidth = Source.GetWidth(),
01558     dstwidth = Destination.GetWidth();
01559   const int srcwidth3 = 3*srcwidth, dstwidth3=3*dstwidth;
01560 
01561 #ifdef BIAS_DEBUG
01562   if (srcwidth*2!=dstwidth) {
01563     BIASERR("wrong image width :"<<srcwidth<<"*2 != "<<dstwidth);
01564     BIASABORT;
01565   }
01566   if (Source.GetChannelCount()!=3) {
01567     BIASERR("wrong channel count in source");
01568     BIASABORT;
01569   }
01570   if (Destination.GetChannelCount()!=3) {
01571     BIASERR("wrong channel count in dest");
01572     BIASABORT;
01573   }
01574   if (Destination.IsPlanar() || Source.IsPlanar()) {
01575     BIASERR("only for interleaved images");
01576     BIASABORT;
01577   }
01578 #endif
01579   dstu = dstl = Destination.GetImageData();
01580   src = Source.GetImageData();
01581   dstl += dstwidth3;
01582 
01583   srclend = src + srcwidth3;
01584 
01585   switch (interpolation_method) {
01586   case NearestNeighbour: {
01587     register const InputImageType* srcend = src + Source.GetPixelCount()*3;
01588     while (src < srcend){
01589       while (src < srclend){
01590         *dstu++ = *dstl++ =  src[0];
01591         *dstu++ = *dstl++ =  src[1];
01592         *dstu++ = *dstl++ =  src[2];
01593         *dstu++ = *dstl++ = *src++;
01594         *dstu++ = *dstl++ = *src++;
01595         *dstu++ = *dstl++ = *src++;
01596       }
01597       dstu += dstwidth3;
01598       dstl += dstwidth3;
01599       srclend += srcwidth3;
01600     }
01601     return 0;
01602   }
01603   case Bilinear: {
01604 #ifdef BIAS_DEBUG
01605     if (srcwidth<2) {
01606       BIASERR("need image of more than one line !");
01607       BIASABORT;
01608     }
01609 #endif
01610 
01611     register const InputImageType* srcend =
01612       src + Source.GetPixelCount()*3 - srcwidth3;
01613     // main part
01614     while (src < srcend){
01615       const InputImageType* c_srclend = srclend-3;
01616       while (src < c_srclend){
01617         *dstl++ = OutputImageType(0.5 * (src[0] + src[srcwidth3]));
01618         *dstl++ = OutputImageType(0.5 * (src[1] + src[srcwidth3+1]));
01619         *dstl++ = OutputImageType(0.5 * (src[2] + src[srcwidth3+2]));
01620 
01621         *dstu++ = src[0];
01622         *dstu++ = src[1];
01623         *dstu++ = src[2];
01624 
01625         *dstu++ = OutputImageType(0.5  * (src[0] + src[3]));
01626         *dstu++ = OutputImageType(0.5  * (src[1] + src[4]));
01627         *dstu++ = OutputImageType(0.5  * (src[2] + src[5]));
01628 
01629         *dstl++ = OutputImageType(0.25 * (src[0] + src[3]+
01630                                           src[srcwidth3] + src[srcwidth3+1]));
01631         src++;
01632         *dstl++ = OutputImageType(0.25 * (src[0] + src[3]+
01633                                           src[srcwidth3] + src[srcwidth3+1]));
01634         src++;
01635         *dstl++ = OutputImageType(0.25 * (src[0] + src[3]+
01636                                           src[srcwidth3] + src[srcwidth3+1]));
01637         src++;
01638 
01639       }
01640       // last column
01641       *dstl++ = OutputImageType(0.5 * (src[0] + src[srcwidth3]));
01642       *dstl++ = OutputImageType(0.5 * (src[1] + src[srcwidth3+1]));
01643       *dstl++ = OutputImageType(0.5 * (src[2] + src[srcwidth3+2]));
01644 
01645       *dstu++ = src[0];
01646       *dstu++ = src[1];
01647       *dstu++ = src[2];
01648 
01649       *dstu++ = src[0];
01650       *dstu++ = src[1];
01651       *dstu++ = src[2];
01652 
01653       *dstl = dstl[-3]; dstl++;
01654       *dstl = dstl[-3]; dstl++;
01655       *dstl = dstl[-3]; dstl++;
01656       src+=3;
01657       // dst: skip already completed line
01658       dstu += dstwidth3;
01659       dstl += dstwidth3;
01660       srclend += srcwidth3;
01661     }
01662     srclend-=3;
01663     // last line only copied
01664     while (src < srclend){
01665       *dstu++ = *dstl++ = *src++;
01666       *dstu++ = *dstl++ = *src++;
01667       *dstu++ = *dstl++ = *src++;
01668 
01669       *dstu++ = *dstl++ = OutputImageType(0.5*(src[-3]+src[0]));
01670       *dstu++ = *dstl++ = OutputImageType(0.5*(src[-2]+src[1]));
01671       *dstu++ = *dstl++ = OutputImageType(0.5*(src[-1]+src[2]));
01672     }
01673     // last pixel (nearest neighbour)
01674     *dstu++ = *dstl++ = src[0];
01675     *dstu++ = *dstl++ = src[1];
01676     *dstu++ = *dstl++ = src[2];
01677     *dstu++ = *dstl++ = src[0];
01678     *dstu++ = *dstl++ = src[1];
01679     *dstu++ = *dstl++ = src[2];
01680 
01681     return 0;
01682   }
01683   case Bicubic: {
01684 #ifdef BIAS_DEBUG
01685     if ((srcwidth<4) || Source.GetHeight()<4) {
01686       BIASERR("need image of more than three lines !");
01687       BIASABORT;
01688     }
01689 #endif
01690 
01691     // first line bilinear
01692     const InputImageType* c_firstsrclend = srclend-3;
01693     while (src < c_firstsrclend){
01694       *dstl++ = OutputImageType(0.5 * (src[0] + src[srcwidth3]));
01695       *dstl++ = OutputImageType(0.5 * (src[1] + src[srcwidth3+1]));
01696       *dstl++ = OutputImageType(0.5 * (src[2] + src[srcwidth3+2]));
01697 
01698       *dstu++ = src[0];
01699       *dstu++ = src[1];
01700       *dstu++ = src[2];
01701 
01702       *dstu++ = OutputImageType(0.5 * (src[0] + src[3]));
01703       *dstu++ = OutputImageType(0.5 * (src[1] + src[4]));
01704       *dstu++ = OutputImageType(0.5 * (src[2] + src[5]));
01705 
01706       *dstl++ = OutputImageType(0.25 * (src[0] + src[3] +
01707                                         src[srcwidth3] + src[srcwidth3+3]));
01708       src++;
01709       *dstl++ = OutputImageType(0.25 * (src[0] + src[3] +
01710                                         src[srcwidth3] + src[srcwidth3+3]));
01711       src++;
01712       *dstl++ = OutputImageType(0.25 * (src[0] + src[3] +
01713                                         src[srcwidth3] + src[srcwidth3+3]));
01714       src++;
01715     }
01716     // last column
01717     *dstl++ = OutputImageType(0.5 * (src[0] + src[srcwidth3]));
01718     *dstl++ = OutputImageType(0.5 * (src[1] + src[srcwidth3+1]));
01719     *dstl++ = OutputImageType(0.5 * (src[2] + src[srcwidth3+2]));
01720     *dstu++ = src[0];
01721     *dstu++ = src[1];
01722     *dstu++ = src[2];
01723     *dstu++ = *src++;
01724     *dstu++ = *src++;
01725     *dstu++ = *src++;
01726     *dstl = dstl[-3]; dstl++;
01727     *dstl = dstl[-3]; dstl++;
01728     *dstl = dstl[-3]; dstl++;
01729 
01730     dstu += dstwidth3;
01731     dstl += dstwidth3;
01732     srclend += srcwidth3;
01733 
01734     register const InputImageType* srcend =
01735       Source.GetImageData() + Source.GetPixelCount()*3 - (srcwidth3*2+6);
01736 
01737     const int width2 = 2*srcwidth3;
01738     const int mwidth = -srcwidth3;
01739     // main part
01740     while (src < srcend){
01741       // first column bilinear
01742       *dstl++ = OutputImageType(0.5 * (*src + src[srcwidth3]));
01743       *dstl++ = OutputImageType(0.5 * (src[1] + src[srcwidth3+1]));
01744       *dstl++ = OutputImageType(0.5 * (src[2] + src[srcwidth3+2]));
01745 
01746       *dstu++ = src[0];
01747       *dstu++ = src[0];
01748       *dstu++ = src[0];
01749 
01750       *dstu++ = OutputImageType(0.5  * (src[0] + src[3]));
01751       *dstu++ = OutputImageType(0.5  * (src[1] + src[4]));
01752       *dstu++ = OutputImageType(0.5  * (src[2] + src[5]));
01753 
01754       *dstl++ = OutputImageType(0.25 * (src[0] + src[3] +
01755                                         src[srcwidth3] + src[srcwidth3+1]));
01756       src++;
01757       *dstl++ = OutputImageType(0.25 * (src[0] + src[3] +
01758                                         src[srcwidth3] + src[srcwidth3+1]));
01759       src++;
01760       *dstl++ = OutputImageType(0.25 * (src[0] + src[3] +
01761                                         src[srcwidth3] + src[srcwidth3+1]));
01762       src++;
01763 
01764       const InputImageType* c_srclend = srclend - 6;
01765       while (src < c_srclend){
01766         // lower left
01767         register double val = 2.0*((WY0*(src[mwidth    ]+src[width2    ])) +
01768                                    (WY1*(src[ 0        ]+src[srcwidth3  ])));
01769         if (val>=valmax) *dstl++ = OutputImageType(valmax);
01770         else if (val<=valmin)  *dstl++ = OutputImageType(valmin);
01771         else  *dstl++ = OutputImageType(val);
01772 
01773         val = 2.0*((WY0*(src[mwidth+1    ]+src[width2+1    ])) +
01774                    (WY1*(src[ 1        ]+src[srcwidth3+1  ])));
01775         if (val>=valmax) *dstl++ = OutputImageType(valmax);
01776         else if (val<=valmin)  *dstl++ = OutputImageType(valmin);
01777         else  *dstl++ = OutputImageType(val);
01778 
01779         val = 2.0*((WY0*(src[mwidth+2    ]+src[width2+2    ])) +
01780                    (WY1*(src[ 2        ]+src[srcwidth3+2  ])));
01781         if (val>=valmax) *dstl++ = OutputImageType(valmax);
01782         else if (val<=valmin)  *dstl++ = OutputImageType(valmin);
01783         else  *dstl++ = OutputImageType(val);
01784 
01785         // upper left is original pixel
01786         *dstu++ = OutputImageType(src[0]);
01787         *dstu++ = OutputImageType(src[1]);
01788         *dstu++ = OutputImageType(src[2]);
01789 
01790         // upper right
01791         val = (WX0*(src[-3] + src[6])+WX1*(src[0]+src[3]))*0.5;
01792         if (val>=valmax) *dstu++ = OutputImageType(valmax);
01793         else if (val<=valmin)  *dstu++ = OutputImageType(valmin);
01794         else  *dstu++ = OutputImageType(val);
01795 
01796         val = (WX0*(src[-2] + src[7])+WX1*(src[1]+src[4]))*0.5;
01797         if (val>=valmax) *dstu++ = OutputImageType(valmax);
01798         else if (val<=valmin)  *dstu++ = OutputImageType(valmin);
01799         else  *dstu++ = OutputImageType(val);
01800 
01801         val = (WX0*(src[-1] + src[8])+WX1*(src[2]+src[5]))*0.5;
01802         if (val>=valmax) *dstu++ = OutputImageType(valmax);
01803         else if (val<=valmin)  *dstu++ = OutputImageType(valmin);
01804         else  *dstu++ = OutputImageType(val);
01805 
01806         // lower right
01807         val = ((WX0*(src[mwidth - 3] + src[mwidth + 6]) +
01808                 WX1*(src[mwidth    ] + src[mwidth + 3])) +
01809                (WX0*(src[width2 - 3] + src[width2 + 6]) +
01810                 WX1*(src[width2    ] + src[width2 + 3])))*WY0 +
01811           ((WX0*(src[-3        ] + src[ 6        ]) +
01812             WX1*(src[ 0        ] + src[ 3        ])) +
01813            (WX0*(src[srcwidth3  - 3] + src[srcwidth3  + 6]) +
01814             WX1*(src[srcwidth3     ] + src[srcwidth3  + 3])))*WY1;
01815         if (val>=valmax) *dstl++ = OutputImageType(valmax);
01816         else if (val<=valmin)  *dstl++ = OutputImageType(valmin);
01817         else  *dstl++ = OutputImageType(val);
01818         src++;
01819         val = ((WX0*(src[mwidth - 3] + src[mwidth + 6]) +
01820                 WX1*(src[mwidth    ] + src[mwidth + 3])) +
01821                (WX0*(src[width2 - 3] + src[width2 + 6]) +
01822                 WX1*(src[width2    ] + src[width2 + 3])))*WY0 +
01823           ((WX0*(src[-3        ] + src[ 6        ]) +
01824             WX1*(src[ 0        ] + src[ 3        ])) +
01825            (WX0*(src[srcwidth3  - 3] + src[srcwidth3  + 6]) +
01826             WX1*(src[srcwidth3     ] + src[srcwidth3  + 3])))*WY1;
01827         if (val>=valmax) *dstl++ = OutputImageType(valmax);
01828         else if (val<=valmin)  *dstl++ = OutputImageType(valmin);
01829         else  *dstl++ = OutputImageType(val);
01830         src++;
01831         val = ((WX0*(src[mwidth - 3] + src[mwidth + 6]) +
01832                 WX1*(src[mwidth    ] + src[mwidth + 3])) +
01833                (WX0*(src[width2 - 3] + src[width2 + 6]) +
01834                 WX1*(src[width2    ] + src[width2 + 3])))*WY0 +
01835           ((WX0*(src[-3        ] + src[ 6        ]) +
01836             WX1*(src[ 0        ] + src[ 3        ])) +
01837            (WX0*(src[srcwidth3  - 3] + src[srcwidth3  + 6]) +
01838             WX1*(src[srcwidth3     ] + src[srcwidth3  + 3])))*WY1;
01839         if (val>=valmax) *dstl++ = OutputImageType(valmax);
01840         else if (val<=valmin)  *dstl++ = OutputImageType(valmin);
01841         else  *dstl++ = OutputImageType(val);
01842         src++;
01843       }
01844 
01845       // last columns
01846       // bilinear
01847       *dstl++ = OutputImageType(0.5 * (src[0] + src[srcwidth3]));
01848       *dstl++ = OutputImageType(0.5 * (src[1] + src[srcwidth3+1]));
01849       *dstl++ = OutputImageType(0.5 * (src[2] + src[srcwidth3+2]));
01850 
01851       *dstu++ = src[0];
01852       *dstu++ = src[1];
01853       *dstu++ = src[2];
01854 
01855       *dstu++ = OutputImageType(0.5 * (src[0] + src[3]));
01856       *dstu++ = OutputImageType(0.5 * (src[1] + src[4]));
01857       *dstu++ = OutputImageType(0.5 * (src[2] + src[5]));
01858 
01859       *dstl++ = OutputImageType(0.25 * (src[0] + src[3] +
01860                                         src[srcwidth] + src[srcwidth+1]));
01861       src++;
01862       *dstl++ = OutputImageType(0.25 * (src[0] + src[3] +
01863                                         src[srcwidth] + src[srcwidth+1]));
01864       src++;
01865       *dstl++ = OutputImageType(0.25 * (src[0] + src[3] +
01866                                         src[srcwidth] + src[srcwidth+1]));
01867       src++;
01868 
01869 
01870       // copy
01871       *dstl++ = OutputImageType(0.5 * (src[0] + src[srcwidth3]));
01872       *dstl++ = OutputImageType(0.5 * (src[1] + src[srcwidth3+1]));
01873       *dstl++ = OutputImageType(0.5 * (src[2] + src[srcwidth3+2]));
01874 
01875       *dstu++ = *src;
01876       *dstu++ = src[1];
01877       *dstu++ = src[2];
01878 
01879       *dstu++ = src[0];
01880       *dstu++ = src[1];
01881       *dstu++ = src[2];
01882 
01883       *dstl = dstl[-3]; dstl++;
01884       src++;
01885       *dstl = dstl[-3]; dstl++;
01886       src++;
01887       *dstl = dstl[-3]; dstl++;
01888       src++;
01889 
01890       dstu += dstwidth3;
01891       dstl += dstwidth3;
01892       srclend += srcwidth3;
01893     }
01894 
01895     // bilinear
01896     register const InputImageType* srcendlin =
01897       Source.GetImageData() + Source.GetPixelCount()- srcwidth3;
01898     // main part
01899     while (src < srcendlin){
01900       const InputImageType* c_srclend = srclend - 3;
01901       while (src < c_srclend){
01902         *dstl++ = OutputImageType(0.5 * (*src + src[srcwidth3]));
01903         *dstl++ = OutputImageType(0.5 * (src[1] + src[srcwidth3+1]));
01904         *dstl++ = OutputImageType(0.5 * (src[2] + src[srcwidth3+2]));
01905 
01906         *dstu++ = *src++;
01907         *dstu++ = *src++;
01908         *dstu++ = *src++;
01909 
01910         *dstu++ = OutputImageType(0.5 * (src[-3] + src[0]));
01911         *dstu++ = OutputImageType(0.5 * (src[-2] + src[1]));
01912         *dstu++ = OutputImageType(0.5 * (src[-1] + src[2]));
01913 
01914         *dstl++ = OutputImageType(0.25 * (src[-3] + src[0] +
01915                                           src[srcwidth3-3] + src[srcwidth3]));
01916         *dstl++ = OutputImageType(0.25 * (src[-2] + src[1] +
01917                                           src[srcwidth3-2] + src[srcwidth3]));
01918         *dstl++ = OutputImageType(0.25 * (src[-1] + src[2] +
01919                                           src[srcwidth3-1] + src[srcwidth3]));
01920       }
01921       *dstl++ = OutputImageType(0.5 * (src[0] + src[srcwidth3]));
01922       *dstl++ = OutputImageType(0.5 * (src[1] + src[srcwidth3+1]));
01923       *dstl++ = OutputImageType(0.5 * (src[2] + src[srcwidth3+2]));
01924 
01925       *dstu++ = *src;
01926       *dstu++ = src[1];
01927       *dstu++ = src[2];
01928 
01929       *dstu++ = *src;
01930       *dstu++ = src[1];
01931       *dstu++ = src[2];
01932 
01933       *dstl = dstl[-3]; dstl++;
01934       src++;
01935       *dstl = dstl[-3]; dstl++;
01936       src++;
01937       *dstl = dstl[-3]; dstl++;
01938       src++;
01939 
01940       dstu += dstwidth3;
01941       dstl += dstwidth3;
01942       srclend += srcwidth3;
01943     }
01944 
01945     srclend-=3;
01946     // last line only copied
01947     while (src < srclend){
01948       *dstu++ = *dstl++ = *src++;
01949       *dstu++ = *dstl++ = *src++;
01950       *dstu++ = *dstl++ = *src++;
01951 
01952       *dstu++ = *dstl++ = OutputImageType(0.5*(src[-3]+src[0]));
01953       *dstu++ = *dstl++ = OutputImageType(0.5*(src[-2]+src[1]));
01954       *dstu++ = *dstl++ = OutputImageType(0.5*(src[-1]+src[2]));
01955     }
01956     // last pixel (nearest neighbour)
01957     *dstu++ = *dstl++ = *src;
01958     *dstu++ = *dstl++ = src[1];
01959     *dstu++ = *dstl++ = src[2];
01960 
01961     *dstu++ = *dstl++ = *src;
01962     *dstu++ = *dstl++ = src[1];
01963     *dstu++ = *dstl++ = src[2];
01964 
01965 
01966 
01967     return 0;
01968   }
01969   default:
01970     BIASERR("interpolation not implemented in rescale");
01971     BIASABORT;
01972   }
01973   return 0;
01974 }
01975 
01976 
01977 template <class InputImageType, class OutputImageType>
01978 int Rescale<InputImageType, OutputImageType>::
01979 UpsampleBy2(const Image<InputImageType>& Source,
01980             Image<OutputImageType>& Destination)
01981 {
01982   BIASDOUT(D_FILTERBASE_CALLSTACK, "Rescale::UpsampleBy2\n");
01983 #ifdef BIAS_DEBUG
01984   if (!Source.IsInterleaved()){
01985     BIASERR("only for interleaved images");
01986     //BIASASSERT(false);
01987     BIASABORT;
01988   }
01989 #endif
01990   register OutputImageType *dstu=NULL, *dstl=NULL;
01991   register const InputImageType *src=NULL, *srcend=NULL,  *srclend=NULL;
01992   register int srcwidth=Source.GetWidthStep(),
01993     dstwidth=Destination.GetWidthStep(),
01994     cc=Source.GetChannelCount();
01995   dstu = dstl = Destination.GetImageData();
01996   src = Source.GetImageData();
01997   dstl += dstwidth;
01998   srcend = src + Source.GetPixelCount() * cc - cc;
01999   srclend = src + srcwidth - cc;
02000 
02001   while (src < srcend){
02002     while (src < srclend){
02003       for (int i=0; i<cc; i++){
02004         *dstu = *dstl = dstu[cc] = dstl[cc] = *src;
02005         src++; dstu++; dstl++;
02006       }
02007       dstu+=cc;
02008       dstl+=cc;
02009     }
02010     dstu += dstwidth;
02011     dstl += dstwidth;
02012     srclend += srcwidth;
02013   }
02014   return 0;
02015 }
02016 
02017 
02018 //////////////////////////////////////////////////////////////////////////
02019 //               protected functions
02020 //////////////////////////////////////////////////////////////////////////
02021 
02022 template <class InputImageType, class OutputImageType>
02023 void Rescale<InputImageType, OutputImageType>::
02024 UpdateLowpass(double factor)
02025 {
02026   Gauss<InputImageType, OutputImageType>* gauss =
02027     dynamic_cast<Gauss<InputImageType, OutputImageType>*>(_lowpass);
02028   if (gauss)
02029   {
02030     gauss->SetSigma(_RescaleMeanSigmaFactor * factor);
02031   }
02032 
02033   Binomial<InputImageType, OutputImageType>* binomial =
02034     dynamic_cast<Binomial<InputImageType, OutputImageType>*>(_lowpass);
02035   if (binomial)
02036   {
02037     binomial->SetHalfWinSize((int)ceil(1.0*factor));
02038   }
02039 
02040   Mean<InputImageType, OutputImageType>* mean =
02041     dynamic_cast<Mean<InputImageType, OutputImageType>*>(_lowpass);
02042   if (mean)
02043   {
02044     mean->SetSize((int)ceil(1.0*factor));
02045   }
02046 }
02047 
02048 template <class InputImageType, class OutputImageType>
02049 void Rescale<InputImageType, OutputImageType>::
02050 SetLowPassFilter(const FilterNToN<InputImageType, OutputImageType>& lowpass)
02051 {
02052   if (_lowpass)
02053     delete _lowpass;
02054   _lowpass = lowpass.Clone();
02055   UpdateLowpass(_RescaleFactor);
02056 }
02057 
02058 template <class InputImageType, class OutputImageType>
02059 void Rescale<InputImageType, OutputImageType>::
02060 SetLowPassType(int lpt)
02061 {
02062   if (_lowpass)
02063     delete _lowpass;
02064 
02065   switch (lpt)
02066   {
02067   case LPT_Mean:
02068     _lowpass = new Mean<InputImageType,OutputImageType>;
02069     break;
02070   case LPT_Gauss:
02071     _lowpass = new Gauss<InputImageType,OutputImageType>;
02072     break;
02073   case LPT_Binomial:
02074     _lowpass = new Binomial<InputImageType,OutputImageType>;
02075     break;
02076   case LPT_GaussThreshold:
02077     _lowpass = new GaussThreshold<InputImageType,OutputImageType>;
02078     break;
02079   default:
02080       BIASERR("Invalid LowPassType: " << lpt);
02081       break;
02082   }
02083 
02084   UpdateLowpass(_RescaleFactor) ;
02085 }
02086 
02087 template <class InputImageType, class OutputImageType>
02088 int Rescale<InputImageType, OutputImageType>::
02089 _ApplyMeanFilter(const Image<InputImageType>& src,
02090                  Image<OutputImageType>& dst, const double factor)
02091 {
02092   BIASDOUT(D_FILTERBASE_CALLSTACK, "Rescale::_ApplyMeanFilter\n");
02093   UpdateLowpass(factor);
02094   //_lowpass->AddDebugLevel(D_CONV_KERNEL);
02095 //  ImageIO::Save("beforeDownsamplingSrc.mip",src);
02096 //  cout << "---------------" << endl;
02097   int res=_lowpass->Filter(src, dst);
02098 //  ImageIO::Save("beforeDownsamplingDst.mip",dst);
02099 //  cout << "---------------" << endl;
02100 #ifdef BIAS_DEBUG
02101   if (this->DebugLevelIsSet(D_RESCALE_WRITE_IMG)) {
02102     //ImageIO::Save("rescale-lowpass-filtered.mip", dst);
02103     ImageIO::Save("rescale-lowpass-filtered.mip", dst);
02104   }
02105 #endif
02106   return res;
02107 }
02108 
02109 
02110 template <class InputImageType, class OutputImageType>
02111 int Rescale<InputImageType, OutputImageType>::
02112 _FillInterpolatedGrey(const Image<OutputImageType>& src,
02113                       Image<OutputImageType>& dst)
02114 {
02115   BIASDOUT(D_FILTERBASE_CALLSTACK, "Rescale::_FillInterpolatedGrey\n");
02116 #ifdef BIAS_DEBUG
02117   BIASASSERT(src.GetChannelCount()==1 &&
02118              src.GetColorModel()==ImageBase::CM_Grey);
02119   BIASASSERT(!dst.IsEmpty());
02120 
02121   if (this->DebugLevelIsSet(D_RESCALE_WRITE_IMG)) {
02122     ImageBase im=src;
02123     //ImageIO::Save("filtered.mip", im);
02124     ImageIO::Save("filtered.mip", im);
02125   }
02126 #endif
02127   int tlx, tly, brx, bry;
02128   dst.GetROI()->GetCorners(tlx, tly, brx, bry);
02129   const int minx=tlx, maxx=brx, miny=tly, maxy=bry;
02130   OutputImageType **idst=dst.GetImageDataArray();
02131   const double facx = (double)src.GetWidth() / (double)dst.GetWidth();
02132   const double facy = (double)src.GetHeight() / (double)dst.GetHeight();
02133   register int x, y;
02134   for (y=miny; y<maxy; y++){
02135     for (x=minx; x<maxx; x++){
02136       idst[y][x]=(OutputImageType)
02137         src.FastBilinearInterpolationGrey(x*facx, y*facy);
02138 
02139       //      cerr << "( "<<x<<", "<<y<<")  ->  ( "<<x*facx<<", "<<y*facy<<")\n";
02140     }
02141   }
02142   return 0;
02143 }
02144 
02145 template <>
02146 int Rescale<unsigned char, unsigned char>::
02147 _FillInterpolatedGrey(const Image<unsigned char>& src,
02148                       Image<unsigned char>& dst)
02149 {
02150   BIASDOUT(D_FILTERBASE_CALLSTACK, "Rescale::_FillInterpolatedGrey\n");
02151 #ifdef BIAS_DEBUG
02152   BIASASSERT(src.GetChannelCount()==1 &&
02153              src.GetColorModel()==ImageBase::CM_Grey);
02154   BIASASSERT(!dst.IsEmpty());
02155 
02156   if (this->DebugLevelIsSet(D_RESCALE_WRITE_IMG)) {
02157     ImageBase im=src;
02158     //ImageIO::Save("filtered.mip", im);
02159     ImageIO::Save("filtered.mip", im);
02160   }
02161 #endif
02162   int tlx, tly, brx, bry;
02163   dst.GetROI()->GetCorners(tlx, tly, brx, bry);
02164   const int minx=tlx, maxx=brx, miny=tly, maxy=bry;
02165   unsigned char **idst=dst.GetImageDataArray();
02166   const double facx = (double)src.GetWidth() / (double)dst.GetWidth();
02167   const double facy = (double)src.GetHeight() / (double)dst.GetHeight();
02168   register int x, y;
02169   for (y=miny; y<maxy; y++){
02170     for (x=minx; x<maxx; x++){
02171       idst[y][x]=(unsigned char)
02172         rint(src.FastBilinearInterpolationGrey(x*facx, y*facy));
02173     }
02174   }
02175   return 0;
02176 }
02177 
02178 #ifdef BUILD_IMAGE_USHORT
02179 template <>
02180 int Rescale<unsigned short, unsigned short>::
02181 _FillInterpolatedGrey(const Image<unsigned short>& src,
02182                       Image<unsigned short>& dst)
02183 {
02184   BIASDOUT(D_FILTERBASE_CALLSTACK, "Rescale::_FillInterpolatedGrey\n");
02185 #ifdef BIAS_DEBUG
02186   BIASASSERT(src.GetChannelCount()==1 &&
02187              src.GetColorModel()==ImageBase::CM_Grey);
02188   BIASASSERT(!dst.IsEmpty());
02189 
02190 //   if (this->DebugLevelIsSet(D_RESCALE_WRITE_IMG)) {
02191 //     ImageIO::Save("filtered.mip", src);
02192 //   }
02193 #endif
02194   int tlx, tly, brx, bry;
02195   dst.GetROI()->GetCorners(tlx, tly, brx, bry);
02196   const int minx=tlx, maxx=brx, miny=tly, maxy=bry;
02197   unsigned short **idst=dst.GetImageDataArray();
02198   const double facx = (double)src.GetWidth() / (double)dst.GetWidth();
02199   const double facy = (double)src.GetHeight() / (double)dst.GetHeight();
02200   register int x, y;
02201   for (y=miny; y<maxy; y++){
02202     for (x=minx; x<maxx; x++){
02203       idst[y][x]=(unsigned short)
02204         rint(src.FastBilinearInterpolationGrey(x*facx, y*facy));
02205     }
02206   }
02207   return 0;
02208 }
02209 #endif
02210 
02211 template <class InputImageType, class OutputImageType>
02212 int Rescale<InputImageType, OutputImageType>::
02213 _FillInterpolatedColor(const Image<OutputImageType>& src,
02214                        Image<OutputImageType>& dst)
02215 {
02216   BIASDOUT(D_FILTERBASE_CALLSTACK, "Rescale::_FillInterpolatedColor\n");
02217 #ifdef BIAS_DEBUG
02218   BIASASSERT(!dst.IsEmpty());
02219   BIASASSERT(src.IsInterleaved());
02220 
02221   if (this->DebugLevelIsSet(D_RESCALE_WRITE_IMG)) {
02222     ImageBase im=src;
02223     //ImageIO::Save("filtered.mip", im);
02224     ImageIO::Save("filtered.mip", im);
02225   }
02226 #endif
02227   int tlx, tly, brx, bry;
02228   dst.GetROI()->GetCorners(tlx, tly, brx, bry);
02229   const int minx=tlx, maxx=brx, miny=tly, maxy=bry;
02230   OutputImageType **idst=dst.GetImageDataArray();
02231   const double facx = (double)src.GetWidth() / (double)dst.GetWidth();
02232   const double facy = (double)src.GetHeight() / (double)dst.GetHeight();
02233   register int x, y;
02234   const unsigned int cc = src.GetChannelCount();
02235   unsigned int c;
02236   for (y=miny; y<maxy; y++){
02237     for (x=minx; x<maxx; x++){
02238       for (c=0;c<cc; c++)
02239       idst[y][x*cc+c]=(OutputImageType)
02240         src.BilinearInterpolationRGBInterleaved(x*facx, y*facy, c);
02241     }
02242   }
02243   return 0;
02244 }
02245 
02246 template <>
02247 int Rescale<unsigned char, unsigned char>::
02248 _FillInterpolatedColor(const Image<unsigned char>& src,
02249                        Image<unsigned char>& dst)
02250 {
02251   BIASDOUT(D_FILTERBASE_CALLSTACK, "Rescale::_FillInterpolatedColor\n");
02252 #ifdef BIAS_DEBUG
02253   BIASASSERT(!dst.IsEmpty());
02254   BIASASSERT(src.IsInterleaved());
02255 
02256   if (this->DebugLevelIsSet(D_RESCALE_WRITE_IMG)) {
02257     ImageBase im=src;
02258     //ImageIO::Save("filtered.mip", im);
02259     ImageIO::Save("filtered.mip", im);
02260   }
02261 #endif
02262   int tlx, tly, brx, bry;
02263   dst.GetROI()->GetCorners(tlx, tly, brx, bry);
02264   const int minx=tlx, maxx=brx, miny=tly, maxy=bry;
02265   unsigned char **idst=dst.GetImageDataArray();
02266   const double facx = (double)src.GetWidth() / (double)dst.GetWidth();
02267   const double facy = (double)src.GetHeight() / (double)dst.GetHeight();
02268   register int x, y;
02269   const unsigned int cc = src.GetChannelCount();
02270   unsigned int c;
02271   for (y=miny; y<maxy; y++){
02272     for (x=minx; x<maxx; x++){
02273       for (c=0;c<cc; c++)
02274         idst[y][x*cc+c]=(unsigned char)rint
02275           (src.BilinearInterpolationRGBInterleaved(x*facx, y*facy, c));
02276     }
02277   }
02278   return 0;
02279 }
02280 
02281 #ifdef BUILD_IMAGE_USHORT
02282 // exactly the same as <unsigned char, int>, but partial specializations
02283 // are not supported by many compilers
02284 template <>
02285 int Rescale<unsigned short,unsigned short>::
02286 _FillInterpolatedColor(const Image<unsigned short>& src,
02287                        Image<unsigned short>& dst)
02288 {
02289   BIASDOUT(D_FILTERBASE_CALLSTACK, "Rescale::_FillInterpolatedColor\n");
02290 #ifdef BIAS_DEBUG
02291   BIASASSERT(!dst.IsEmpty());
02292   BIASASSERT(src.IsInterleaved());
02293 
02294 //   if (this->DebugLevelIsSet(D_RESCALE_WRITE_IMG)) {
02295 //     ImageIO::Save("filtered.mip", src);
02296 //   }
02297 #endif
02298 
02299   int tlx, tly, brx, bry;
02300   dst.GetROI()->GetCorners(tlx, tly, brx, bry);
02301   const int minx=tlx, maxx=brx, miny=tly, maxy=bry;
02302   unsigned short **idst=dst.GetImageDataArray();
02303   const double facx = (double)src.GetWidth() / (double)dst.GetWidth();
02304   const double facy = (double)src.GetHeight() / (double)dst.GetHeight();
02305   register int x, y;
02306   const unsigned int cc = src.GetChannelCount();
02307   unsigned int c;
02308   for (y=miny; y<maxy; y++){
02309     for (x=minx; x<maxx; x++){
02310       for (c=0;c<cc; c++)
02311         idst[y][x*cc+c]=(unsigned char)rint
02312           (src.BilinearInterpolationRGBInterleaved(x*facx, y*facy, c));
02313     }
02314   }
02315 
02316 
02317   return 0;
02318 }
02319 #endif
02320 
02321 template <class InputImageType, class OutputImageType>
02322 void Rescale<InputImageType, OutputImageType>::
02323 GetBordersValid_(int &border_x, int &border_y) const
02324 {
02325   if (_RescaleFactor<=1.0){
02326     border_x = border_y = 0;
02327   } else {
02328     _lowpass->GetBorders(border_x, border_y);
02329   }
02330 }
02331 
02332 
02333 template <class InputImageType, class OutputImageType>
02334 int Rescale<InputImageType, OutputImageType>::
02335 _FillInterpolated(const Image<OutputImageType>& src,
02336                   const double &factor, Image<OutputImageType>& dst)
02337 {
02338   BIASDOUT(D_FILTERBASE_CALLSTACK, "Rescale::_FillInterpolatedColor\n");
02339 #ifdef BIAS_DEBUG
02340   BIASASSERT(!dst.IsEmpty());
02341   BIASASSERT(src.IsInterleaved());
02342 
02343   if (this->DebugLevelIsSet(D_RESCALE_WRITE_IMG)) {
02344     ImageBase im=src;
02345     //ImageIO::Save("filtered.mip", im);
02346     ImageIO::Save("filtered.mip", im);
02347   }
02348 #endif
02349   int tlx, tly, brx, bry;
02350   dst.GetROI()->GetCorners(tlx, tly, brx, bry);
02351   const int minx=tlx, maxx=brx, miny=tly, maxy=bry;
02352   OutputImageType **idst=dst.GetImageDataArray();
02353   register int x, xc, y;
02354   const unsigned cc = src.GetChannelCount();
02355   register unsigned c;
02356   double src_x, src_y;
02357   for (y=miny; y<maxy; y++){
02358     src_y = (y+0.5)*factor-0.5;
02359     for (x=minx; x<maxx; x++){
02360       src_x = (x+0.5)*factor-0.5;
02361       xc = x*cc;
02362       for (c=0; c<cc; c++){
02363         idst[y][xc+c] = Cast<OutputImageType>
02364           (src.BilinearInterpolationRGBInterleaved(src_x, src_y, c));
02365       }
02366     }
02367   }
02368   _LastPositionOffset = 0.5*factor-0.5;
02369   return 0;
02370 }
02371 
02372 
02373 //////////////////////////////////////////////////////////////////////////
02374 //                 instantiation
02375 //////////////////////////////////////////////////////////////////////////
02376 #define FILTER_INSTANTIATION_CLASS Rescale
02377 #include "Filterinst.hh"
02378 
02379 } // namespace BIAS
 All Classes Functions Variables Typedefs Enumerations Enumerator Friends