aboutsummaryrefslogtreecommitdiff
path: root/fly-tools/FlyTrackingFilter.cpp
diff options
context:
space:
mode:
Diffstat (limited to 'fly-tools/FlyTrackingFilter.cpp')
-rw-r--r--fly-tools/FlyTrackingFilter.cpp666
1 files changed, 0 insertions, 666 deletions
diff --git a/fly-tools/FlyTrackingFilter.cpp b/fly-tools/FlyTrackingFilter.cpp
deleted file mode 100644
index 9020522..0000000
--- a/fly-tools/FlyTrackingFilter.cpp
+++ /dev/null
@@ -1,666 +0,0 @@
-#include <iostream>
-#include <string>
-#include <sstream>
-#include <cmath>
-#include <vector>
-#include <list>
-#include <fstream>
-#include <cassert>
-#include <cstdlib>
-
-#include <ImageMagick/Magick++.h>
-
-#include <gsl/gsl_matrix.h>
-#include <gsl/gsl_vector.h>
-#include <gsl/gsl_blas.h>
-#include <gsl/gsl_eigen.h>
-
-#include "FrameInfo.h"
-
-using namespace Magick;
-using namespace std;
-void findObj(Image* img, int x, int y, vector<pair<int,int> > & shape ,bool eightCon=true, bool colorLookingFor=true);
-void eightConnObj(Image* img, int x, int y, vector<pair<int, int> > & obj, bool color=true);
-void fourConnObj(Image* img, int x, int y, vector<pair<int, int> > & obj, bool color=true);
-vector<double> covariantDecomposition(vector<pair<int,int> > & points);
-pair<int,int> getCentroid(vector<pair<int,int> > & points);
-bool isInterface(Image* orig, unsigned int x, unsigned int y);
-void writeFrameImage(int fn, string imS);
-int roundT(double v) {return int(v+0.5);}
-
-
-const double PI = atan(1.0)*4.0;
-const double FACTOR_EIGEN = 100;
-
-Image* residual;
-
-ostream &operator<<(ostream &out, FlyObject & fO) {
- fO.output(out);
- return out;
-}
-
-ostream &operator<<(ostream &out, FrameInfo & fI) {
- fI.output(out);
- return out;
-}
-
-vector<vector<pair<int, int> > > shapeVectors;
-vector<pair<int,int> > shape;
-vector<pair<int, int> > sizeNIndexVector;
-
-void bubbleSort() {
-
- for(int i=1; i<sizeNIndexVector.size(); i++) {
- for(int j=0; j<sizeNIndexVector.size()-i; j++) {
- pair<int, int> a = sizeNIndexVector[j];
- pair<int, int> b = sizeNIndexVector[j+1];
-
- if (a.first < b.first) {
- pair<int, int> c = sizeNIndexVector[j];
- sizeNIndexVector[j] = sizeNIndexVector[j+1];
- sizeNIndexVector[j+1] = c;
- }
- }
- }
-
-}
-
-void fillResidualWithObj(vector<pair<int, int> > & obj, ColorRGB c)
-{
- for (unsigned int i = 0; i<obj.size(); i++)
- residual->pixelColor(obj[i].first, obj[i].second, c);
-}
-
-void writeHist(const char* filename, map<unsigned int, unsigned int> & len)
-{
- map<unsigned int,unsigned int>::iterator front = len.begin(),
- back = len.end();
- back--;
-
-
- unsigned int first = front->first, last = back->first;
- /*if (cutoff != -1 && cutoff < int(last))
- last = cutoff;
- */
- cout << "Min: " << first << endl
- << "Max: " << last << endl
- << "Count: " << last-first << endl;
- //vector<unsigned int> hist(last-first, 0);
- vector<unsigned int> hist(last+1, 0);
-
- cout << "hist size: " << hist.size() << endl;
- try{
- for(unsigned int j = 0; j<first; j++) {
- hist[j] = 0;
- }
- for (unsigned int j = first; j<=last; j++)
- {
- /*if ( roundT(j-first) >= int(hist.size()) )
- hist.resize(j-first,0);
- hist[roundT(j-first)] = len[j];
- */
-
- /*if ( roundT(j) >= int(hist.size()) )
- hist.resize(j,0);
- hist[roundT(j)] = len[j];
- */
- hist[j] = len[j];
- }
- }
- catch (...)
- { cerr << "Bad histogram bucketing" << endl; }
-
- /*if ( (cutoff >= 0) && (cutoff<int(hist.size())) )
- hist.resize(cutoff);
- */
- len.clear();
- try
- {
- ofstream fout(filename);
- for (unsigned int i = 0; i<hist.size(); i++) {
- fout << hist[i] << endl;
- }
- fout << first << " " << last << " " << hist.size() << endl;
- fout.close();
- }
- catch (...)
- { cerr << "Bad memory loc for opening file" << endl; }
-}
-
-int main(int argc, char* argv[])
-{
- if (argc < 5)
- {
- cerr << "Usage: executablename <inputFile.txt> <ratio_largest_to_second_largest> <InputLocationOfMaskImage> <outputFolderName>" << endl; // input file contains name of the
- // input image files
- return -1;
- }
-
- MagickCore::SetMagickResourceLimit(MagickCore::MemoryResource, 1536);
- MagickCore::SetMagickResourceLimit(MagickCore::MapResource, 2048);
-
- string fileName;
- ifstream inputFile(argv[1]);
- if (inputFile.fail() ) {
- cout << "cannot open the input file that contains name of the input images\n";
- exit(1);
- }
-
- // get the output file name along with the location from argv[4]
- string outputFileLocation(argv[4]);
- string outputFileName = outputFileLocation + "final/outputFile.txt";
- ofstream outputFile(outputFileName.c_str());
-
- // get the location of the input mask file
- string inputMaskFileLocation(argv[3]);
-
- int frameCounter = 0;
- // use 15, 20 or 10
- // 15 found to be working correct
- double ratioSecondLargestToLargest = atof(argv[2]);
-
- cout << "Ratio is 1/"<<ratioSecondLargestToLargest<< " = "<<(1/ratioSecondLargestToLargest)<<endl;
-
- outputFile<<"Ratio given 1/"<<ratioSecondLargestToLargest<<" = "<<(1/ratioSecondLargestToLargest)<<endl;
-
- ratioSecondLargestToLargest = 1/ratioSecondLargestToLargest;
- // to find the largest object avg size
-
- long totalSize = 0;
-
- char buffer[100];
-
- while (inputFile>>fileName) {
-
- string savedFileName = fileName;
-
- // get the input mask file
-// fileName = "input/"+fileName;
- fileName = inputMaskFileLocation + fileName;
-
- Image* img = new Image(fileName.c_str());
- int width = img->columns(),height = img->rows();
- Image* imgForFilter;
- sprintf(buffer,"%ix%i",width,height);
-
- // residual image is initialized with black representing not visited.
- residual = new Image(buffer, "black");
-
- imgForFilter = new Image(buffer, "white");
- cout << "Reading "<<savedFileName<<endl;
- cout << "Filter wxh "<<width<<","<<height<<endl;
- shape.clear();
- // find the black background from location (0,0)
- findObj(img, 0, 0, shape, false, false);
- int s = shape.size();
- if (s > 0)
- cout << "black object size is "<<s;
- for (int i=0; i<s; i++) {
- imgForFilter->pixelColor(shape[i].first, shape[i].second, "black");
- }
-
- // store the intermediate file in temp folder under the output location
-// string oFilteredFileName = "output/filtered/temp/"+outputFile;
- string oFilteredFileName = outputFileLocation +"temp/"+savedFileName;
-
- cout << "Saving the filtered image "<< oFilteredFileName<<endl;
- imgForFilter->write(oFilteredFileName.c_str());
-
- delete residual;
- /*
- residual = new Image(buffer, "black");
-
- shapeVectors.clear();
- sizeNIndexVector.clear();
-
- // find the two largest object
- int objectCounter = 0;
- for (int x=0; x<width; x++) {
- for (int y=0; y<height; y++) {
- //find the white object now using eight connected neighbours
- shape.clear();
- findObj(imgForFilter, x, y, shape, true, true);
- int s = shape.size();
- if (s>0) {
- shapeVectors.push_back(shape);
- cout << "new object pushed back at position "<<x<<","<<y<<" of size "<<s<<endl;
- pair<int, int> si(s, objectCounter);
- sizeNIndexVector.push_back(si);
- objectCounter++;
- }
-
- }
- }
-
- // sort the sizes and take the largest to find average largest
- cout<<"shapeVectors size = "<<shapeVectors.size()<<endl;
-
- bubbleSort();
-
- cout <<"Largest object size is "<<sizeNIndexVector[0].first<<endl;
-
- // add the largest size to sum
- totalSize += static_cast<long>(sizeNIndexVector[0].first);
-
- cout << "Current total size is "<<totalSize<<" for object "<<(frameCounter+1)<<endl;
-
- cout << "-----------------------------------------------------------"<<endl;
- */
-
- frameCounter++;
-
- delete imgForFilter;
-
-// delete residual;
-
- delete img;
-
-
- }
-
- inputFile.close();
-
- /*
- avgLargestSize = static_cast<long> (totalSize/frameCounter);
- cout << "Average largest size is "<<avgLargestSize<<endl;
- */
-
- // previous loop calculates the average largest size
-
- inputFile.open(argv[1]);
-
- if (inputFile.fail() == true) {
- cout << "Cannot open the input file again"<<endl;
- exit(1);
- }
-
-
- while (inputFile>>fileName) {
-
- //string inputFileName = "output/filtered/temp/"+fileName;
- string inputFileName = outputFileLocation + "temp/"+fileName;
-
- Image* img = new Image(inputFileName.c_str());
- int width = img->columns();
- int height = img->rows();
-
- outputFile<<"File name is "<<inputFileName<<endl;
- outputFile<<"----------------------------------------------\n";
-
- sprintf(buffer,"%ix%i",width,height);
-
- // residual image is initialized with black representing not visited.
- residual = new Image(buffer, "black");
-
- cout << "Reading "<<inputFileName<<endl;
- cout << "Filter wxh "<<width<<","<<height<<endl;
-
- shapeVectors.clear();
- sizeNIndexVector.clear();
-
- // find the objects and sort according to size
- int objectCounter = 0;
- for (int x=0; x<width; x++) {
- for (int y=0; y<height; y++) {
- // find the white object using eight connected
- shape.clear();
- findObj(img, x, y, shape, true, true);
- int s = shape.size();
-
- if (s > 0) {
-
- shapeVectors.push_back(shape);
- cout << "New object found at position ("<<x<<","<<y<<") of size "<<s<<endl;
- pair<int, int> si(s, objectCounter);
- sizeNIndexVector.push_back(si);
- objectCounter++;
-
- }
-
- }
- }
-
- cout << "Total object found "<<sizeNIndexVector.size()<<endl;
-
- outputFile<<"Total objects found "<<sizeNIndexVector.size()<<endl;
-
- bubbleSort();
-
- // take the largest object
- double currentLargestSize = static_cast<double>(sizeNIndexVector[0].first);
-
- double secondLargest = 0;
- double ratio = 0;
- if (sizeNIndexVector.size() > 1) {
-
- secondLargest = static_cast<double>(sizeNIndexVector[1].first);
- ratio = secondLargest/currentLargestSize;
- cout << "Ratio is "<<secondLargest<<"/"<<currentLargestSize<<" = "<<ratio<<endl;
- outputFile<<"secondLargest = "<<secondLargest<<"\ncurrentLargest = "<<currentLargestSize<<"\nRatio = "<<ratio<<endl;
-
- }
-
- // find the largest to second largest ratio if it is less than the defined ratio then
- // the objects are single object
-
- int numberOfObjects = 0;
-
- if (sizeNIndexVector.size() == 1) {
- cout << "Frame contains one object\n";
- outputFile << "Frame contains one object\n";
- numberOfObjects = 1;
- }
- else if (ratio <= ratioSecondLargestToLargest ) {
- cout << "Single object calculated in the frame because current ratio = "<<ratio<<" is less than defined ratio = "<<ratioSecondLargestToLargest<<endl;
- outputFile<< "Single object calculated in the frame because current ratio = "<<ratio<<" is less than defined ratio = "<<ratioSecondLargestToLargest<<endl;
- numberOfObjects = 1;
- } else {
- cout << "Two objects in the frame\n";
- outputFile<<"Two objects in the frame\n";
- numberOfObjects = 2;
- }
-
- cout << "Total object found "<<numberOfObjects<<endl;
- outputFile << "Total object found "<<numberOfObjects<<endl;
-
- Image* imgFinal = new Image(buffer, "black");
-
- for (int n=0; n<numberOfObjects; n++) {
-
- int totalPoints = sizeNIndexVector[n].first;
-
- for (int i=0; i<totalPoints; i++) {
-
- imgFinal->pixelColor(shapeVectors[ sizeNIndexVector[n].second ][i].first, shapeVectors[ sizeNIndexVector[n].second ][i].second, "white");
- }
-
-
- }
-
- //string finalImageName = "output/filtered/final/"+fileName;
-
- string finalImageName = outputFileLocation + "final/"+fileName;
-
- imgFinal->write( finalImageName.c_str() );
-
-
-
-
- // writing the single in red
- if (numberOfObjects == 1) {
-
- Image* singleObjectFinal = new Image(buffer, "black");
-
- int totalPoints = sizeNIndexVector[0].first;
-
- cout << "Output the single object of size = "<<totalPoints<<endl;
-
- for (int i=0; i<totalPoints; i++) {
-
- singleObjectFinal->pixelColor(shapeVectors[ sizeNIndexVector[0].second ][i].first, shapeVectors[ sizeNIndexVector[0].second ][i].second, "red");
- }
-
- //string singleImageName = "output/filtered/single/"+fileName;
- string singleImageName = outputFileLocation + "single/"+fileName;
-
- singleObjectFinal->write(singleImageName.c_str());
-
- delete singleObjectFinal;
-
- }
-
- outputFile<<"----------------------------------------------------\n";
-
-
- delete img;
-
- delete residual;
-
- delete imgFinal;
-
-
- }
-
- inputFile.close();
- outputFile.close();
-
- return 0;
-}
-
-
-void findObj(Image* img, int x, int y, vector<pair<int,int> > & shape ,bool eightCon, bool colorLookingFor)
-{
- assert(residual != NULL);
-
- if (eightCon == true)
- eightConnObj(img, x, y, shape, colorLookingFor);
- else {
- fourConnObj(img, x, y, shape, colorLookingFor);
- }
-}
-
-int barrier = 1000;
-void fourConnObj(Image* img, int x, int y, vector<pair<int, int> > & obj, bool color)
-{
- int width = img->columns(),height = img->rows();
-
- // boundary violation check
- if ( (x >= (width)) || (x < 0) || (y >= (height) ) || (y < 0) )
- return;
-
- // residualpixel.mono() == true implies it is visited. Otherwise not visited.
- ColorMono residualpixel = ColorMono(residual->pixelColor(x,y));
- // originalpixel.mono() == true implies it is an object pixel. Otherwise it is blank region pixel.
- ColorMono originalpixel = ColorMono(img->pixelColor(x,y));
-
- // If the current pixel is already visited then return
- if (residualpixel.mono() == true)
- return;
-
- // Else if current pixel is not visited and it is black, which means it is not an object pixel; so return
- else if (residualpixel.mono() == false && originalpixel.mono() != color)
- return;
- // If current pixel is not visited and its value is white, which means a new object is found.
- else if (residualpixel.mono() == false && originalpixel.mono() == color) {
- // Save the coordinates of the current pixel into the vector and make the pixel visited in the residual image
- pair<int,int> p;
- p.first = x;
- p.second = y;
- obj.push_back(p);
-
-// if (obj.size() > barrier) {
-// cout<<obj.size()<<endl;
-// barrier = barrier + 1000;
-// }
- // setting the residual image at pixel(x,y) to white.
- residual->pixelColor(x,y, ColorMono(true));
-
- // Recursively call all of it's eight neighbours.
- fourConnObj(img, x+1, y, obj, color);
- fourConnObj(img, x, y-1, obj, color);
-
- fourConnObj(img, x-1, y, obj, color);
- fourConnObj(img, x, y+1, obj, color);
- }
-
-}
-
-void eightConnObj(Image* img, int x, int y, vector<pair<int, int> > & obj, bool color)
-{
- int width = img->columns(),height = img->rows();
-
- // boundary violation check
- if ( (x >= (width)) || (x < 0) || (y >= (height) ) || (y < 0) )
- return;
-
- // residualpixel.mono() == true implies it is visited. Otherwise not visited.
- ColorMono residualpixel = ColorMono(residual->pixelColor(x,y));
- // originalpixel.mono() == true implies it is an object pixel. Otherwise it is blank region pixel.
- ColorMono originalpixel = ColorMono(img->pixelColor(x,y));
-
- // If the current pixel is already visited then return
- if (residualpixel.mono() == true)
- return;
-
- // Else if current pixel is not visited and it is black, which means it is not an object pixel; so return
- else if (residualpixel.mono() == false && originalpixel.mono() != color)
- return;
- // If current pixel is not visited and its value is white, which means a new object is found.
- else if (residualpixel.mono() == false && originalpixel.mono() == color) {
- // Save the coordinates of the current pixel into the vector and make the pixel visited in the residual image
- pair<int,int> p;
- p.first = x;
- p.second = y;
- obj.push_back(p);
-
-// if (obj.size() > barrier) {
-// //cout<<obj.size()<<endl;
-// barrier = barrier + 1000;
-// }
- // setting the residual image at pixel(x,y) to white.
- residual->pixelColor(x,y, ColorMono(true));
-
- // Recursively call all of it's eight neighbours.
- eightConnObj(img, x+1, y, obj, color);
- eightConnObj(img, x+1, y-1, obj, color);
- eightConnObj(img, x, y-1, obj, color);
- eightConnObj(img, x-1, y-1, obj, color);
-
- eightConnObj(img, x-1, y, obj, color);
- eightConnObj(img, x-1, y+1, obj, color);
- eightConnObj(img, x, y+1, obj, color);
- eightConnObj(img, x+1, y+1, obj, color);
-
- }
-
-}
-
-// Aspect Ratio
-pair<int,int> getCentroid(vector<pair<int,int> > & points)
-{
- pair<int,int> centroid;
- centroid.first = 0;
- centroid.second = 0;
-
- for (unsigned int i = 0; i<points.size(); i++)
- {
- centroid.first += points[i].first;
- centroid.second += points[i].second;
- }
-
- centroid.first = roundT(double(centroid.first)/points.size());
- centroid.second = roundT(double(centroid.second)/points.size());
-
- return centroid;
-}
-
-
-vector<double> covariantDecomposition(vector<pair<int,int> > & points)
-{
- unsigned int i,j,k;
- pair<int,int> centroid = getCentroid(points);
- vector<double> retval;
-
- gsl_matrix* matrice = gsl_matrix_alloc(2, 2);
-
- double sumX2 = 0, sumXY = 0, sumY2 = 0;
- for (k = 0; k<points.size(); k++)
- {
- sumX2 += pow(double(points[k].first - centroid.first),2.0);
- sumY2 += pow(double(points[k].second - centroid.second),2.0);
- // should we take the absolute value of X*Y
- sumXY += (points[k].first - centroid.first) * (points[k].second - centroid.second);
- }
- gsl_matrix_set(matrice, 0, 0, roundT(sumX2/points.size()));
- gsl_matrix_set(matrice, 0, 1, roundT(sumXY/points.size()));
- gsl_matrix_set(matrice, 1, 0, roundT(sumXY/points.size()));
- gsl_matrix_set(matrice, 1, 1, roundT(sumY2/points.size()));
-
- // outputMatrix("Covariant", matrice);
-
- // This function allocates a workspace for computing eigenvalues of n-by-n
- // real symmetric matrices. The size of the workspace is O(2n).
- gsl_eigen_symmv_workspace* eigenSpace = gsl_eigen_symmv_alloc(2);
- gsl_vector* eigenVal = gsl_vector_alloc(2);
- gsl_matrix* eigenVec = gsl_matrix_alloc(2, 2);
- // This function computes the eigenvalues and eigenvectors of the real
- // symmetric matrix A. Additional workspace of the appropriate size must be
- // provided in w. The diagonal and lower triangular part of A are destroyed
- // during the computation, but the strict upper triangular part is not
- // referenced. The eigenvalues are stored in the vector eval and are unordered.
- // The corresponding eigenvectors are stored in the columns of the matrix evec.
- // For example, the eigenvector in the first column corresponds to the first
- // eigenvalue. The eigenvectors are guaranteed to be mutually orthogonal and
- // normalised to unit magnitude.
- gsl_eigen_symmv (matrice, eigenVal, eigenVec, eigenSpace);
- gsl_eigen_symmv_free (eigenSpace);
-
- gsl_eigen_symmv_sort(eigenVal, eigenVec, GSL_EIGEN_SORT_VAL_ASC);
-
- for (i = 0; i<eigenVal->size; i++)
- retval.push_back(gsl_vector_get(eigenVal, i));
-
- for (j = 0; j<eigenVec->size2; j++)
- for (i = 0; i<eigenVec->size1; i++)
- retval.push_back(gsl_matrix_get(eigenVec, i, j));
-
- retval.push_back(static_cast<double>(centroid.first));
- retval.push_back(static_cast<double> (centroid.second));
-
-// for (i=0; i<2; i++) {
-// gsl_vector_view evec_i = gsl_matrix_column (eigenVec, i);
-// //printf ("eigenvalue = %g\n", eval_i);
-// cout<<"eigenvector = \n";
-// gsl_vector_fprintf (stdout, &evec_i.vector, "%g");
-// }
-
- gsl_vector_free(eigenVal);
- gsl_matrix_free(matrice);
- gsl_matrix_free(eigenVec);
-
- return retval;
-}
-
-// isInterface for binary image
-bool isInterface(Image* orig, unsigned int x, unsigned int y)
-{
- // Get the current pixel's color
- ColorMono currentpixel = (ColorMono)orig->pixelColor(x,y);
- // If the current pixel is black pixel then it is not boundary pixel
- // error check
- if (currentpixel.mono() == false)
- return false;
-
- // If the current pixel is not black then it is white. So, now we need
- // to check whether any four of its neighbor pixels (left, top, right,
- // bottom ) is black. If any of this neighbor is black then current
- // pixel is a neighbor pixel. Otherwise current pixel is not neighbor
- // pixel.
-
- ColorMono leftneighborpixel = (ColorMono)orig->pixelColor(x-1,y);
- ColorMono topneighborpixel = (ColorMono)orig->pixelColor(x,y-1);
- ColorMono rightneighborpixel = (ColorMono)orig->pixelColor(x+1,y);
- ColorMono bottomneighborpixel = (ColorMono)orig->pixelColor(x,y+1);
-
- // If leftneighborpixel is black and currentpixel is white then it is
- // boundary pixel
- if ( leftneighborpixel.mono() != currentpixel.mono())
- return true;
- // If topneighborpixel is black and currentpixel is white then it is
- // boundary pixel
- else if (topneighborpixel.mono() != currentpixel.mono())
- return true;
- // If rightneighborpixel is black and currentpixel is white then it
- // is boundary pixel
- else if (rightneighborpixel.mono() != currentpixel.mono())
- return true;
- // If bottomneighborpixel is black and currentpixel is white then it
- // is boundary pixel
- else if (bottomneighborpixel.mono() != currentpixel.mono())
- return true;
- // Else all of its neighbor pixels are white so it can not be a
- // boundary pixel
- else
- return false;
-
-} \ No newline at end of file