PanoDetectorLogic.cpp

Go to the documentation of this file.
00001 // -*- c-basic-offset: 4 ; tab-width: 4 -*-
00002 /*
00003 * Copyright (C) 2007-2008 Anael Orlinski
00004 *
00005 * This file is part of Panomatic.
00006 *
00007 * Panomatic is free software; you can redistribute it and/or modify
00008 * it under the terms of the GNU General Public License as published by
00009 * the Free Software Foundation; either version 2 of the License, or
00010 * (at your option) any later version.
00011 *
00012 * Panomatic is distributed in the hope that it will be useful,
00013 * but WITHOUT ANY WARRANTY; without even the implied warranty of
00014 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
00015 * GNU General Public License for more details.
00016 *
00017 * You should have received a copy of the GNU General Public License
00018 * along with Panomatic; if not, write to the Free Software
00019 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
00020 */
00021 
00022 #include "ImageImport.h"
00023 
00024 #include "PanoDetector.h"
00025 #include <iostream>
00026 #include <fstream>
00027 #include <boost/foreach.hpp>
00028 #include <vigra/distancetransform.hxx>
00029 #include "vigra/impex.hxx"   // debug image save
00030 
00031 #include <localfeatures/Sieve.h>
00032 #include <localfeatures/PointMatch.h>
00033 #include <localfeatures/RansacFiltering.h>
00034 #include <localfeatures/KeyPointIO.h>
00035 #include <localfeatures/CircularKeyPointDescriptor.h>
00036 
00037 /*
00038 #include "KDTree.h"
00039 #include "KDTreeImpl.h"
00040 */
00041 #include "Utils.h"
00042 #include "Tracer.h"
00043 
00044 #include <algorithms/nona/ComputeImageROI.h>
00045 #include <algorithms/optimizer/PTOptimizer.h>
00046 #include <nona/RemappedPanoImage.h>
00047 #include <nona/ImageRemapper.h>
00048 
00049 #include <time.h>
00050 
00051 #define TRACE_IMG(X) {if (iPanoDetector.getVerbose() > 1) { TRACE_INFO("i" << ioImgInfo._number << " : " << X << endl);} }
00052 #define TRACE_PAIR(X) {if (iPanoDetector.getVerbose() > 1){ TRACE_INFO("i" << ioMatchData._i1->_number << " <> " \
00053                 "i" << ioMatchData._i2->_number << " : " << X << endl)}}
00054 
00055 using namespace std;
00056 using namespace lfeat;
00057 using namespace HuginBase;
00058 using namespace AppBase;
00059 using namespace HuginBase::Nona;
00060 using namespace hugin_utils;
00061 
00062 
00063 static ZThread::FastMutex aPanoToolsMutex;
00064 
00065 // define a Keypoint insertor
00066 class KeyPointVectInsertor : public lfeat::KeyPointInsertor
00067 {
00068 public:
00069     KeyPointVectInsertor(KeyPointVect_t& iVect) : _v(iVect) {};
00070     inline virtual void operator()(const lfeat::KeyPoint& k)
00071     {
00072         _v.push_back(KeyPointPtr(new lfeat::KeyPoint(k)));
00073     }
00074 
00075 private:
00076     KeyPointVect_t& _v;
00077 
00078 };
00079 
00080 
00081 // define a sieve extractor
00082 class SieveExtractorKP : public lfeat::SieveExtractor<KeyPointPtr>
00083 {
00084 public:
00085     SieveExtractorKP(KeyPointVect_t& iV) : _v(iV) {};
00086     inline virtual void operator()(const KeyPointPtr& k)
00087     {
00088         _v.push_back(k);
00089     }
00090 private:
00091     KeyPointVect_t& _v;
00092 };
00093 
00094 class SieveExtractorMatch : public lfeat::SieveExtractor<lfeat::PointMatchPtr>
00095 {
00096 public:
00097     SieveExtractorMatch(lfeat::PointMatchVector_t& iM) : _m(iM) {};
00098     inline virtual void operator()(const lfeat::PointMatchPtr& m)
00099     {
00100         _m.push_back(m);
00101     }
00102 private:
00103     lfeat::PointMatchVector_t& _m;
00104 };
00105 
00106 bool PanoDetector::LoadKeypoints(ImgData& ioImgInfo, const PanoDetector& iPanoDetector)
00107 {
00108     TRACE_IMG("Loading keypoints...");
00109 
00110     ImageInfo info = lfeat::loadKeypoints(ioImgInfo._keyfilename, ioImgInfo._kp);
00111     ioImgInfo._loadFail = (info.filename.size() == 0);
00112 
00113     // update ImgData
00114     if(ioImgInfo._needsremap)
00115     {
00116         ioImgInfo._detectWidth = max(info.width,info.height);
00117         ioImgInfo._detectHeight = max(info.width,info.height);
00118         ioImgInfo._projOpts.setWidth(ioImgInfo._detectWidth);
00119         ioImgInfo._projOpts.setHeight(ioImgInfo._detectHeight);
00120     }
00121     else
00122     {
00123         ioImgInfo._detectWidth = info.width;
00124         ioImgInfo._detectHeight = info.height;
00125     };
00126     ioImgInfo._descLength = info.dimensions;
00127 
00128     return true;
00129 }
00130 
00131 vigra::RGBValue<float> gray2RGB(float const& v)
00132 {
00133     return vigra::RGBValue<float>(v,v,v);
00134 }
00135 
00136 template <class SrcImageIterator, class SrcAccessor>
00137 void applyMaskAndCrop(vigra::triple<SrcImageIterator, SrcImageIterator, SrcAccessor> img, const HuginBase::SrcPanoImage& SrcImg)
00138 {
00139     vigra::Diff2D imgSize = img.second - img.first;
00140 
00141     // create dest y iterator
00142     SrcImageIterator yd(img.first);
00143     // loop over the image and transform
00144     for(int y=0; y < imgSize.y; ++y, ++yd.y)
00145     {
00146         // create x iterators
00147         SrcImageIterator xd(yd);
00148         for(int x=0; x < imgSize.x; ++x, ++xd.x)
00149         {
00150             if(!SrcImg.isInside(vigra::Point2D(x,y)))
00151             {
00152                 *xd=0;
00153             };
00154         }
00155     }
00156 }
00157 
00158 // save some intermediate images to disc if defined
00159 //#define DEBUG_LOADING_REMAPPING
00160 bool PanoDetector::AnalyzeImage(ImgData& ioImgInfo, const PanoDetector& iPanoDetector)
00161 {
00162     vigra::DImage final_img;
00163     vigra::BImage final_mask;
00164 
00165     try
00166     {
00167         ioImgInfo._loadFail=false;
00168 
00169         TRACE_IMG("Load image...");
00170         vigra::ImageImportInfo aImageInfo(ioImgInfo._name.c_str());
00171         vigra::FRGBImage RGBimg(aImageInfo.width(), aImageInfo.height());
00172         vigra::BImage mask;
00173 
00174         if (aImageInfo.numBands() == 1)
00175         {
00176             // grayscale image, convert to RGB. This is kind of stupid, but celeste wants RGB images...
00177             vigra::FImage tmpImg(aImageInfo.width(), aImageInfo.height());
00178             if (aImageInfo.numExtraBands() == 0)
00179             {
00180                 vigra::importImage(aImageInfo, destImage(tmpImg));
00181             }
00182             else if (aImageInfo.numExtraBands() == 1)
00183             {
00184                 mask.resize(aImageInfo.size());
00185                 importImageAlpha(aImageInfo, destImage(tmpImg), destImage(mask));
00186             }
00187             else
00188             {
00189                 TRACE_INFO("Image with multiple alpha channels are not supported");
00190                 ioImgInfo._loadFail = true;
00191                 return false;
00192             }
00193             //vigra::GrayToRGBAccessor<vigra::RGBValue<float> > ga;
00194             RGBimg.resize(aImageInfo.size());
00195             vigra::transformImage(srcImageRange(tmpImg), destImage(RGBimg), &gray2RGB);
00196         }
00197         else
00198         {
00199             if(aImageInfo.numExtraBands() == 1)
00200             {
00201                 mask.resize(aImageInfo.size());
00202                 importImageAlpha(aImageInfo, destImage(RGBimg), destImage(mask));
00203             }
00204             else
00205             {
00206                 if (aImageInfo.numExtraBands() == 0)
00207                 {
00208                     vigra::importImage(aImageInfo, destImage(RGBimg));
00209                 }
00210                 else
00211                 {
00212                     TRACE_INFO("Image with multiple alpha channels are not supported");
00213                     ioImgInfo._loadFail = true;
00214                     return false;
00215                 };
00216             };
00217         }
00218 
00219         //convert image, so that all (rgb) values are between 0 and 1
00220         if(aImageInfo.getPixelType() == std::string("FLOAT") || aImageInfo.getPixelType() == std::string("DOUBLE"))
00221         {
00222             vigra::RGBToGrayAccessor<vigra::RGBValue<float> > ga;
00223             vigra::FindMinMax<float> minmax;   // init functor
00224             vigra::inspectImage(srcImageRange(RGBimg, ga),minmax);
00225             double minVal = minmax.min;
00226             double maxVal = minmax.max;
00227             vigra_ext::applyMapping(srcImageRange(RGBimg), destImage(RGBimg), minVal, maxVal, 0);
00228         }
00229         else
00230         {
00231             vigra::transformImage(srcImageRange(RGBimg), destImage(RGBimg),
00232                                   vigra::functor::Arg1()/vigra::functor::Param(vigra_ext::getMaxValForPixelType(aImageInfo.getPixelType())));
00233         };
00234 
00235         if(ioImgInfo._needsremap)
00236         {
00237             TRACE_IMG("Remap image...");
00238             RemappedPanoImage<vigra::FRGBImage,vigra::BImage>* remapped=new RemappedPanoImage<vigra::FRGBImage,vigra::BImage>;
00239             vigra::FImage ffImg;
00240             MultiProgressDisplay* progress=new DummyMultiProgressDisplay();
00241             remapped->setPanoImage(iPanoDetector._panoramaInfoCopy.getImage(ioImgInfo._number),
00242                                    ioImgInfo._projOpts, ioImgInfo._projOpts.getROI());
00243             if(mask.width()>0)
00244             {
00245                 remapped->remapImage(vigra::srcImageRange(RGBimg),vigra::srcImage(mask),vigra_ext::INTERP_CUBIC,*progress);
00246             }
00247             else
00248             {
00249                 remapped->remapImage(vigra::srcImageRange(RGBimg),vigra_ext::INTERP_CUBIC,*progress);
00250             };
00251             RGBimg.resize(0,0);
00252             mask.resize(0,0);
00253             RGBimg=remapped->m_image;
00254             mask=remapped->m_mask;
00255             delete remapped;
00256             delete progress;
00257         }
00258         else
00259         {
00260             const SrcPanoImage& SrcImg=iPanoDetector._panoramaInfoCopy.getImage(ioImgInfo._number);
00261             if(SrcImg.hasActiveMasks() || (SrcImg.getCropMode()!=SrcPanoImage::NO_CROP && !SrcImg.getCropRect().isEmpty()))
00262             {
00263                 if(mask.width()!=aImageInfo.width() || mask.height()!=aImageInfo.height())
00264                 {
00265                     mask.resize(aImageInfo.size().width(),aImageInfo.size().height(),255);
00266                 };
00267                 //copy mask and crop from pto file into alpha layer
00268                 applyMaskAndCrop(vigra::destImageRange(mask), SrcImg);
00269             };
00270         };
00271 
00272         if(iPanoDetector.getCeleste())
00273         {
00274             vigra::FRGBImage scaled(ioImgInfo._detectWidth,ioImgInfo._detectHeight);
00275             if(iPanoDetector._downscale && (!ioImgInfo._needsremap))
00276             {
00277                 TRACE_IMG("Resize image...");
00278                 vigra::resizeImageNoInterpolation(srcImageRange(RGBimg),destImageRange(scaled));
00279                 if(mask.width()>0 && mask.height()>0)
00280                 {
00281                     final_mask.resize(ioImgInfo._detectWidth, ioImgInfo._detectHeight);
00282                     vigra::resizeImageNoInterpolation(srcImageRange(mask),destImageRange(final_mask));
00283                 };
00284             }
00285             else
00286             {
00287                 vigra::copyImage(srcImageRange(RGBimg),destImage(scaled));
00288                 if(mask.width()>0 && mask.height()>0)
00289                 {
00290                     final_mask.resize(ioImgInfo._detectWidth, ioImgInfo._detectHeight);
00291                     vigra::copyImage(srcImageRange(mask),destImage(final_mask));
00292                 };
00293             };
00294             RGBimg.resize(0,0);
00295             mask.resize(0,0);
00296             TRACE_IMG("Mask areas with clouds...");
00297             vigra::UInt16RGBImage image16(scaled.size());
00298             vigra::transformImage(srcImageRange(scaled),destImage(image16),vigra::functor::Arg1()*vigra::functor::Param(65535));
00299             int radius=iPanoDetector.getCelesteRadius();
00300             if(iPanoDetector._downscale)
00301             {
00302                 radius>>= 1;
00303             };
00304             if(radius<2)
00305             {
00306                 radius=2;
00307             };
00308             vigra::BImage celeste_mask=celeste::getCelesteMask(iPanoDetector.svmModel,image16,radius,iPanoDetector.getCelesteThreshold(),800,true,false);
00309 #ifdef DEBUG_LOADING_REMAPPING
00310             // DEBUG: export celeste mask
00311             std::ostringstream maskfilename;
00312             maskfilename << ioImgInfo._name << "_celeste_mask.JPG";
00313             vigra::ImageExportInfo maskexinfo(maskfilename.str().c_str());
00314             vigra::exportImage(srcImageRange(celeste_mask), maskexinfo);
00315 #endif
00316             image16.resize(0,0);
00317             if(final_mask.width()>0)
00318             {
00319                 vigra::copyImageIf(srcImageRange(celeste_mask),srcImage(final_mask),destImage(final_mask));
00320             }
00321             else
00322             {
00323                 final_mask.resize(ioImgInfo._detectWidth,ioImgInfo._detectHeight);
00324                 vigra::copyImage(srcImageRange(celeste_mask),destImage(final_mask));
00325             };
00326             celeste_mask.resize(0,0);
00327             TRACE_IMG("Convert to greyscale double...");
00328             final_img.resize(ioImgInfo._detectWidth,ioImgInfo._detectHeight);
00329             vigra::copyImage(scaled.upperLeft(), scaled.lowerRight(), vigra::RGBToGrayAccessor<vigra::RGBValue<double> >(),
00330                              final_img.upperLeft(), vigra::DImage::Accessor());
00331             scaled.resize(0,0);
00332         }
00333         else
00334         {
00335             //without celeste
00336             final_img.resize(ioImgInfo._detectWidth,ioImgInfo._detectHeight);
00337             if (iPanoDetector._downscale && !ioImgInfo._needsremap)
00338             {
00339                 // Downscale and convert to grayscale double format
00340                 TRACE_IMG("Resize to greyscale double...");
00341                 vigra::resizeImageNoInterpolation(
00342                     RGBimg.upperLeft(),
00343                     RGBimg.upperLeft() + vigra::Diff2D(aImageInfo.width(), aImageInfo.height()),
00344                     vigra::RGBToGrayAccessor<vigra::RGBValue<double> >(),
00345                     final_img.upperLeft(),
00346                     final_img.lowerRight(),
00347                     vigra::DImage::Accessor());
00348                 RGBimg.resize(0,0);
00349                 //downscale mask
00350                 if(mask.width()>0 && mask.height()>0)
00351                 {
00352                     final_mask.resize(ioImgInfo._detectWidth, ioImgInfo._detectHeight);
00353                     vigra::resizeImageNoInterpolation(srcImageRange(mask),destImageRange(final_mask));
00354                     mask.resize(0,0);
00355                 };
00356             }
00357             else
00358             {
00359                 // convert to grayscale
00360                 TRACE_IMG("Convert to greyscale double...");
00361                 vigra::copyImage(RGBimg.upperLeft(), RGBimg.lowerRight(), vigra::RGBToGrayAccessor<vigra::RGBValue<double> >(),
00362                                  final_img.upperLeft(), vigra::DImage::Accessor());
00363                 RGBimg.resize(0,0);
00364                 if(mask.width()>0 && mask.height()>0)
00365                 {
00366                     final_mask.resize(ioImgInfo._detectWidth, ioImgInfo._detectHeight);
00367                     vigra::copyImage(srcImageRange(mask),destImage(final_mask));
00368                     mask.resize(0,0);
00369                 };
00370             };
00371         };
00372 
00373         //now scale image from 0..1 to 0..255
00374         vigra::transformImage(srcImageRange(final_img),destImage(final_img),vigra::functor::Arg1()*vigra::functor::Param(255));
00375 
00376 #ifdef DEBUG_LOADING_REMAPPING
00377         // DEBUG: export remapped
00378         std::ostringstream filename;
00379         filename << ioImgInfo._name << "_grey.JPG";
00380         vigra::ImageExportInfo exinfo(filename.str().c_str());
00381         vigra::exportImage(srcImageRange(final_img), exinfo);
00382 #endif
00383 
00384         // Build integral image
00385         TRACE_IMG("Build integral image...");
00386         ioImgInfo._ii.init(final_img);
00387         final_img.resize(0,0);
00388 
00389         // compute distance map
00390         if(final_mask.width()>0 && final_mask.height()>0)
00391         {
00392             TRACE_IMG("Build distance map...");
00393             //apply threshold, in case loaded mask contains other values than 0 and 255
00394             vigra::transformImage(srcImageRange(final_mask), destImage(final_mask),
00395                                   vigra::Threshold<vigra::BImage::PixelType, vigra::BImage::PixelType>(1, 255, 0, 255));
00396             ioImgInfo._distancemap.resize(final_mask.width(),final_mask.height(),0);
00397             vigra::distanceTransform(srcImageRange(final_mask), destImage(ioImgInfo._distancemap), 255, 2);
00398 #ifdef DEBUG_LOADING_REMAPPING
00399             std::ostringstream maskfilename;
00400             maskfilename << ioImgInfo._name << "_mask.JPG";
00401             vigra::ImageExportInfo maskexinfo(maskfilename.str().c_str());
00402             vigra::exportImage(srcImageRange(final_mask), maskexinfo);
00403             std::ostringstream distfilename;
00404             distfilename << ioImgInfo._name << "_distancemap.JPG";
00405             vigra::ImageExportInfo distexinfo(distfilename.str().c_str());
00406             vigra::exportImage(srcImageRange(ioImgInfo._distancemap), distexinfo);
00407 #endif
00408             final_mask.resize(0,0);
00409         };
00410     }
00411     catch (std::exception& e)
00412     {
00413         TRACE_INFO("An error happened while loading image : caught exception: " << e.what() << endl);
00414         ioImgInfo._loadFail=true;
00415         return false;
00416     }
00417 
00418     return true;
00419 }
00420 
00421 
00422 bool PanoDetector::FindKeyPointsInImage(ImgData& ioImgInfo, const PanoDetector& iPanoDetector)
00423 {
00424     TRACE_IMG("Find keypoints...");
00425 
00426     // setup the detector
00427     KeyPointDetector aKP;
00428 
00429     // detect the keypoints
00430     KeyPointVectInsertor aInsertor(ioImgInfo._kp);
00431     aKP.detectKeypoints(ioImgInfo._ii, aInsertor);
00432 
00433     TRACE_IMG("Found "<< ioImgInfo._kp.size() << " interest points.");
00434 
00435     return true;
00436 }
00437 
00438 bool PanoDetector::FilterKeyPointsInImage(ImgData& ioImgInfo, const PanoDetector& iPanoDetector)
00439 {
00440     TRACE_IMG("Filtering keypoints...");
00441 
00442     lfeat::Sieve<lfeat::KeyPointPtr, lfeat::KeyPointPtrSort > aSieve(iPanoDetector.getSieve1Width(),
00443             iPanoDetector.getSieve1Height(),
00444             iPanoDetector.getSieve1Size());
00445 
00446     // insert the points in the Sieve if they are not masked
00447     double aXF = (double)iPanoDetector.getSieve1Width() / (double)ioImgInfo._detectWidth;
00448     double aYF = (double)iPanoDetector.getSieve1Height() / (double)ioImgInfo._detectHeight;
00449 
00450     bool distmap_valid=(ioImgInfo._distancemap.width()>0 && ioImgInfo._distancemap.height()>0);
00451     BOOST_FOREACH(KeyPointPtr& aK, ioImgInfo._kp)
00452     {
00453         if(distmap_valid)
00454         {
00455             if(aK->_x > 0 && aK->_x < ioImgInfo._distancemap.width() && aK->_y > 0 && aK->_y < ioImgInfo._distancemap.height()
00456                     && ioImgInfo._distancemap((int)(aK->_x),(int)(aK->_y)) >aK->_scale*8)
00457             {
00458                 //cout << " dist from border:" << ioImgInfo._distancemap((int)(aK->_x),(int)(aK->_y)) << " required dist: " << aK->_scale*12 << std::endl;
00459                 aSieve.insert(aK, (int)(aK->_x * aXF), (int)(aK->_y * aYF));
00460             }
00461         }
00462         else
00463         {
00464             aSieve.insert(aK, (int)(aK->_x * aXF), (int)(aK->_y * aYF));
00465         };
00466     }
00467 
00468     // pull remaining values from the sieve
00469     ioImgInfo._kp.clear();
00470 
00471     // make an extractor and pull the points
00472     SieveExtractorKP aSieveExt(ioImgInfo._kp);
00473     aSieve.extract(aSieveExt);
00474 
00475     TRACE_IMG("Kept " << ioImgInfo._kp.size() << " interest points.");
00476 
00477     return true;
00478 
00479 }
00480 
00481 bool PanoDetector::MakeKeyPointDescriptorsInImage(ImgData& ioImgInfo, const PanoDetector& iPanoDetector)
00482 {
00483     TRACE_IMG("Make keypoint descriptors...");
00484 
00485     // build a keypoint descriptor
00486     CircularKeyPointDescriptor aKPD(ioImgInfo._ii);
00487 
00488     // vector for keypoints with more than one orientation
00489     KeyPointVect_t kp_new_ori;
00490     BOOST_FOREACH(KeyPointPtr& aK, ioImgInfo._kp)
00491     {
00492         double angles[4];
00493         int nAngles = aKPD.assignOrientation(*aK, angles);
00494         for (int i=0; i < nAngles; i++)
00495         {
00496             // duplicate Keypoint with additional angles
00497             KeyPointPtr aKn = KeyPointPtr ( new lfeat::KeyPoint ( *aK ) );
00498             aKn->_ori = angles[i];
00499             kp_new_ori.push_back(aKn);
00500         }
00501     }
00502     ioImgInfo._kp.insert(ioImgInfo._kp.end(), kp_new_ori.begin(), kp_new_ori.end());
00503 
00504     BOOST_FOREACH(KeyPointPtr& aK, ioImgInfo._kp)
00505     {
00506         aKPD.makeDescriptor(*aK);
00507     }
00508     // store the descriptor length
00509     ioImgInfo._descLength = aKPD.getDescriptorLength();
00510     return true;
00511 }
00512 
00513 bool PanoDetector::RemapBackKeypoints(ImgData& ioImgInfo, const PanoDetector& iPanoDetector)
00514 {
00515 
00516     double scale=iPanoDetector._downscale ? 2.0:1.0;
00517 
00518     if (!ioImgInfo._needsremap)
00519     {
00520         if(scale != 1.0)
00521         {
00522             BOOST_FOREACH(KeyPointPtr& aK, ioImgInfo._kp)
00523             {
00524                 aK->_x *= scale;
00525                 aK->_y *= scale;
00526                 aK->_scale *= scale;
00527             }
00528         };
00529     }
00530     else
00531     {
00532         TRACE_IMG("Remapping back keypoints...");
00533         HuginBase::PTools::Transform trafo1;
00534         trafo1.createTransform(iPanoDetector._panoramaInfoCopy.getSrcImage(ioImgInfo._number),
00535                                ioImgInfo._projOpts);
00536 
00537         int dx1 = ioImgInfo._projOpts.getROI().left();
00538         int dy1 = ioImgInfo._projOpts.getROI().top();
00539 
00540         BOOST_FOREACH(KeyPointPtr& aK, ioImgInfo._kp)
00541         {
00542             double xout, yout;
00543             if(trafo1.transformImgCoord(xout, yout, aK->_x + dx1, aK->_y+ dy1))
00544             {
00545                 // downscaling is take care of by the remapping transform
00546                 // no need for multiplying the scale factor...
00547                 aK->_x=xout;
00548                 aK->_y=yout;
00549                 aK->_scale *= scale;
00550             }
00551         }
00552     }
00553     return true;
00554 }
00555 
00556 bool PanoDetector::BuildKDTreesInImage(ImgData& ioImgInfo, const PanoDetector& iPanoDetector)
00557 {
00558     TRACE_IMG("Build KDTree...");
00559 
00560     if(ioImgInfo._kp.size()==0)
00561     {
00562         return false;
00563     };
00564     // build a vector of KDElemKeyPointPtr
00565 
00566     // create feature vector matrix for flann
00567     ioImgInfo._flann_descriptors = flann::Matrix<double>(new double[ioImgInfo._kp.size()*ioImgInfo._descLength],
00568                                    ioImgInfo._kp.size(), ioImgInfo._descLength);
00569     int i = 0;
00570     BOOST_FOREACH(KeyPointPtr& aK, ioImgInfo._kp)
00571     {
00572         memcpy(ioImgInfo._flann_descriptors[i++], aK->_vec, sizeof(double)*ioImgInfo._descLength);
00573     }
00574 
00575     // build query structure
00576     ioImgInfo._flann_index = new flann::Index<flann::L2<double> > (ioImgInfo._flann_descriptors, flann::KDTreeIndexParams(4));
00577     ioImgInfo._flann_index->buildIndex();
00578 
00579     return true;
00580 }
00581 
00582 bool PanoDetector::FreeMemoryInImage(ImgData& ioImgInfo, const PanoDetector& iPanoDetector)
00583 {
00584     TRACE_IMG("Freeing memory...");
00585 
00586     ioImgInfo._ii.clean();
00587     ioImgInfo._distancemap.resize(0,0);
00588 
00589     return true;
00590 }
00591 
00592 
00593 bool PanoDetector::FindMatchesInPair(MatchData& ioMatchData, const PanoDetector& iPanoDetector)
00594 {
00595     TRACE_PAIR("Find Matches...");
00596 
00597     // retrieve the KDTree of image 2
00598     flann::Index<flann::L2<double> > * index2 = ioMatchData._i2->_flann_index;
00599 
00600     // retrieve query points from image 1
00601     flann::Matrix<double> & query = ioMatchData._i1->_flann_descriptors;
00602 
00603     // storage for sorted 2 best matches
00604     int nn = 2;
00605     flann::Matrix<int> indices(new int[query.rows*nn], query.rows, nn);
00606     flann::Matrix<double> dists(new double[query.rows*nn], query.rows, nn);
00607 
00608     // perform matching using flann
00609     index2->knnSearch(query, indices, dists, nn, flann::SearchParams(iPanoDetector.getKDTreeSearchSteps()));
00610 
00611     //typedef KDTreeSpace::BestMatch<KDElemKeyPoint>            BM_t;
00612     //std::set<BM_t, std::greater<BM_t> >       aBestMatches;
00613 
00614     // store the matches already found to avoid 2 points in image1
00615     // match the same point in image2
00616     // both matches will be removed.
00617     set<int> aAlreadyMatched;
00618     set<int> aBadMatch;
00619 
00620     // unfiltered vector of matches
00621     typedef std::pair<KeyPointPtr, int> TmpPair_t;
00622     std::vector<TmpPair_t>      aUnfilteredMatches;
00623 
00624     //PointMatchVector_t aMatches;
00625 
00626     // go through all the keypoints of image 1
00627     for (unsigned aKIt = 0; aKIt < query.rows; ++aKIt)
00628     {
00629         // accept the match if the second match is far enough
00630         // put a lower value for stronger matching default 0.15
00631         if (dists[aKIt][0] > iPanoDetector.getKDTreeSecondDistance()  * dists[aKIt][1])
00632         {
00633             continue;
00634         }
00635 
00636         // check if the kdtree match number is already in the already matched set
00637         if (aAlreadyMatched.find(indices[aKIt][0]) != aAlreadyMatched.end())
00638         {
00639             // add to delete list and continue
00640             aBadMatch.insert(indices[aKIt][0]);
00641             continue;
00642         }
00643 
00644         // TODO: add check for duplicate matches (can happen if a keypoint gets multiple orientations)
00645 
00646         // add the match number in already matched set
00647         aAlreadyMatched.insert(indices[aKIt][0]);
00648 
00649         // add the match to the unfiltered list
00650         aUnfilteredMatches.push_back(TmpPair_t(ioMatchData._i1->_kp[aKIt], indices[aKIt][0]));
00651     }
00652 
00653     // now filter and fill the vector of matches
00654     BOOST_FOREACH(TmpPair_t& aP, aUnfilteredMatches)
00655     {
00656         // if the image2 match number is in the badmatch set, skip it.
00657         if (aBadMatch.find(aP.second) != aBadMatch.end())
00658         {
00659             continue;
00660         }
00661 
00662         // add the match in the output vector
00663         ioMatchData._matches.push_back(lfeat::PointMatchPtr( new lfeat::PointMatch(aP.first, ioMatchData._i2->_kp[aP.second])));
00664     }
00665 
00666     TRACE_PAIR("Found " << ioMatchData._matches.size() << " matches.");
00667     return true;
00668 }
00669 
00670 bool PanoDetector::RansacMatchesInPair(MatchData& ioMatchData, const PanoDetector& iPanoDetector)
00671 {
00672     // Use panotools model for wide angle lenses
00673     RANSACOptimizer::Mode rmode = iPanoDetector._ransacMode;
00674     if (rmode == RANSACOptimizer::HOMOGRAPHY ||
00675             (rmode == RANSACOptimizer::AUTO && iPanoDetector._panoramaInfo->getImage(ioMatchData._i1->_number).getHFOV() < 65 &&
00676              iPanoDetector._panoramaInfo->getImage(ioMatchData._i2->_number).getHFOV() < 65))
00677     {
00678         return RansacMatchesInPairHomography(ioMatchData, iPanoDetector);
00679     }
00680     else
00681     {
00682         return RansacMatchesInPairCam(ioMatchData, iPanoDetector);
00683     }
00684 }
00685 
00686 // new code with fisheye aware ransac
00687 bool PanoDetector::RansacMatchesInPairCam(MatchData& ioMatchData, const PanoDetector& iPanoDetector)
00688 {
00689     TRACE_PAIR("RANSAC Filtering with Panorama model...");
00690 
00691     if (ioMatchData._matches.size() < (unsigned int)iPanoDetector.getMinimumMatches())
00692     {
00693         TRACE_PAIR("Too few matches ... removing all of them.");
00694         ioMatchData._matches.clear();
00695         return true;
00696     }
00697 
00698     if (ioMatchData._matches.size() < 6)
00699     {
00700         TRACE_PAIR("Not enough matches for RANSAC filtering.");
00701         return true;
00702     }
00703 
00704     // setup a panorama project with the two images.
00705     // is this threadsafe (is this read only access?)
00706     UIntSet imgs;
00707     int pano_i1 = ioMatchData._i1->_number;
00708     int pano_i2 = ioMatchData._i2->_number;
00709     imgs.insert(pano_i1);
00710     imgs.insert(pano_i2);
00711     int pano_local_i1 = 0;
00712     int pano_local_i2 = 1;
00713     if (pano_i1 > pano_i2)
00714     {
00715         pano_local_i1 = 1;
00716         pano_local_i2 = 0;
00717     }
00718 
00719     // perform ransac matching.
00720     // ARGH the panotools optimizer uses global variables is not reentrant
00721     std::vector<int> inliers;
00722     {
00723         ZThread::Guard<ZThread::FastMutex> g(aPanoToolsMutex);
00724 
00725         PanoramaData* panoSubset = iPanoDetector._panoramaInfo->getNewSubset(imgs);
00726 
00727         // create control point vector
00728         CPVector controlPoints(ioMatchData._matches.size());
00729         int i=0;
00730         BOOST_FOREACH(PointMatchPtr& aM, ioMatchData._matches)
00731         {
00732             controlPoints[i] = ControlPoint(pano_local_i1, aM->_img1_x, aM->_img1_y,
00733                                             pano_local_i2, aM->_img2_x, aM->_img2_y);
00734             i++;
00735         }
00736         panoSubset->setCtrlPoints(controlPoints);
00737 
00738 
00739         PT_setProgressFcn(ptProgress);
00740         PT_setInfoDlgFcn(ptinfoDlg);
00741 
00742         RANSACOptimizer::Mode rmode = iPanoDetector._ransacMode;
00743         if (rmode == RANSACOptimizer::AUTO)
00744         {
00745             rmode = RANSACOptimizer::RPY;
00746         }
00747         inliers = HuginBase::RANSACOptimizer::findInliers(*panoSubset, pano_local_i1, pano_local_i2,
00748                   iPanoDetector.getRansacDistanceThreshold(), rmode);
00749         PT_setProgressFcn(NULL);
00750         PT_setInfoDlgFcn(NULL);
00751         delete panoSubset;
00752 
00753         TRACE_PAIR("Removed " << controlPoints.size() - inliers.size() << " matches. " << inliers.size() << " remaining.");
00754         if (inliers.size() < 0.5 * controlPoints.size())
00755         {
00756             // more than 50% of matches were removed, ignore complete pair...
00757             TRACE_PAIR("RANSAC found more than 50% outliers, removing all matches");
00758             ioMatchData._matches.clear();
00759             return true;
00760         }
00761     }
00762 
00763 
00764     if (inliers.size() < (unsigned int)iPanoDetector.getMinimumMatches())
00765     {
00766         TRACE_PAIR("Too few matches ... removing all of them.");
00767         ioMatchData._matches.clear();
00768         return true;
00769     }
00770 
00771     // keep only inlier matches
00772     PointMatchVector_t aInlierMatches;
00773     aInlierMatches.reserve(inliers.size());
00774 
00775     BOOST_FOREACH(int idx, inliers)
00776     {
00777         aInlierMatches.push_back(ioMatchData._matches[idx]);
00778     }
00779     ioMatchData._matches = aInlierMatches;
00780 
00781     /*
00782     if (iPanoDetector.getTest())
00783         TestCode::drawRansacMatches(ioMatchData._i1->_name, ioMatchData._i2->_name, ioMatchData._matches,
00784                                                                 aRemovedMatches, aRansacFilter, iPanoDetector.getDownscale());
00785     */
00786 
00787     return true;
00788 }
00789 
00790 // homography based ransac matching
00791 bool PanoDetector::RansacMatchesInPairHomography(MatchData& ioMatchData, const PanoDetector& iPanoDetector)
00792 {
00793     TRACE_PAIR("RANSAC Filtering...");
00794 
00795     if (ioMatchData._matches.size() < (unsigned int)iPanoDetector.getMinimumMatches())
00796     {
00797         TRACE_PAIR("Too few matches ... removing all of them.");
00798         ioMatchData._matches.clear();
00799         return true;
00800     }
00801 
00802     if (ioMatchData._matches.size() < 6)
00803     {
00804         TRACE_PAIR("Not enough matches for RANSAC filtering.");
00805         return true;
00806     }
00807 
00808     PointMatchVector_t aRemovedMatches;
00809 
00810     Ransac aRansacFilter;
00811     aRansacFilter.setIterations(iPanoDetector.getRansacIterations());
00812     int thresholdDistance=iPanoDetector.getRansacDistanceThreshold();
00813     //increase RANSAC distance if the image were remapped to not exclude
00814     //too much points in this case
00815     if(ioMatchData._i1->_needsremap || ioMatchData._i2->_needsremap)
00816     {
00817         thresholdDistance*=5;
00818     }
00819     aRansacFilter.setDistanceThreshold(thresholdDistance);
00820     aRansacFilter.filter(ioMatchData._matches, aRemovedMatches);
00821 
00822 
00823     TRACE_PAIR("Removed " << aRemovedMatches.size() << " matches. " << ioMatchData._matches.size() << " remaining.");
00824 
00825     if (aRemovedMatches.size() > ioMatchData._matches.size())
00826     {
00827         // more than 50% of matches were removed, ignore complete pair...
00828         TRACE_PAIR("More than 50% outliers, removing all matches");
00829         ioMatchData._matches.clear();
00830         return true;
00831     }
00832 
00833     if (iPanoDetector.getTest())
00834         TestCode::drawRansacMatches(ioMatchData._i1->_name, ioMatchData._i2->_name, ioMatchData._matches,
00835                                     aRemovedMatches, aRansacFilter, iPanoDetector.getDownscale());
00836 
00837     return true;
00838 
00839 }
00840 
00841 
00842 bool PanoDetector::FilterMatchesInPair(MatchData& ioMatchData, const PanoDetector& iPanoDetector)
00843 {
00844     TRACE_PAIR("Clustering matches...");
00845 
00846     if (ioMatchData._matches.size() < 2)
00847     {
00848         return true;
00849     }
00850 
00851     // compute min,max of x,y for image1
00852 
00853     double aMinX = numeric_limits<double>::max();
00854     double aMinY = numeric_limits<double>::max();
00855     double aMaxX = -numeric_limits<double>::max();
00856     double aMaxY = -numeric_limits<double>::max();
00857 
00858     BOOST_FOREACH(PointMatchPtr& aM, ioMatchData._matches)
00859     {
00860         if (aM->_img1_x < aMinX)
00861         {
00862             aMinX = aM->_img1_x;
00863         }
00864         if (aM->_img1_x > aMaxX)
00865         {
00866             aMaxX = aM->_img1_x;
00867         }
00868 
00869         if (aM->_img1_y < aMinY)
00870         {
00871             aMinY = aM->_img1_y;
00872         }
00873         if (aM->_img1_y > aMaxY)
00874         {
00875             aMaxY = aM->_img1_y;
00876         }
00877     }
00878 
00879     double aSizeX = aMaxX - aMinX + 2; // add 2 so max/aSize is strict < 1
00880     double aSizeY = aMaxY - aMinY + 2;
00881 
00882     //
00883 
00884     Sieve<PointMatchPtr, PointMatchPtrSort> aSieve(iPanoDetector.getSieve2Width(),
00885             iPanoDetector.getSieve2Height(),
00886             iPanoDetector.getSieve2Size());
00887 
00888     // insert the points in the Sieve
00889     double aXF = (double)iPanoDetector.getSieve2Width() / aSizeX;
00890     double aYF = (double)iPanoDetector.getSieve2Height() / aSizeY;
00891     int aCount = 0;
00892     BOOST_FOREACH(PointMatchPtr& aM, ioMatchData._matches)
00893     {
00894         aSieve.insert(aM, (int)((aM->_img1_x - aMinX) * aXF), (int)((aM->_img1_y - aMinY) * aYF));
00895         aCount++;
00896     }
00897 
00898     // pull remaining values from the sieve
00899     ioMatchData._matches.clear();
00900 
00901     // make an extractor and pull the points
00902     SieveExtractorMatch aSieveExt(ioMatchData._matches);
00903     aSieve.extract(aSieveExt);
00904 
00905     TRACE_PAIR("Kept " << ioMatchData._matches.size() << " matches.");
00906     return true;
00907 }
00908 
00909 void PanoDetector::writeOutput()
00910 {
00911     // Write output pto file
00912 
00913     ofstream aOut(_outputFile.c_str(), ios_base::trunc);
00914     if( !aOut )
00915     {
00916         cerr << "ERROR : "
00917              << "Couldn't open file '" << _outputFile << "'!" << endl; //STS
00918         return;
00919     }
00920 
00921     aOut << "# pto project file generated by Hugin's cpfind" << endl << endl;
00922 
00923     _panoramaInfo->removeDuplicateCtrlPoints();
00924     AppBase::DocumentData::ReadWriteError err = _panoramaInfo->writeData(aOut);
00925     if (err != AppBase::DocumentData::SUCCESSFUL)
00926     {
00927         cerr << "ERROR couldn't write to output file '" << _outputFile << "'!" << endl;
00928         return;
00929     }
00930 }
00931 
00932 void PanoDetector::writeKeyfile(ImgData& imgInfo)
00933 {
00934     // Write output keyfile
00935 
00936     ofstream aOut(imgInfo._keyfilename.c_str(), ios_base::trunc);
00937 
00938     SIFTFormatWriter writer(aOut);
00939 
00940     int origImgWidth =  _panoramaInfo->getImage(imgInfo._number).getSize().width();
00941     int origImgHeight =  _panoramaInfo->getImage(imgInfo._number).getSize().height();
00942 
00943     ImageInfo img_info(imgInfo._name, origImgWidth, origImgHeight);
00944 
00945     writer.writeHeader ( img_info, imgInfo._kp.size(), imgInfo._descLength );
00946 
00947     BOOST_FOREACH ( KeyPointPtr& aK, imgInfo._kp )
00948     {
00949         writer.writeKeypoint ( aK->_x, aK->_y, aK->_scale, aK->_ori, aK->_score,
00950                                imgInfo._descLength, aK->_vec );
00951     }
00952     writer.writeFooter();
00953 }
00954 

Generated on 24 Nov 2014 for Hugintrunk by  doxygen 1.4.7