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 <algorithms/PanoramaAlgorithm.h>
00031 
00032 #include <boost/random.hpp>
00033 #include <vigra/impex.hxx>
00034 #include <vigra_ext/utils.h>
00035 #include <vigra_ext/Pyramid.h>
00036 #include <appbase/ProgressReporterOld.h>
00037 #include <panodata/PanoramaData.h>
00038 
00039 namespace HuginBase
00040 {
00041 
00042     class IMPEX PointSampler : public TimeConsumingPanoramaAlgorithm
00043     {
00044         protected:
00046             PointSampler(PanoramaData& panorama, AppBase::ProgressDisplay* progressDisplay,
00047                          std::vector<vigra::FRGBImage*> images,
00048                          int nPoints)
00049                 : TimeConsumingPanoramaAlgorithm(panorama, progressDisplay),
00050                   o_images(images), o_numPoints(nPoints)
00051             {};
00052         
00053         public:        
00055             virtual ~PointSampler() {};
00056         
00057         
00058         public:
00060             static void extractPoints(PanoramaData& pano, std::vector<vigra::FRGBImage*> images, int nPoints,
00061                                       bool randomPoints, AppBase::ProgressReporter& progress,
00062                                       std::vector<vigra_ext::PointPairRGB>& points);
00063         
00064         protected:
00066             typedef vigra_ext::ImageInterpolator<vigra::FRGBImage::const_traverser,
00067                                                  vigra::FRGBImage::ConstAccessor,
00068                                                  vigra_ext::interp_cubic           > InterpolImg;
00069             
00071             void sampleAndExtractPoints(AppBase::ProgressReporter& progress);
00072             
00074             virtual void samplePoints(const std::vector<InterpolImg>& imgs,
00075                                       const std::vector<vigra::FImage*>& voteImgs,
00076                                       const std::vector<SrcPanoImage>& src,
00077                                       const PanoramaOptions& dest,
00078                                       float minI,
00079                                       float maxI,
00080                                       std::vector<std::multimap<double,vigra_ext::PointPairRGB> >& radiusHist,
00081                                       unsigned& nGoodPoints,
00082                                       unsigned& nBadPoints,
00083                                       AppBase::ProgressReporter& progress) =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::ProgressReporter& progress);
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 std::vector<SrcPanoImage> & src,
00154                                      const PanoramaOptions & dest,
00155                                      int nPoints,
00156                                      float minI,
00157                                      float maxI,
00158                                      std::vector<std::multimap<double, PP > > & radiusHist,
00159                                      unsigned & nGoodPoints,
00160                                      unsigned & nBadPoints,
00161                                      AppBase::ProgressReporter & progress);
00162             
00163         protected:
00165             virtual void samplePoints(const std::vector<InterpolImg>& imgs,
00166                                       const std::vector<vigra::FImage*>& voteImgs,
00167                                       const std::vector<SrcPanoImage>& src,
00168                                       const PanoramaOptions& dest,
00169                                       float minI,
00170                                       float maxI,
00171                                       std::vector<std::multimap<double,vigra_ext::PointPairRGB> >& radiusHist,
00172                                       unsigned& nGoodPoints,
00173                                       unsigned& nBadPoints,
00174                                       AppBase::ProgressReporter& progress)
00175             {
00176                 sampleAllPanoPoints(imgs,
00177                                     voteImgs,
00178                                     src,
00179                                     dest,
00180                                      o_numPoints,
00181                                     minI,
00182                                     maxI,
00183                                     radiusHist,
00184                                     nGoodPoints,
00185                                     nBadPoints,
00186                                     progress);
00187             }
00188             
00189     };
00190 
00191 
00195     class RandomPointSampler : public  PointSampler
00196     {
00197         public:
00199             RandomPointSampler(PanoramaData& panorama, AppBase::ProgressDisplay* progressDisplay,
00200                                std::vector<vigra::FRGBImage*> images,
00201                                int nPoints)
00202              : PointSampler(panorama, progressDisplay, images, nPoints)
00203             {};
00204 
00206             virtual ~RandomPointSampler() {};
00207             
00208             
00209         public:
00210             template <class Img, class VoteImg, class PP>
00211             static void sampleRandomPanoPoints(const std::vector<Img> imgs,
00212                                                const std::vector<VoteImg *> &voteImgs,
00213                                                const std::vector<SrcPanoImage> & src,
00214                                                const PanoramaOptions & dest,
00215                                                int nPoints,
00216                                                float minI,
00217                                                float maxI,
00218                                                std::vector<std::multimap<double, PP > > & radiusHist,
00219                                                unsigned & nBadPoints,
00220                                                AppBase::ProgressReporter & progress);
00221             
00222         protected:
00224             virtual void samplePoints(const std::vector<InterpolImg>& imgs,
00225                                       const std::vector<vigra::FImage*>& voteImgs,
00226                                       const std::vector<SrcPanoImage>& src,
00227                                       const PanoramaOptions& dest,
00228                                       float minI,
00229                                       float maxI,
00230                                       std::vector<std::multimap<double,vigra_ext::PointPairRGB> >& radiusHist,
00231                                       unsigned& nGoodPoints,
00232                                       unsigned& nBadPoints,
00233                                       AppBase::ProgressReporter& progress)
00234             {
00235                  sampleRandomPanoPoints(imgs,
00236                                         voteImgs,
00237                                         src,
00238                                         dest,
00239                                          5*o_numPoints,
00240                                         minI,
00241                                         maxI,
00242                                         radiusHist,
00243                                         nBadPoints,
00244                                         progress);
00245             }
00246             
00247     };
00248 
00249     
00250     
00251 } // namespace
00252     
00253     
00254 
00255 
00256     
00257 //==============================================================================
00258 //  templated methods
00259 
00260 
00261 #include <boost/random.hpp>
00262 #include <panotools/PanoToolsInterface.h>
00263 
00264 
00265 namespace HuginBase {
00266 
00267 
00268 template<class PointPairClass>
00269 void PointSampler::sampleRadiusUniform(const std::vector<std::multimap<double, PointPairClass> >& radiusHist,
00270                                        unsigned nPoints,
00271                                        std::vector<PointPairClass> &selectedPoints,
00272                                        AppBase::ProgressReporter & progress)
00273 {
00274     typedef std::multimap<double,PointPairClass> PointPairMap;
00275     
00276     // reserve selected points..
00277     int drawsPerBin = nPoints / radiusHist.size();
00278     // reallocate output vector.
00279     selectedPoints.reserve(drawsPerBin*radiusHist.size());
00280     for (typename std::vector<PointPairMap>::const_iterator bin = radiusHist.begin();
00281          bin != radiusHist.end(); ++bin) 
00282     {
00283         unsigned i=drawsPerBin;
00284         // copy into vector
00285         for (typename PointPairMap::const_iterator it= (*bin).begin();
00286              it != (*bin).end(); ++it) 
00287         {
00288             selectedPoints.push_back(it->second);
00289             // do not extract more than drawsPerBin Point pairs.
00290             --i;
00291             if (i == 0)
00292                 break;
00293         }
00294         progress.increaseProgress(1.0/radiusHist.size());
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::ProgressReporter& 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->size() == 0) {
00407                                 destMap = map1;
00408                             } else if (map2->size() == 0) {
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         int h = dest.getROI().bottom() - dest.getROI().top();
00434         if ((y-dest.getROI().top())%(h/10) == 0) {
00435             progress.increaseProgress(0.1);
00436         }
00437     }
00438 
00439     for(unsigned i=0; i < imgs.size(); i++) {
00440         delete transf[i];
00441     }
00442 }
00443 
00444 
00445 
00446 template <class Img, class VoteImg, class PP>
00447 void RandomPointSampler::sampleRandomPanoPoints(const std::vector<Img> imgs,
00448                                                 const std::vector<VoteImg *> &voteImgs,
00449                                                 const std::vector<SrcPanoImage> & src,
00450                                                 const PanoramaOptions & dest,
00451                                                 int nPoints,
00452                                                 float minI,
00453                                                 float maxI,
00454                                                 //std::vector<PP> &points,
00455                                                 std::vector<std::multimap<double, PP > > & radiusHist,
00456                                                 unsigned & nBadPoints,
00457                                                 AppBase::ProgressReporter & progress)
00458 {
00459     typedef typename Img::PixelType PixelType;
00460 
00461     vigra_precondition(imgs.size() > 1, "sampleRandomPanoPoints: At least two images required");
00462     vigra_precondition(imgs.size() == src.size(), "number of src images doesn't match");
00463     
00464     unsigned nImg = imgs.size();
00465 
00466     unsigned nBins = radiusHist.size();
00467     unsigned pairsPerBin = nPoints / nBins;
00468 
00469     int allPoints = nPoints;
00470 
00471     // create an array of transforms.
00472     //std::vector<SpaceTransform> transf(imgs.size());
00473     std::vector<PTools::Transform *> transf(imgs.size());
00474     std::vector<double> maxr(imgs.size());
00475 
00476     // initialize transforms, and interpolating accessors
00477     for(unsigned i=0; i < imgs.size(); i++) {
00478         // same size is not needed?
00479 //        vigra_precondition(src[i].getSize() == srcSize, "images need to have the same size");
00480         transf[i] = new PTools::Transform;
00481         transf[i]->createTransform(src[i], dest);
00482         vigra::Size2D srcSize = src[i].getSize();
00483         maxr[i] = sqrt(((double)srcSize.x)*srcSize.x + ((double)srcSize.y)*srcSize.y) / 2.0;
00484     }
00485     // init random number generator
00486     boost::mt19937 rng;
00487     // start with a different seed every time.
00488     rng.seed(static_cast<unsigned int>(std::time(0)));
00489     // randomly sample points.
00490     boost::uniform_int<> distribx(dest.getROI().left(), dest.getROI().right()-1);
00491     boost::uniform_int<> distriby(dest.getROI().top(), dest.getROI().bottom()-1);
00492 
00493     boost::variate_generator<boost::mt19937&, boost::uniform_int<> >
00494             randX(rng, distribx);             // glues randomness with mapping
00495     boost::variate_generator<boost::mt19937&, boost::uniform_int<> >
00496             randY(rng, distriby);             // glues randomness with mapping
00497 
00498     double percentReported = 0.0;
00499     
00500     for (unsigned maxTry = nPoints*5; nPoints > 0 && maxTry > 0; maxTry--) {
00501         unsigned x = randX();
00502         unsigned y = randY();
00503         hugin_utils::FDiff2D panoPnt(x,y);
00504         for (unsigned i=0; i< nImg; i++) {
00505             // transform pixel
00506             PixelType i1;
00507             hugin_utils::FDiff2D p1;
00508             if(!transf[i]->transformImgCoord(p1, panoPnt))
00509                 continue;
00510             vigra::Point2D p1Int(p1.toDiff2D());
00511             // check if pixel is valid
00512             if (!src[i].isInside(p1Int)) {
00513                 // point is outside image
00514                 continue;
00515             }
00516             vigra::UInt8 maskI;
00517             if ( imgs[i](p1.x,p1.y, i1, maskI)){
00518                 float im1 = vigra_ext::getMaxComponent(i1);
00519                 if (minI > im1 || maxI < im1) {
00520                     // ignore pixels that are too dark or bright
00521                     continue;
00522                 }
00523                 double r1 = hugin_utils::norm((p1 - src[i].getRadialVigCorrCenter())/maxr[i]);
00524                 for (unsigned j=i+1; j < nImg; j++) {
00525                     PixelType i2;
00526                     hugin_utils::FDiff2D p2;
00527                     if(!transf[j]->transformImgCoord(p2, panoPnt))
00528                         continue;
00529                     // check if a pixel is inside the source image
00530                     vigra::Point2D p2Int(p2.toDiff2D());
00531                     if (!src[j].isInside(p2Int)) {
00532                         // point is outside image
00533                         continue;
00534                     }
00535                     vigra::UInt8 maskI2;
00536                     if (imgs[j](p2.x, p2.y, i2, maskI2)){
00537                         float im2 = vigra_ext::getMaxComponent(i2);
00538                         if (minI > im2 || maxI < im2) {
00539                             // ignore pixels that are too dark or bright
00540                             continue;
00541                         }
00542                         // TODO: add check for gradient radius.
00543                         double r2 = hugin_utils::norm((p2 - src[j].getRadialVigCorrCenter())/maxr[j]);
00544 #if 0
00545                         // add pixel
00546                         if (im1 <= im2) {
00547                             points.push_back(PP(i, i1, p1, r1,   j, i2, p2, r2) );
00548                         } else {
00549                             points.push_back(PP(j, i2, p2, r2,   i, i1, p1, r1) );
00550                         }
00551 #else
00552                             // add pixel
00553                             const VoteImg & vimg1 =  *voteImgs[i];
00554                             const VoteImg & vimg2 =  *voteImgs[j];
00555                             double laplace = hugin_utils::sqr(vimg1[p1Int]) + hugin_utils::sqr(vimg2[p2Int]);
00556                             size_t bin1 = (size_t)(r1*nBins);
00557                             size_t bin2 = (size_t)(r2*nBins);
00558                             // a center shift might lead to radi > 1.
00559                             if (bin1+1 > nBins) bin1 = nBins-1;
00560                             if (bin2+1 > nBins) bin2 = nBins-1;
00561 
00562                             PP pp;
00563                             if (im1 <= im2) {
00564                                 // choose i1 to be smaller than i2
00565                                 pp = PP(i, i1, p1, r1,   j, i2, p2, r2);
00566                             } else {
00567                                 pp = PP(j, i2, p2, r2,   i, i1, p1, r1);
00568                             }
00569 
00570                             // decide which bin should be used.
00571                             std::multimap<double, PP> * map1 = &radiusHist[bin1];
00572                             std::multimap<double, PP> * map2 = &radiusHist[bin2];
00573                             std::multimap<double, PP> * destMap;
00574                             if (map1->size() == 0) {
00575                                 destMap = map1;
00576                             } else if (map2->size() == 0) {
00577                                 destMap = map2;
00578                             } else if (map1->size() < map2->size()) {
00579                                 destMap = map1;
00580                             } else if (map1->size() > map2->size()) {
00581                                 destMap = map2;
00582                                                         } else if (map1->rbegin()->first > map2->rbegin()->first) {
00583                                 // heuristic: insert into bin with higher maximum laplacian filter response
00584                                 // (higher probablity of misregistration).
00585                                 destMap = map1;
00586                             } else {
00587                                 destMap = map2;
00588                             }
00589                             // insert
00590                             destMap->insert(std::make_pair(laplace,pp));
00591                             // remove last element if too many elements have been gathered
00592                             if (destMap->size() > pairsPerBin) {
00593                                 destMap->erase((--(destMap->end())));
00594                             }
00595 //                            nGoodPoints++;
00596 #endif
00597                         nPoints--;
00598                     }
00599                 }
00600             }
00601         }
00602         double pc = (allPoints - nPoints);
00603         double percentNow = (pc / allPoints) * 100.0;
00604         if (percentNow - percentReported >= 10) {
00605             percentReported = percentNow;
00606             progress.increaseProgress(0.1);
00607         }
00608     }
00609     for(unsigned i=0; i < imgs.size(); i++) {
00610         delete transf[i];
00611     }
00612     
00613     DEBUG_INFO("Point sampled: " << allPoints-nPoints)
00614     progress.increaseProgress(1.0 - percentReported/100.0);
00615 }
00616 
00617 
00618 } // namespace
00619 #endif //_H

Generated on Thu Jul 31 01:25:40 2014 for Hugintrunk by  doxygen 1.3.9.1