source: branches/NewStructure/src/Finder/Searcher.cc @ 1441

#include <iostream>
#include <duchamp/duchamp.hh>
#include <duchamp/Finder/Searcher.hh>
#include <duchamp/Finder/DuchampFinder.hh>
#include <duchamp/DataArrays/AstroImage.hh>
#include <duchamp/Detection/detection.hh>
#include <duchamp/Utils/feedback.hh>

namespace duchamp
{

    Searcher::Searcher()
    {
    }

    Searcher::Searcher(DuchampFinder& finder)
    {
        this->itsParam = finder.par();
        this->itsStats = finder.stats();

        if(finder.par()->getFlagATrous())
            this->itsImage = finder.recon();
        else if(finder.par()->getFlagSmooth())
            this->itsImage = finder.recon();
        else this->itsImage = finder.image();
    }

    Searcher::Searcher(const Searcher& other)
    {
        this->operator=(other);
    }

    Searcher& Searcher::operator= (const Searcher& other)
    {
        if ( this == &other) return *this;
        this->itsImage = other.itsImage;
        this->itsDetectionList = other.itsDetectionList;
        this->itsParam = other.itsParam;
        this->itsStats = other.itsStats;
        return *this;
    }

    void Searcher::search()
    {
        if(this->itsParam->getSearchType() == "spatial")  
            this->spatialSearch();
        else
            this->spectralSearch();
    }
    
    void Searcher::spatialSearch()
    {
    }

    void Searcher::spectralSearch()
    {
        std::vector<size_t> dim = this->itsImage->dim();
        size_t zdim = dim[2];
        size_t xySize = dim[0] * dim[1];
        int num = 0;

        if(zdim>1){
    
            ProgressBar bar;
            if(this->itsParam->isVerbose()) bar.init(xySize);

            bool *doPixel = new bool[xySize];
            for(size_t npix=0; npix<xySize; npix++){
                doPixel[npix] = false;
                for(size_t z=0;z<zdim || doPixel[npix];z++){
                    doPixel[npix] = doPixel[npix] || 
                        (!this->itsParam->isBlank(Array[npix+xySize*z]) && !this->itsParam->isFlaggedChannel(z));
                }
                // doPixel[i] is false only when there are no good pixels in spectrum
                //  of pixel #i.
            }

            size_t *specdim = new size_t[2];
            specdim[0] = zdim; specdim[1]=1;
            Image *spectrum = new Image(specdim);
            delete [] specdim;
            spectrum->saveParam(par);
            spectrum->saveStats(stats);
            //    spectrum->setMinSize(this->itsParam->getMinChannels());
            spectrum->setMinSize(1);

            for(size_t y=0; y<dim[1]; y++){
                for(size_t x=0; x<dim[0]; x++){

                    size_t npix = y*dim[0] + x;
                    if( this->itsParam->isVerbose() ) bar.update(npix+1);
        
                    if(doPixel[npix]){
                        spectrum->extractSpectrum(Array,dim,npix);
                        spectrum->removeFlaggedChannels();
                        std::vector<Scan> objlist = spectrum->findSources1D();
                        std::vector<Scan>::iterator obj;
                        num += objlist.size();
                        for(obj=objlist.begin();obj<objlist.end();obj++){
                            Detection newObject;
                            // Fix up coordinates of each pixel to match original array
                            for(int z=obj->getX();z<=obj->getXmax();z++) {
                                newObject.addPixel(x,y,z);
                            }
                            newObject.setOffsets(par);
                            if(this->itsParam->getFlagTwoStageMerging()) mergeIntoList(newObject,outputList,par);
                            else outputList.push_back(newObject);
                        }
                    }
                }
            }

            delete spectrum;
            delete [] doPixel;
  

            if(this->itsParam->isVerbose()){
                bar.remove();
                std::cout << "Found " << num << ".\n";
            }

        }

        return outputList;

    }

}