contrib/mul/vil3d/vil3d_resample_tricubic.txx

Go to the documentation of this file.
00001 #ifndef vil3d_resample_tricubic_txx_
00002 #define vil3d_resample_tricubic_txx_
00003 //:
00004 // \file
00005 // \brief Resample a 3D image by a tricubic method
00006 // \author Gwenael Guillard
00007 
00008 #include "vil3d_resample_tricubic.h"
00009 
00010 #include <vil/vil_convert.h>
00011 #include <vil3d/vil3d_tricub_interp.h>
00012 #include <vil3d/vil3d_plane.h>
00013 #include <vcl_cassert.h>
00014 
00015 
00016 inline bool vil3d_grid_corner_in_image(double x0, double y0, double z0,
00017                                        const vil3d_image_view_base& image)
00018 {
00019   if (x0<1.0) return false;
00020   if (x0>image.ni()-3.0) return false;
00021   if (y0<1.0) return false;
00022   if (y0>image.nj()-3.0) return false;
00023   if (z0<1.0) return false;
00024   if (z0>image.nk()-3.0) return false;
00025   return true;
00026 }
00027 
00028 
00029 //: Sample grid of points in one image and place in another, using tricubic interpolation.
00030 //  dest_image(i,j,k,p) is sampled from the src_image at
00031 //  (x0+i.dx1+j.dx2+k.dx3, y0+i.dy1+j.dy2+k.dy3, z0+i.dz1+j.dz2+k.dz3),
00032 //  where i=[0..n1-1], j=[0..n2-1], k=[0..n3-1]
00033 //  dest_image resized to (n1,n2,n3,src_image.nplanes())
00034 //  Points outside interpolatable region return zero or \a outval
00035 template <class S, class T>
00036 void vil3d_resample_tricubic(const vil3d_image_view<S>& src_image,
00037                              vil3d_image_view<T>& dest_image,
00038                              double x0, double y0, double z0,
00039                              double dx1, double dy1, double dz1,
00040                              double dx2, double dy2, double dz2,
00041                              double dx3, double dy3, double dz3,
00042                              int n1, int n2, int n3,
00043                              T outval/*=0*/)
00044 {
00045   bool all_in_image =
00046     vil3d_grid_corner_in_image(x0,
00047                                y0,
00048                                z0,
00049                                src_image) &&
00050     vil3d_grid_corner_in_image(x0 + (n1-1)*dx1,
00051                                y0 + (n1-1)*dy1,
00052                                z0 + (n1-1)*dz1,
00053                                src_image) &&
00054     vil3d_grid_corner_in_image(x0 + (n2-1)*dx2,
00055                                y0 + (n2-1)*dy2,
00056                                z0 + (n2-1)*dz2,
00057                                src_image) &&
00058     vil3d_grid_corner_in_image(x0 + (n1-1)*dx1 + (n2-1)*dx2,
00059                                y0 + (n1-1)*dy1 + (n2-1)*dy2,
00060                                z0 + (n1-1)*dz1 + (n2-1)*dz2,
00061                                src_image) &&
00062     vil3d_grid_corner_in_image(x0 + (n3-1)*dx3,
00063                                y0 + (n3-1)*dy3,
00064                                z0 + (n3-1)*dz3,
00065                                src_image) &&
00066     vil3d_grid_corner_in_image(x0 + (n1-1)*dx1 + (n3-1)*dx3,
00067                                y0 + (n1-1)*dy1 + (n3-1)*dy3,
00068                                z0 + (n1-1)*dz1 + (n3-1)*dz3,
00069                                src_image) &&
00070     vil3d_grid_corner_in_image(x0 + (n2-1)*dx2 + (n3-1)*dx3,
00071                                y0 + (n2-1)*dy2 + (n3-1)*dy3,
00072                                z0 + (n2-1)*dz2 + (n3-1)*dz3,
00073                                src_image) &&
00074     vil3d_grid_corner_in_image(x0 + (n1-1)*dx1 + (n2-1)*dx2 + (n3-1)*dx3,
00075                                y0 + (n1-1)*dy1 + (n2-1)*dy2 + (n3-1)*dy3,
00076                                z0 + (n1-1)*dz1 + (n2-1)*dz2 + (n3-1)*dz3,
00077                                src_image);
00078 
00079   vil_convert_round_pixel<double, T> cast_and_possibly_round;
00080 
00081   const unsigned ni = src_image.ni();
00082   const unsigned nj = src_image.nj();
00083   const unsigned nk = src_image.nk();
00084   const unsigned np = src_image.nplanes();
00085   const vcl_ptrdiff_t istep = src_image.istep();
00086   const vcl_ptrdiff_t jstep = src_image.jstep();
00087   const vcl_ptrdiff_t kstep = src_image.kstep();
00088   const vcl_ptrdiff_t pstep = src_image.planestep();
00089   const S* plane0 = src_image.origin_ptr();
00090 
00091   dest_image.set_size(n1,n2,n3,np);
00092   const vcl_ptrdiff_t d_istep = dest_image.istep();
00093   const vcl_ptrdiff_t d_jstep = dest_image.jstep();
00094   const vcl_ptrdiff_t d_kstep = dest_image.kstep();
00095   const vcl_ptrdiff_t d_pstep = dest_image.planestep();
00096   T* d_plane0 = dest_image.origin_ptr();
00097 
00098   if (all_in_image)
00099   {
00100     if (np==1)
00101     {
00102       double xk=x0, yk=y0, zk=z0;
00103       T *slice = d_plane0;
00104       for (int k=0; k<n3; ++k, xk+=dx3, yk+=dy3, zk+=dz3, slice+=d_kstep)
00105       {
00106         double xj=xk, yj=yk, zj=zk;  // Start of k-th slice
00107         T *row = slice;
00108         for (int j=0; j<n2; ++j, xj+=dx2, yj+=dy2, zj+=dz2, row+=d_jstep)
00109         {
00110           double x=xj, y=yj, z=zj;  // Start of j-th row
00111           T *dpt = row;
00112           for (int i=0; i<n1; ++i, x+=dx1, y+=dy1, z+=dz1, dpt+=d_istep)
00113             cast_and_possibly_round( vil3d_tricub_interp_raw( x, y, z,
00114                                                               plane0,
00115                                                               istep, jstep, kstep),
00116                                      *dpt);
00117         }
00118       }
00119     }
00120     else
00121     {
00122       double xk=x0, yk=y0, zk=z0;
00123       T *slice = d_plane0;
00124       for (int k=0; k<n3; ++k, xk+=dx3, yk+=dy3, zk+=dz3, slice+=d_kstep)
00125       {
00126         double xj=xk, yj=yk, zj=zk;  // Start of k-th slice
00127         T *row = slice;
00128         for (int j=0; j<n2; ++j, xj+=dx2, yj+=dy2, zj+=dz2, row+=d_jstep)
00129         {
00130           double x=xj, y=yj, z=zj;  // Start of j-th row
00131           T *dpt = row;
00132           for (int i=0; i<n1; ++i, x+=dx1, y+=dy1, z+=dz1, dpt+=d_istep)
00133           {
00134             for (unsigned int p=0; p<np; ++p)
00135               cast_and_possibly_round( vil3d_tricub_interp_raw( x, y, z,
00136                                                                 plane0+p*pstep,
00137                                                                 istep, jstep, kstep),
00138                                        dpt[p*d_pstep]);
00139           }
00140         }
00141       }
00142     }
00143   }
00144   else
00145   {
00146     // Use safe interpolation
00147     if (np==1)
00148     {
00149       double xk=x0, yk=y0, zk=z0;
00150       T *slice = d_plane0;
00151       for (int k=0; k<n3; ++k, xk+=dx3, yk+=dy3, zk+=dz3, slice+=d_kstep)
00152       {
00153         double xj=xk, yj=yk, zj=zk;  // Start of k-th slice
00154         T *row = slice;
00155         for (int j=0; j<n2; ++j, xj+=dx2, yj+=dy2, zj+=dz2, row+=d_jstep)
00156         {
00157           double x=xj, y=yj, z=zj;  // Start of j-th row
00158           T *dpt = row;
00159           for (int i=0; i<n1; ++i, x+=dx1, y+=dy1, z+=dz1, dpt+=d_istep)
00160             cast_and_possibly_round( vil3d_tricub_interp_safe( x, y, z,
00161                                                                plane0,
00162                                                                ni, nj, nk,
00163                                                                istep, jstep, kstep,
00164                                                                outval),
00165                                      *dpt);
00166         }
00167       }
00168     }
00169     else
00170     {
00171       double xk=x0, yk=y0, zk=z0;
00172       T *slice = d_plane0;
00173       for (int k=0; k<n3; ++k, xk+=dx3, yk+=dy3, zk+=dz3, slice+=d_kstep)
00174       {
00175         double xj=xk, yj=yk, zj=zk;  // Start of k-th slice
00176         T *row = slice;
00177         for (int j=0; j<n2; ++j, xj+=dx2, yj+=dy2, zj+=dz2, row+=d_jstep)
00178         {
00179           double x=xj, y=yj, z=zj;  // Start of j-th row
00180           T *dpt = row;
00181           for (int i=0; i<n1; ++i, x+=dx1, y+=dy1, z+=dz1, dpt+=d_istep)
00182           {
00183             for (unsigned int p=0; p<np; ++p)
00184               cast_and_possibly_round( vil3d_tricub_interp_safe( x, y, z,
00185                                                                  plane0+p*pstep,
00186                                                                  ni, nj, nk,
00187                                                                  istep, jstep, kstep,
00188                                                                  outval),
00189                                        dpt[p*d_pstep]);
00190           }
00191         }
00192       }
00193     }
00194   }
00195 }
00196 
00197 //: Sample grid of points in one image and place in another, using tricubic interpolation.
00198 //  dest_image(i,j,k,p) is sampled from the src_image at
00199 //  (x0+i.dx1+j.dx2+k.dx3, y0+i.dy1+j.dy2+k.dy3, z0+i.dz1+j.dz2+k.dz3),
00200 //  where i=[0..n1-1], j=[0..n2-1], k=[0..n3-1]
00201 //  dest_image resized to (n1,n2,n3,src_image.nplanes())
00202 //  Points outside interpolatable return the value of the nearest valid pixel.
00203 template <class S, class T>
00204 void vil3d_resample_tricubic_edge_extend(const vil3d_image_view<S>& src_image,
00205                                          vil3d_image_view<T>& dest_image,
00206                                          double x0, double y0, double z0,
00207                                          double dx1, double dy1, double dz1,
00208                                          double dx2, double dy2, double dz2,
00209                                          double dx3, double dy3, double dz3,
00210                                          int n1, int n2, int n3)
00211 {
00212   bool all_in_image =
00213     vil3d_grid_corner_in_image(x0,
00214                                y0,
00215                                z0,
00216                                src_image) &&
00217     vil3d_grid_corner_in_image(x0 + (n1-1)*dx1,
00218                                y0 + (n1-1)*dy1,
00219                                z0 + (n1-1)*dz1,
00220                                src_image) &&
00221     vil3d_grid_corner_in_image(x0 + (n2-1)*dx2,
00222                                y0 + (n2-1)*dy2,
00223                                z0 + (n2-1)*dz2,
00224                                src_image) &&
00225     vil3d_grid_corner_in_image(x0 + (n1-1)*dx1 + (n2-1)*dx2,
00226                                y0 + (n1-1)*dy1 + (n2-1)*dy2,
00227                                z0 + (n1-1)*dz1 + (n2-1)*dz2,
00228                                src_image) &&
00229     vil3d_grid_corner_in_image(x0 + (n3-1)*dx3,
00230                                y0 + (n3-1)*dy3,
00231                                z0 + (n3-1)*dz3,
00232                                src_image) &&
00233     vil3d_grid_corner_in_image(x0 + (n1-1)*dx1 + (n3-1)*dx3,
00234                                y0 + (n1-1)*dy1 + (n3-1)*dy3,
00235                                z0 + (n1-1)*dz1 + (n3-1)*dz3,
00236                                src_image) &&
00237     vil3d_grid_corner_in_image(x0 + (n2-1)*dx2 + (n3-1)*dx3,
00238                                y0 + (n2-1)*dy2 + (n3-1)*dy3,
00239                                z0 + (n2-1)*dz2 + (n3-1)*dz3,
00240                                src_image) &&
00241     vil3d_grid_corner_in_image(x0 + (n1-1)*dx1 + (n2-1)*dx2 + (n3-1)*dx3,
00242                                y0 + (n1-1)*dy1 + (n2-1)*dy2 + (n3-1)*dy3,
00243                                z0 + (n1-1)*dz1 + (n2-1)*dz2 + (n3-1)*dz3,
00244                                src_image);
00245 
00246   vil_convert_round_pixel<double, T> cast_and_possibly_round;
00247 
00248   const unsigned ni = src_image.ni();
00249   const unsigned nj = src_image.nj();
00250   const unsigned nk = src_image.nk();
00251   const unsigned np = src_image.nplanes();
00252   const vcl_ptrdiff_t istep = src_image.istep();
00253   const vcl_ptrdiff_t jstep = src_image.jstep();
00254   const vcl_ptrdiff_t kstep = src_image.kstep();
00255   const vcl_ptrdiff_t pstep = src_image.planestep();
00256   const S* plane0 = src_image.origin_ptr();
00257 
00258   dest_image.set_size(n1,n2,n3,np);
00259   const vcl_ptrdiff_t d_istep = dest_image.istep();
00260   const vcl_ptrdiff_t d_jstep = dest_image.jstep();
00261   const vcl_ptrdiff_t d_kstep = dest_image.kstep();
00262   const vcl_ptrdiff_t d_pstep = dest_image.planestep();
00263   T* d_plane0 = dest_image.origin_ptr();
00264 
00265   if ( all_in_image )
00266   {
00267     if (np==1)
00268     {
00269       double xk=x0, yk=y0, zk=z0;
00270       T *slice = d_plane0;
00271       for (int k=0; k<n3; ++k, xk+=dx3, yk+=dy3, zk+=dz3, slice+=d_kstep)
00272       {
00273         double xj=xk, yj=yk, zj=zk;  // Start of k-th slice
00274         T *row = slice;
00275         for (int j=0; j<n2; ++j, xj+=dx2, yj+=dy2, zj+=dz2, row+=d_jstep)
00276         {
00277           double x=xj, y=yj, z=zj;  // Start of j-th row
00278           T *dpt = row;
00279           for (int i=0; i<n1; ++i, x+=dx1, y+=dy1, z+=dz1, dpt+=d_istep)
00280             cast_and_possibly_round( vil3d_tricub_interp_raw(x,y,z,
00281                                                              plane0,
00282                                                              istep,jstep,kstep),
00283                                      *dpt );
00284         }
00285       }
00286     }
00287     else
00288     {
00289       double xk=x0, yk=y0, zk=z0;
00290       T *slice = d_plane0;
00291       for (int k=0; k<n3; ++k, xk+=dx3, yk+=dy3, zk+=dz3, slice+=d_kstep)
00292       {
00293         double xj=xk, yj=yk, zj=zk;  // Start of k-th slice
00294         T *row = slice;
00295         for (int j=0; j<n2; ++j, xj+=dx2, yj+=dy2, zj+=dz2, row+=d_jstep)
00296         {
00297           double x=xj, y=yj, z=zj;  // Start of j-th row
00298           T *dpt = row;
00299           for (int i=0; i<n1; ++i, x+=dx1, y+=dy1, z+=dz1, dpt+=d_istep)
00300           {
00301             for (unsigned int p=0; p<np; ++p)
00302               cast_and_possibly_round( vil3d_tricub_interp_raw(x,y,z,
00303                                                                plane0+p*pstep,
00304                                                                istep,jstep,kstep),
00305                                        dpt[p*d_pstep] );
00306           }
00307         }
00308       }
00309     }
00310   }
00311   else
00312   {
00313     if (np==1)
00314     {
00315       double xk=x0, yk=y0, zk=z0;
00316       T *slice = d_plane0;
00317       for (int k=0; k<n3; ++k, xk+=dx3, yk+=dy3, zk+=dz3, slice+=d_kstep)
00318       {
00319         double xj=xk, yj=yk, zj=zk;  // Start of k-th slice
00320         T *row = slice;
00321         for (int j=0; j<n2; ++j, xj+=dx2, yj+=dy2, zj+=dz2, row+=d_jstep)
00322         {
00323           double x=xj, y=yj, z=zj;  // Start of j-th row
00324           T *dpt = row;
00325           for (int i=0; i<n1; ++i, x+=dx1, y+=dy1, z+=dz1, dpt+=d_istep)
00326             cast_and_possibly_round( vil3d_tricub_interp_safe_extend(x,y,z,
00327                                                                      plane0,
00328                                                                      ni,nj,nk,
00329                                                                      istep,jstep,kstep),
00330                                      *dpt );
00331         }
00332       }
00333     }
00334     else
00335     {
00336       double xk=x0, yk=y0, zk=z0;
00337       T *slice = d_plane0;
00338       for (int k=0; k<n3; ++k, xk+=dx3, yk+=dy3, zk+=dz3, slice+=d_kstep)
00339       {
00340         double xj=xk, yj=yk, zj=zk;  // Start of k-th slice
00341         T *row = slice;
00342         for (int j=0; j<n2; ++j, xj+=dx2, yj+=dy2, zj+=dz2, row+=d_jstep)
00343         {
00344           double x=xj, y=yj, z=zj;  // Start of j-th row
00345           T *dpt = row;
00346           for (int i=0; i<n1; ++i, x+=dx1, y+=dy1, z+=dz1, dpt+=d_istep)
00347           {
00348             for (unsigned int p=0; p<np; ++p)
00349               cast_and_possibly_round( vil3d_tricub_interp_safe_extend(x,y,z,
00350                                                                        plane0+p*pstep,
00351                                                                        ni,nj,nk,
00352                                                                        istep,jstep,kstep),
00353                                        dpt[p*d_pstep] );
00354           }
00355         }
00356       }
00357     }
00358   }
00359 }
00360 
00361 
00362 //: Sample grid of points in one image and place in another, using tricubic interpolation and edge extension.
00363 //  dest_image(i,j,k,p) is sampled from the src_image at
00364 //  (x0+i.dx1+j.dx2+k.dx3, y0+i.dy1+j.dy2+k.dy3, z0+i.dz1+j.dz2+k.dz3),
00365 //  where i=[0..n1-1], j=[0..n2-1], k=[0..n3-1].
00366 //  dest_image resized to (n1,n2,n3,src_image.nplanes())
00367 //  Points outside interpolatable return the value of the nearest valid pixel.
00368 template <class S, class T>
00369 void vil3d_resample_tricubic_edge_extend(const vil3d_image_view<S>& src_image,
00370                                          vil3d_image_view<T>& dest_image,
00371                                          int n1, int n2, int n3)
00372 {
00373   double f= 0.9999999; // so sampler doesn't go off edge of image
00374   double x0=0;
00375   double y0=0;
00376   double z0=0;
00377 
00378   double dx1=f*(src_image.ni()-1)*1.0/(n1-1);
00379   double dy1=0;
00380   double dz1=0;
00381 
00382   double dx2=0;
00383   double dy2=f*(src_image.nj()-1)*1.0/(n2-1);
00384   double dz2=0;
00385 
00386   double dx3=0;
00387   double dy3=0;
00388   double dz3=f*(src_image.nk()-1)*1.0/(n3-1);
00389 
00390   vil3d_resample_tricubic_edge_extend(src_image, dest_image,
00391                                       x0, y0, z0,
00392                                       dx1, dy1, dz1, dx2, dy2, dz2, dx3, dy3, dz3,
00393                                       n1, n2, n3);
00394 }
00395 
00396 
00397 
00398 //: Sample grid of points in one image and place in another, using tricubic interpolation.
00399 //  dest_image(i,j,k,p) is sampled from the src_image at
00400 //  (x0+i.dx1+j.dx2+k.dx3, y0+i.dy1+j.dy2+k.dy3, z0+i.dz1+j.dz2+k.dz3),
00401 //  where i=[0..n1-1], j=[0..n2-1], k=[0..n3-1]
00402 //  dest_image resized to (n1,n2,n3,src_image.nplanes())
00403 //  Points outside interpolatable return the trilinear interpolated value of the nearest valid pixels.
00404 template <class S, class T>
00405 void vil3d_resample_tricubic_edge_trilin_extend(const vil3d_image_view<S>& src_image,
00406                                                 vil3d_image_view<T>& dest_image,
00407                                                 double x0, double y0, double z0,
00408                                                 double dx1, double dy1, double dz1,
00409                                                 double dx2, double dy2, double dz2,
00410                                                 double dx3, double dy3, double dz3,
00411                                                 int n1, int n2, int n3)
00412 {
00413   vil_convert_round_pixel<double, T> cast_and_possibly_round;
00414 
00415   const unsigned ni = src_image.ni();
00416   const unsigned nj = src_image.nj();
00417   const unsigned nk = src_image.nk();
00418   const unsigned np = src_image.nplanes();
00419   const vcl_ptrdiff_t istep = src_image.istep();
00420   const vcl_ptrdiff_t jstep = src_image.jstep();
00421   const vcl_ptrdiff_t kstep = src_image.kstep();
00422   const vcl_ptrdiff_t pstep = src_image.planestep();
00423   const S* plane0 = src_image.origin_ptr();
00424 
00425   dest_image.set_size(n1,n2,n3,np);
00426   const vcl_ptrdiff_t d_istep = dest_image.istep();
00427   const vcl_ptrdiff_t d_jstep = dest_image.jstep();
00428   const vcl_ptrdiff_t d_kstep = dest_image.kstep();
00429   const vcl_ptrdiff_t d_pstep = dest_image.planestep();
00430   T* d_plane0 = dest_image.origin_ptr();
00431   // Use safe interpolation with edge-extension
00432   if (np==1)
00433   {
00434     double xk=x0, yk=y0, zk=z0;
00435     T *slice = d_plane0;
00436     for (int k=0; k<n3; ++k, xk+=dx3, yk+=dy3, zk+=dz3, slice+=d_kstep)
00437     {
00438       double xj=xk, yj=yk, zj=zk;  // Start of k-th slice
00439       T *row = slice;
00440       for (int j=0; j<n2; ++j, xj+=dx2, yj+=dy2, zj+=dz2, row+=d_jstep)
00441       {
00442         double x=xj, y=yj, z=zj;  // Start of j-th row
00443         T *dpt = row;
00444         for (int i=0; i<n1; ++i, x+=dx1, y+=dy1, z+=dz1, dpt+=d_istep)
00445           cast_and_possibly_round( vil3d_tricub_interp_safe_trilinear_extend(x,y,z,
00446                                                                              plane0,
00447                                                                              ni,nj,nk,
00448                                                                              istep,jstep,kstep),
00449                                    *dpt );
00450       }
00451     }
00452   }
00453   else
00454   {
00455     double xk=x0, yk=y0, zk=z0;
00456     T *slice = d_plane0;
00457     for (int k=0; k<n3; ++k, xk+=dx3, yk+=dy3, zk+=dz3, slice+=d_kstep)
00458     {
00459       double xj=xk, yj=yk, zj=zk;  // Start of k-th slice
00460       T *row = slice;
00461       for (int j=0; j<n2; ++j, xj+=dx2, yj+=dy2, zj+=dz2, row+=d_jstep)
00462       {
00463         double x=xj, y=yj, z=zj;  // Start of j-th row
00464         T *dpt = row;
00465         for (int i=0; i<n1; ++i, x+=dx1, y+=dy1, z+=dz1, dpt+=d_istep)
00466         {
00467           for (unsigned int p=0; p<np; ++p)
00468             cast_and_possibly_round( vil3d_tricub_interp_safe_trilinear_extend(x,y,z,
00469                                                                                plane0+p*pstep,
00470                                                                                ni,nj,nk,
00471                                                                                istep,jstep,kstep),
00472                                      dpt[p*d_pstep] );
00473         }
00474       }
00475     }
00476   }
00477 }
00478 
00479 //: Sample grid of points in one image and place in another, using tricubic interpolation.
00480 //  dest_image(i,j,k,p) is sampled from the src_image at
00481 //  (x0+i.dx1+j.dx2+k.dx3, y0+i.dy1+j.dy2+k.dy3, z0+i.dz1+j.dz2+k.dz3),
00482 //  where i=[0..n1-1], j=[0..n2-1], k=[0..n3-1]
00483 //  dest_image resized to (n1,n2,n3,src_image.nplanes())
00484 //  Points outside interpolatable return the trilinear interpolated value of the nearest valid pixels.
00485 template <class S, class T>
00486 void vil3d_resample_tricubic_edge_trilin_extend(const vil3d_image_view<S>& src_image,
00487                                                 vil3d_image_view<T>& dest_image,
00488                                                 int n1, int n2, int n3)
00489 {
00490   double f= 0.9999999; // so sampler doesn't go off edge of image
00491   double x0=0;
00492   double y0=0;
00493   double z0=0;
00494 
00495   double dx1=f*(src_image.ni()-1)*1.0/(n1-1);
00496   double dy1=0;
00497   double dz1=0;
00498 
00499   double dx2=0;
00500   double dy2=f*(src_image.nj()-1)*1.0/(n2-1);
00501   double dz2=0;
00502 
00503   double dx3=0;
00504   double dy3=0;
00505   double dz3=f*(src_image.nk()-1)*1.0/(n3-1);
00506 
00507   vil3d_resample_tricubic_edge_trilin_extend(src_image, dest_image,
00508                                              x0, y0, z0,
00509                                              dx1, dy1, dz1, dx2, dy2, dz2, dx3, dy3, dz3,
00510                                              n1, n2, n3);
00511 }
00512 
00513 
00514 //: Resample image to a specified dimensions (n1 * n2 * n3)
00515 template <class S, class T>
00516 void vil3d_resample_tricubic(const vil3d_image_view<S>& src_image,
00517                              vil3d_image_view<T>& dest_image,
00518                              int n1, int n2, int n3)
00519 {
00520   double f= 0.9999999; // so sampler doesn't go off edge of image
00521   double x0=0;
00522   double y0=0;
00523   double z0=0;
00524 
00525   double dx1=f*(src_image.ni()-1)*1.0/(n1-1);
00526   double dy1=0;
00527   double dz1=0;
00528 
00529   double dx2=0;
00530   double dy2=f*(src_image.nj()-1)*1.0/(n2-1);
00531   double dz2=0;
00532 
00533   double dx3=0;
00534   double dy3=0;
00535   double dz3=f*(src_image.nk()-1)*1.0/(n3-1);
00536 
00537   vil3d_resample_tricubic(src_image, dest_image,
00538                           x0, y0, z0,
00539                           dx1, dy1, dz1,
00540                           dx2, dy2, dz2,
00541                           dx3, dy3, dz3,
00542                           n1, n2, n3);
00543 }
00544 
00545 
00546 #define VIL3D_RESAMPLE_TRICUBIC_INSTANTIATE( S, T ) \
00547 template void vil3d_resample_tricubic(const vil3d_image_view< S >& src_image, \
00548                                       vil3d_image_view< T >& dest_image, \
00549                                       double x0, double y0, double z0, \
00550                                       double dx1, double dy1, double dz1, \
00551                                       double dx2, double dy2, double dz2, \
00552                                       double dx3, double dy3, double dz3, \
00553                                       int n1, int n2, int n3, \
00554                                       T outval); \
00555 template void vil3d_resample_tricubic_edge_extend(const vil3d_image_view< S >& src_image, \
00556                                                   vil3d_image_view< T >& dest_image, \
00557                                                   double x0, double y0, double z0, \
00558                                                   double dx1, double dy1, double dz1, \
00559                                                   double dx2, double dy2, double dz2, \
00560                                                   double dx3, double dy3, double dz3, \
00561                                                   int n1, int n2, int n3); \
00562 template void vil3d_resample_tricubic_edge_extend(const vil3d_image_view< S >& src_image, \
00563                                                   vil3d_image_view< T >& dest_image, \
00564                                                   int n1, int n2, int n3); \
00565 template void vil3d_resample_tricubic_edge_trilin_extend(const vil3d_image_view< S >& src_image, \
00566                                                          vil3d_image_view< T >& dest_image, \
00567                                                          double x0, double y0, double z0, \
00568                                                          double dx1, double dy1, double dz1, \
00569                                                          double dx2, double dy2, double dz2, \
00570                                                          double dx3, double dy3, double dz3, \
00571                                                          int n1, int n2, int n3); \
00572 template void vil3d_resample_tricubic_edge_trilin_extend(const vil3d_image_view< S >& src_image, \
00573                                                          vil3d_image_view< T >& dest_image, \
00574                                                          int n1, int n2, int n3); \
00575 template void vil3d_resample_tricubic(const vil3d_image_view< S >& src_image, \
00576                                       vil3d_image_view< T >& dest_image, \
00577                                       int n1, int n2, int n3)
00578 #if 0
00579 template void vil3d_resample_trilinear(const vil3d_image_view< T >& src_image, \
00580                                        vil3d_image_view< T >& dst_image, \
00581                                        const double dx, \
00582                                        const double dy, \
00583                                        const double dz); \
00584 template void vil3d_resample_trilinear_scale_2(const vil3d_image_view< T >& src_image, \
00585                                                vil3d_image_view< T >& dst_image)
00586 #endif
00587 
00588 #endif // vil3d_resample_tricubic_txx_

Generated on Mon Nov 23 05:15:07 2009 for contrib/mul/vil3d by  doxygen 1.5.1