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

Generated on 29 Jul 2015 for Hugintrunk by  doxygen 1.4.7