PointSampler.h

Go to the documentation of this file.
00001 // -*- c-basic-offset: 4 -*-
00024 #ifndef _POINTSAMPLER_H
00025 #define _POINTSAMPLER_H
00026 
00027 #include <ctime>
00028 
00029 #include <hugin_shared.h>
00030 #include <hugin_config.h>
00031 #include <vigra/stdimage.hxx>
00032 #include <algorithms/PanoramaAlgorithm.h>
00033 
00034 #ifdef HAVE_CXX11
00035 #include <random>
00036 #include <functional>
00037 #else
00038 #include <boost/random.hpp>
00039 #endif
00040 #include <vigra_ext/utils.h>
00041 #include <appbase/ProgressDisplay.h>
00042 #include <panodata/PanoramaData.h>
00043 
00044 namespace HuginBase
00045 {
00046 
00047     class IMPEX PointSampler : public TimeConsumingPanoramaAlgorithm
00048     {
00049         protected:
00051             PointSampler(PanoramaData& panorama, AppBase::ProgressDisplay* progressDisplay,
00052                          std::vector<vigra::FRGBImage*> images,
00053                          int nPoints)
00054                 : TimeConsumingPanoramaAlgorithm(panorama, progressDisplay),
00055                   o_images(images), o_numPoints(nPoints)
00056             {};
00057         
00058         public:        
00060             virtual ~PointSampler() {};
00061         
00062         
00063         public:
00065             static void extractPoints(PanoramaData& pano, std::vector<vigra::FRGBImage*> images, int nPoints,
00066                                       bool randomPoints, AppBase::ProgressDisplay* progress,
00067                                       std::vector<vigra_ext::PointPairRGB>& points);
00068         
00069         protected:
00071             typedef vigra_ext::ImageInterpolator<vigra::FRGBImage::const_traverser,
00072                                                  vigra::FRGBImage::ConstAccessor,
00073                                                  vigra_ext::interp_cubic           > InterpolImg;
00074             
00076             void sampleAndExtractPoints(AppBase::ProgressDisplay* progress);
00077             
00079             virtual void samplePoints(const std::vector<InterpolImg>& imgs,
00080                                       const std::vector<vigra::FImage*>& voteImgs,
00081                                       const PanoramaData& pano,
00082                                       float minI,
00083                                       float maxI,
00084                                       std::vector<std::multimap<double,vigra_ext::PointPairRGB> >& radiusHist,
00085                                       unsigned& nGoodPoints,
00086                                       unsigned& nBadPoints,
00087                                       AppBase::ProgressDisplay*) = 0;
00088             
00093             template<class PointPairClass>
00094             static void sampleRadiusUniform(const std::vector<std::multimap<double,PointPairClass> >& radiusHist,
00095                                             unsigned nPoints,
00096                                             std::vector<PointPairClass>& selectedPoints,
00097                                             AppBase::ProgressDisplay*);
00098             
00099         public:
00101             virtual bool modifiesPanoramaData() const
00102                 { return false; }
00103             
00105             virtual bool runAlgorithm();
00106 
00107             PointSampler & execute()
00108             {
00109                                 run();
00110                 return *this;
00111             }
00112             
00113         public:
00115             typedef std::vector<vigra_ext::PointPairRGB> PointPairs;
00116             
00118             PointPairs getResultPoints()
00119             {
00120                 //[TODO] debug unsuccessful
00121                 
00122                 return o_resultPoints;
00123             }
00124             
00125             
00126         protected:
00127             std::vector<vigra::FRGBImage*> o_images;
00128             int o_numPoints;
00129             PointPairs o_resultPoints;
00130     };
00131 
00132 
00136     class AllPointSampler : public  PointSampler
00137     {
00138         public:
00140             AllPointSampler(PanoramaData& panorama, AppBase::ProgressDisplay* progressDisplay,
00141                                std::vector<vigra::FRGBImage*> images,
00142                                int nPoints)
00143              : PointSampler(panorama, progressDisplay, images, nPoints)
00144             {};
00145 
00147             virtual ~AllPointSampler() {};
00148             
00149             
00150         public:
00154             template <class Img, class VoteImg, class PP>
00155             static void sampleAllPanoPoints(const std::vector<Img> &imgs,
00156                                      const std::vector<VoteImg *> &voteImgs,
00157                                      const PanoramaData& pano,
00158                                      int nPoints,
00159                                      float minI,
00160                                      float maxI,
00161                                      std::vector<std::multimap<double, PP > > & radiusHist,
00162                                      unsigned & nGoodPoints,
00163                                      unsigned & nBadPoints,
00164                                      AppBase::ProgressDisplay* progress);
00165             
00166         protected:
00168             virtual void samplePoints(const std::vector<InterpolImg>& imgs,
00169                                       const std::vector<vigra::FImage*>& voteImgs,
00170                                       const PanoramaData& pano,
00171                                       float minI,
00172                                       float maxI,
00173                                       std::vector<std::multimap<double,vigra_ext::PointPairRGB> >& radiusHist,
00174                                       unsigned& nGoodPoints,
00175                                       unsigned& nBadPoints,
00176                                       AppBase::ProgressDisplay* progress)
00177             {
00178                 sampleAllPanoPoints(imgs,
00179                                     voteImgs,
00180                                     pano,
00181                                     o_numPoints,
00182                                     minI,
00183                                     maxI,
00184                                     radiusHist,
00185                                     nGoodPoints,
00186                                     nBadPoints,
00187                                     progress);
00188             }
00189             
00190     };
00191 
00192 
00196     class RandomPointSampler : public  PointSampler
00197     {
00198         public:
00200             RandomPointSampler(PanoramaData& panorama, AppBase::ProgressDisplay* progressDisplay,
00201                                std::vector<vigra::FRGBImage*> images,
00202                                int nPoints)
00203              : PointSampler(panorama, progressDisplay, images, nPoints)
00204             {};
00205 
00207             virtual ~RandomPointSampler() {};
00208             
00209             
00210         public:
00211             template <class Img, class VoteImg, class PP>
00212             static void sampleRandomPanoPoints(const std::vector<Img>& imgs,
00213                                                const std::vector<VoteImg *> &voteImgs,
00214                                                const PanoramaData& pano,
00215                                                int nPoints,
00216                                                float minI,
00217                                                float maxI,
00218                                                std::vector<std::multimap<double, PP > > & radiusHist,
00219                                                unsigned & nBadPoints,
00220                                                AppBase::ProgressDisplay* progress);
00221             
00222         protected:
00224             virtual void samplePoints(const std::vector<InterpolImg>& imgs,
00225                                       const std::vector<vigra::FImage*>& voteImgs,
00226                                       const PanoramaData& pano,
00227                                       float minI,
00228                                       float maxI,
00229                                       std::vector<std::multimap<double,vigra_ext::PointPairRGB> >& radiusHist,
00230                                       unsigned& nGoodPoints,
00231                                       unsigned& nBadPoints,
00232                                       AppBase::ProgressDisplay* progress)
00233             {
00234                  sampleRandomPanoPoints(imgs,
00235                                         voteImgs,
00236                                         pano,
00237                                         5*o_numPoints,
00238                                         minI,
00239                                         maxI,
00240                                         radiusHist,
00241                                         nBadPoints,
00242                                         progress);
00243             }
00244             
00245     };
00246 
00247     
00248     
00249 } // namespace
00250     
00251     
00252 
00253 
00254     
00255 //==============================================================================
00256 //  templated methods
00257 
00258 #include <panotools/PanoToolsInterface.h>
00259 
00260 namespace HuginBase {
00261 
00262 template<class PointPairClass>
00263 void PointSampler::sampleRadiusUniform(const std::vector<std::multimap<double, PointPairClass> >& radiusHist,
00264                                        unsigned nPoints,
00265                                        std::vector<PointPairClass> &selectedPoints,
00266                                        AppBase::ProgressDisplay* progress)
00267 {
00268     typedef std::multimap<double,PointPairClass> PointPairMap;
00269     
00270     // reserve selected points..
00271     int drawsPerBin = nPoints / radiusHist.size();
00272     // reallocate output vector.
00273     selectedPoints.reserve(drawsPerBin*radiusHist.size());
00274     for (typename std::vector<PointPairMap>::const_iterator bin = radiusHist.begin();
00275          bin != radiusHist.end(); ++bin) 
00276     {
00277         unsigned i=drawsPerBin;
00278         // copy into vector
00279         for (typename PointPairMap::const_iterator it= (*bin).begin();
00280              it != (*bin).end(); ++it) 
00281         {
00282             selectedPoints.push_back(it->second);
00283             // do not extract more than drawsPerBin Point pairs.
00284             --i;
00285             if (i == 0)
00286                 break;
00287         }
00288     }
00289 }
00290     
00291     
00292     
00293 template <class Img, class VoteImg, class PP>
00294 void AllPointSampler::sampleAllPanoPoints(const std::vector<Img> &imgs,
00295                                           const std::vector<VoteImg *> &voteImgs,
00296                                           const PanoramaData& pano,
00297                                           int nPoints,
00298                                           float minI,
00299                                           float maxI,
00300                                           //std::vector<vigra_ext::PointPair> &points,
00301                                           std::vector<std::multimap<double, PP > > & radiusHist,
00302                                           unsigned & nGoodPoints,
00303                                           unsigned & nBadPoints,
00304                                           AppBase::ProgressDisplay* progress)
00305 {
00306     typedef typename Img::PixelType PixelType;
00307 
00308     // use 10 bins
00309     radiusHist.resize(10);
00310     const unsigned pairsPerBin = nPoints / radiusHist.size();
00311 
00312     nGoodPoints = 0;
00313     nBadPoints = 0;
00314     vigra_precondition(imgs.size() > 1, "sampleAllPanoPoints: At least two images required");
00315     
00316     const unsigned nImg = imgs.size();
00317 
00318     vigra::Size2D srcSize = pano.getImage(0).getSize();
00319     double maxr = sqrt(((double)srcSize.x)*srcSize.x + ((double)srcSize.y)*srcSize.y) / 2.0;
00320 
00321     // create an array of transforms.
00322     //std::vector<SpaceTransform> transf(imgs.size());
00323     std::vector<PTools::Transform*> transf(imgs.size());
00324 
00325     // initialize transforms, and interpolating accessors
00326     for(unsigned i=0; i < nImg; i++) {
00327         vigra_precondition(pano.getImage(i).getSize() == srcSize, "images need to have the same size");
00328         transf[i] = new PTools::Transform;
00329         transf[i]->createTransform(pano.getImage(i), pano.getOptions());
00330     }
00331 
00332     const vigra::Rect2D roi = pano.getOptions().getROI();
00333     for (int y=roi.top(); y < roi.bottom(); ++y) {
00334         for (int x=roi.left(); x < roi.right(); ++x) {
00335             hugin_utils::FDiff2D panoPnt(x,y);
00336             for (unsigned i=0; i< nImg-1; i++) {
00337             // transform pixel
00338                 hugin_utils::FDiff2D p1;
00339                 if(!transf[i]->transformImgCoord(p1, panoPnt))
00340                     continue;
00341                 vigra::Point2D p1Int(p1.toDiff2D());
00342                 // is inside:
00343                 if (!pano.getImage(i).isInside(p1Int)) {
00344                     // point is outside image
00345                     continue;
00346                 }
00347                 PixelType i1;
00348                 vigra::UInt8 maskI;
00349                 if (imgs[i](p1.x,p1.y, i1, maskI)){
00350                     float im1 = vigra_ext::getMaxComponent(i1);
00351                     if (minI > im1 || maxI < im1 || maskI == 0) {
00352                         // ignore pixels that are too dark or bright
00353                         continue;
00354                     }
00355                     double r1 = hugin_utils::norm((p1 - pano.getImage(i).getRadialVigCorrCenter()) / maxr);
00356 
00357                     // check inner image
00358                     for (unsigned j=i+1; j < nImg; j++) {
00359                         hugin_utils::FDiff2D p2;
00360                         if(!transf[j]->transformImgCoord(p2, panoPnt))
00361                             continue;
00362                         vigra::Point2D p2Int(p2.toDiff2D());
00363                         if (!pano.getImage(j).isInside(p2Int)) {
00364                             // point is outside image
00365                             continue;
00366                         }
00367                         PixelType i2;
00368                         vigra::UInt8 maskI2;
00369                         if (imgs[j](p2.x, p2.y, i2, maskI2)){
00370                             float im2 = vigra_ext::getMaxComponent(i2);
00371                             if (minI > im2 || maxI < im2 || maskI2 == 0) {
00372                                 // ignore pixels that are too dark or bright
00373                                 continue;
00374                             }
00375                             double r2 = hugin_utils::norm((p2 - pano.getImage(j).getRadialVigCorrCenter()) / maxr);
00376                             // add pixel
00377                             const VoteImg & vimg1 =  *voteImgs[i];
00378                             const VoteImg & vimg2 =  *voteImgs[j];
00379                             double laplace = hugin_utils::sqr(vimg1[p1Int]) + hugin_utils::sqr(vimg2[p2Int]);
00380                             int bin1 = (int)(r1*10);
00381                             int bin2 = (int)(r2*10);
00382                             // a center shift might lead to radi > 1.
00383                             if (bin1 > 9) bin1 = 9;
00384                             if (bin2 > 9) bin2 = 9;
00385 
00386                             PP pp;
00387                             if (im1 <= im2) {
00388                                 // choose i1 to be smaller than i2
00389                                 pp = PP(i, i1, p1, r1,   j, i2, p2, r2);
00390                             } else {
00391                                 pp = PP(j, i2, p2, r2,   i, i1, p1, r1);
00392                             }
00393 
00394                             // decide which bin should be used.
00395                             std::multimap<double, PP> * map1 = &radiusHist[bin1];
00396                             std::multimap<double, PP> * map2 = &radiusHist[bin2];
00397                             std::multimap<double, PP> * destMap;
00398                             if (map1->empty()) {
00399                                 destMap = map1;
00400                             } else if (map2->empty()) {
00401                                 destMap = map2;
00402                             } else if (map1->size() < map2->size()) {
00403                                 destMap = map1;
00404                             } else if (map1->size() > map2->size()) {
00405                                 destMap = map2;
00406                             } else if (map1->rbegin()->first > map2->rbegin()->first) {
00407                                 // heuristic: insert into bin with higher maximum laplacian filter response
00408                                 // (higher probablity of misregistration).
00409                                 destMap = map1;
00410                             } else {
00411                                 destMap = map2;
00412                             }
00413                             // insert
00414                             destMap->insert(std::make_pair(laplace,pp));
00415                             // remove last element if too many elements have been gathered
00416                             if (destMap->size() > pairsPerBin) {
00417                                 destMap->erase((--(destMap->end())));
00418                             }
00419                             nGoodPoints++;
00420                         }
00421                     }
00422                 }
00423             }
00424         }
00425     }
00426 
00427     for(unsigned i=0; i < nImg; i++) {
00428         delete transf[i];
00429     }
00430 }
00431 
00432 
00433 
00434 template <class Img, class VoteImg, class PP>
00435 void RandomPointSampler::sampleRandomPanoPoints(const std::vector<Img>& imgs,
00436                                                 const std::vector<VoteImg *> &voteImgs,
00437                                                 const PanoramaData& pano,
00438                                                 int nPoints,
00439                                                 float minI,
00440                                                 float maxI,
00441                                                 //std::vector<PP> &points,
00442                                                 std::vector<std::multimap<double, PP > > & radiusHist,
00443                                                 unsigned & nBadPoints,
00444                                                 AppBase::ProgressDisplay* progress)
00445 {
00446     typedef typename Img::PixelType PixelType;
00447 
00448     vigra_precondition(imgs.size() > 1, "sampleRandomPanoPoints: At least two images required");
00449     
00450     const unsigned nImg = imgs.size();
00451 
00452     const unsigned nBins = radiusHist.size();
00453     const unsigned pairsPerBin = nPoints / nBins;
00454 
00455     int allPoints = nPoints;
00456 
00457     // create an array of transforms.
00458     //std::vector<SpaceTransform> transf(imgs.size());
00459     std::vector<PTools::Transform *> transf(imgs.size());
00460     std::vector<double> maxr(imgs.size());
00461 
00462     // initialize transforms, and interpolating accessors
00463     for(unsigned i=0; i < nImg; i++) {
00464         // same size is not needed?
00465 //        vigra_precondition(src[i].getSize() == srcSize, "images need to have the same size");
00466         transf[i] = new PTools::Transform;
00467         transf[i]->createTransform(pano.getImage(i), pano.getOptions());
00468         vigra::Size2D srcSize = pano.getImage(i).getSize();
00469         maxr[i] = sqrt(((double)srcSize.x)*srcSize.x + ((double)srcSize.y)*srcSize.y) / 2.0;
00470     }
00471     // init random number generator
00472     const vigra::Rect2D roi = pano.getOptions().getROI();
00473 #ifdef HAVE_CXX11
00474     std::mt19937 rng(static_cast<unsigned int>(std::time(0)));
00475     std::uniform_int_distribution<unsigned int> distribx(roi.left(), roi.right()-1);
00476     std::uniform_int_distribution<unsigned int> distriby(roi.top(), roi.bottom()-1);
00477 
00478     auto randX = std::bind(distribx, std::ref(rng));
00479     auto randY = std::bind(distriby, std::ref(rng));
00480 #else
00481     boost::mt19937 rng;
00482     // start with a different seed every time.
00483     rng.seed(static_cast<unsigned int>(std::time(0)));
00484     // randomly sample points.
00485     boost::uniform_int<> distribx(roi.left(), roi.right()-1);
00486     boost::uniform_int<> distriby(roi.top(), roi.bottom()-1);
00487 
00488     boost::variate_generator<boost::mt19937&, boost::uniform_int<> >
00489             randX(rng, distribx);             // glues randomness with mapping
00490     boost::variate_generator<boost::mt19937&, boost::uniform_int<> >
00491             randY(rng, distriby);             // glues randomness with mapping
00492 #endif
00493 
00494     for (unsigned maxTry = nPoints*5; nPoints > 0 && maxTry > 0; maxTry--) {
00495         unsigned x = randX();
00496         unsigned y = randY();
00497         hugin_utils::FDiff2D panoPnt(x,y);
00498         for (unsigned i=0; i< nImg-1; i++) {
00499             // transform pixel
00500             PixelType i1;
00501             hugin_utils::FDiff2D p1;
00502             if(!transf[i]->transformImgCoord(p1, panoPnt))
00503                 continue;
00504             vigra::Point2D p1Int(p1.toDiff2D());
00505             // check if pixel is valid
00506             if (!pano.getImage(i).isInside(p1Int)) {
00507                 // point is outside image
00508                 continue;
00509             }
00510             vigra::UInt8 maskI;
00511             if ( imgs[i](p1.x,p1.y, i1, maskI)){
00512                 float im1 = vigra_ext::getMaxComponent(i1);
00513                 if (minI > im1 || maxI < im1) {
00514                     // ignore pixels that are too dark or bright
00515                     continue;
00516                 }
00517                 double r1 = hugin_utils::norm((p1 - pano.getImage(i).getRadialVigCorrCenter()) / maxr[i]);
00518                 for (unsigned j=i+1; j < nImg; j++) {
00519                     PixelType i2;
00520                     hugin_utils::FDiff2D p2;
00521                     if(!transf[j]->transformImgCoord(p2, panoPnt))
00522                         continue;
00523                     // check if a pixel is inside the source image
00524                     vigra::Point2D p2Int(p2.toDiff2D());
00525                     if (!pano.getImage(j).isInside(p2Int)) {
00526                         // point is outside image
00527                         continue;
00528                     }
00529                     vigra::UInt8 maskI2;
00530                     if (imgs[j](p2.x, p2.y, i2, maskI2)){
00531                         float im2 = vigra_ext::getMaxComponent(i2);
00532                         if (minI > im2 || maxI < im2) {
00533                             // ignore pixels that are too dark or bright
00534                             continue;
00535                         }
00536                         // TODO: add check for gradient radius.
00537                         double r2 = hugin_utils::norm((p2 - pano.getImage(j).getRadialVigCorrCenter()) / maxr[j]);
00538 #if 0
00539                         // add pixel
00540                         if (im1 <= im2) {
00541                             points.push_back(PP(i, i1, p1, r1,   j, i2, p2, r2) );
00542                         } else {
00543                             points.push_back(PP(j, i2, p2, r2,   i, i1, p1, r1) );
00544                         }
00545 #else
00546                             // add pixel
00547                             const VoteImg & vimg1 =  *voteImgs[i];
00548                             const VoteImg & vimg2 =  *voteImgs[j];
00549                             double laplace = hugin_utils::sqr(vimg1[p1Int]) + hugin_utils::sqr(vimg2[p2Int]);
00550                             size_t bin1 = (size_t)(r1*nBins);
00551                             size_t bin2 = (size_t)(r2*nBins);
00552                             // a center shift might lead to radi > 1.
00553                             if (bin1+1 > nBins) bin1 = nBins-1;
00554                             if (bin2+1 > nBins) bin2 = nBins-1;
00555 
00556                             PP pp;
00557                             if (im1 <= im2) {
00558                                 // choose i1 to be smaller than i2
00559                                 pp = PP(i, i1, p1, r1,   j, i2, p2, r2);
00560                             } else {
00561                                 pp = PP(j, i2, p2, r2,   i, i1, p1, r1);
00562                             }
00563 
00564                             // decide which bin should be used.
00565                             std::multimap<double, PP> * map1 = &radiusHist[bin1];
00566                             std::multimap<double, PP> * map2 = &radiusHist[bin2];
00567                             std::multimap<double, PP> * destMap;
00568                             if (map1->empty()) {
00569                                 destMap = map1;
00570                             } else if (map2->empty()) {
00571                                 destMap = map2;
00572                             } else if (map1->size() < map2->size()) {
00573                                 destMap = map1;
00574                             } else if (map1->size() > map2->size()) {
00575                                 destMap = map2;
00576                                                         } else if (map1->rbegin()->first > map2->rbegin()->first) {
00577                                 // heuristic: insert into bin with higher maximum laplacian filter response
00578                                 // (higher probablity of misregistration).
00579                                 destMap = map1;
00580                             } else {
00581                                 destMap = map2;
00582                             }
00583                             // insert
00584                             destMap->insert(std::make_pair(laplace,pp));
00585                             // remove last element if too many elements have been gathered
00586                             if (destMap->size() > pairsPerBin) {
00587                                 destMap->erase((--(destMap->end())));
00588                             }
00589 //                            nGoodPoints++;
00590 #endif
00591                         nPoints--;
00592                     }
00593                 }
00594             }
00595         }
00596     }
00597     for(unsigned i=0; i < imgs.size(); i++) {
00598         delete transf[i];
00599     }
00600     
00601     DEBUG_INFO("Point sampled: " << allPoints-nPoints)
00602 }
00603 
00604 
00605 } // namespace
00606 #endif //_H

Generated on 9 Feb 2016 for Hugintrunk by  doxygen 1.4.7