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

Generated on 6 May 2016 for Hugintrunk by  doxygen 1.4.7