API tools faq
paste
Advertisement
Guest User

Untitled

a guest
Jan 18th, 2017
74
0
Never
Add comment
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
C++ 21.15 KB | None | 0 0
raw download clone embed print report
  1.  
  2. #include "stdafx.h"
  3. #define _USE_MATH_DEFINES
  4. #include <iostream>
  5. #include <cmath>
  6. #include <string>
  7. #include <vector>
  8. #include <complex>
  9. #include "CImg.h"
  10.  
  11. using namespace std;
  12. using namespace cimg_library;
  13.  
  14. void DFT(const CImg<float>&image, vector<vector<complex<double>>>&factors)
  15. {
  16.     double width = image.width();
  17.     double height = image.height();
  18.     double angle;
  19.     for (int i = 0; i < width; i++)
  20.     {
  21.         for (int j = 0; j < height; j++)
  22.         {
  23.             complex<double> sum(0, 0);
  24.             complex<double> Wn(0, 0);
  25.             complex<double> Wm(0, 0);
  26.             for (int x = 0; x < width; x++)
  27.             {
  28.                 for (int y = 0; y < height; y++)
  29.                 {
  30.                     Wn.real(cos(2 * M_PI * x * i / width));
  31.                     Wn.imag(-sin(2 * M_PI * x * i / width));
  32.                     Wm.real(cos(2 * M_PI * y * j / height));
  33.                     Wm.imag(-sin(2 * M_PI * y * j / height));
  34.                     complex<double> valuePix(double(image(x, y, 0)));
  35.                     sum += (valuePix * Wn * Wm);
  36.                 }
  37.             }
  38.             factors[i][j] = sum / (sqrt(width*height)); //normalization
  39.         }
  40.     }
  41. }
  42.  
  43. void IDFT(CImg<float> &image, const vector<vector<complex<double>>>&factors)
  44. {
  45.     double width = image.width();
  46.     double height = image.height();
  47.  
  48.     for (int j = 0; j < height; j++)
  49.     {
  50.         for (int i = 0; i < width; i++)
  51.         {
  52.             complex<double> suma(0, 0);
  53.             complex<double> Wn(0, 0);
  54.             complex<double> Wm(0, 0);
  55.             for (int y = 0; y < height; y++)
  56.             {
  57.                 for (int x = 0; x < width; x++)
  58.                 {
  59.                     Wn.real(cos(2 * M_PI * i * x / width));
  60.                     Wn.imag(sin(2 * M_PI * i * x / width));
  61.                     Wm.real(cos(2 * M_PI * j * y / height));
  62.                     Wm.imag(sin(2 * M_PI * j * y / height));
  63.                     suma += (factors[x][y] * Wn * Wm);
  64.                 }
  65.             }
  66.             image(i, j) = suma.real() / sqrt(width*height); //normalization!
  67.         }
  68.     }
  69. }
  70.  
  71. void changePosition(CImg<unsigned char> &image)
  72. {
  73.     const int width = image.width();
  74.     const int height = image.height();
  75.     CImg<unsigned char> temp(width, height, 1, 3, 0);
  76.     for (int j = height / 2; j < height; j++)
  77.         for (int i = 0; i < width / 2; i++)
  78.         {
  79.             temp((width / 2) + i, j - (height / 2)) = image(i, j);
  80.         }
  81.     for (int j = 0; j < height / 2; j++)
  82.         for (int i = 0; i < width / 2; i++)
  83.         {
  84.             temp((width / 2) + i, (height / 2) + j) = image(i, j);
  85.         }
  86.     for (int j = 0; j < height / 2; j++)
  87.         for (int i = width / 2; i < width; i++)
  88.         {
  89.             temp(i - (width / 2), (height / 2) + j) = image(i, j);
  90.         }
  91.     for (int j = height / 2; j < height; j++)
  92.         for (int i = width / 2; i < width; i++)
  93.         {
  94.             temp(i - (width / 2), j - (width / 2)) = image(i, j);
  95.         }
  96.     image = temp;
  97. }
  98.  
  99. CImg<unsigned char> getMagnitude(vector<vector<complex<double>>> &factors)
  100. {
  101.     const int width = factors.size();
  102.     const int height = factors[0].size();
  103.     int magnitudeMax = 0;
  104.     int presentMag = 0;
  105.     CImg<unsigned char> magnitude(width, height, 1, 3, 0);
  106.     for (int x = 0; x < width; x++)
  107.         for (int y = 0; y < height; y++)
  108.         {
  109.             float re = factors[x][y].real();
  110.             float im = factors[x][y].imag();
  111.             presentMag = sqrt((re * re) + (im * im));
  112.             if (presentMag > magnitudeMax) magnitudeMax = presentMag;   //seting max
  113.         }
  114.     for (int x = 0; x < width; x++)
  115.         for (int y = 0; y < height; y++)
  116.         {
  117.             float re = factors[x][y].real();
  118.             float im = factors[x][y].imag();
  119.             magnitude(x, y) = (255 / log(1.0 + magnitudeMax)) * log(1.0 + presentMag);
  120.         }
  121.     changePosition(magnitude);
  122.     return magnitude;
  123. }
  124.  
  125. CImg<unsigned char> getPhase(vector<vector<complex<double>>> &factors)
  126. {
  127.     const int width = factors.size();
  128.     const int height = factors[0].size();
  129.     CImg<unsigned char> phase(width, height, 1, 3, 0);
  130.  
  131.     for (int i = 0; i<width; i++)
  132.         for (int j = 0; j<height; j++)
  133.         {
  134.             float re = factors[i][j].real();
  135.             float im = factors[i][j].imag();
  136.             phase(i, j) = atan(im / re);
  137.             phase(i, j) = (phase(i, j)) * (255 / (2 * M_PI));
  138.         }
  139.     changePosition(phase);
  140.     return phase;
  141. }
  142.  
  143. vector< complex<double> > fft1d(vector< complex<double> > input)
  144. {
  145.     const int N = input.size();
  146.     const complex<double> I = (0.0, 1.0);
  147.     if (N == 1) return input;
  148.     else {
  149.  
  150.         vector< complex<double> > left, right, temp1, temp2, temp3, temp4, temp5;
  151.  
  152.         for (int i = 0; i < (N / 2); ++i)
  153.         {
  154.             left.push_back(input[i]);
  155.             right.push_back((input[N / 2 + i]));
  156.         }
  157.         for (int n = 0; n < N / 2; ++n)
  158.         {
  159.             temp1.push_back(left[n] + right[n]);
  160.             temp2.push_back(left[n] - right[n]);
  161.         }
  162.         //fourierShift(temp2, inverse);
  163.         complex<double> alfa;
  164.         complex<double> t;
  165.         for (int k = 0; k < N / 2; ++k)
  166.         {
  167.             alfa = ((2 * M_PI * k) / N);
  168.             t = (cos(alfa) - I*sin(alfa))* temp2[k];
  169.             temp3.push_back(t);
  170.         }
  171.         temp4 = fft1d(temp1);
  172.         temp5 = fft1d(temp3);
  173.         int j = 0;
  174.         for (int k = 0; k<N / 2; ++k)
  175.         {
  176.             input[j] = temp4[k];
  177.             input[j + 1] = temp5[k];
  178.             j += 2;
  179.         }
  180.         return input;
  181.     }
  182. }
  183. void normalization(vector< complex<double> > &c)
  184. {
  185.     for (auto i = c.begin(); i < c.end(); ++i)
  186.     {
  187.         (*i) /= sqrt(c.size());
  188.     }
  189. }
  190.  
  191. vector<vector<complex<double>>> FFT(vector<vector<complex<double>>> image)
  192. {
  193.  
  194.     for (size_t i = 0; i < image[0].size(); i++)
  195.         fft1d(image[i]);
  196.  
  197.     //swap(input);
  198.  
  199.     for (size_t j = 0; j < image.size(); j++)
  200.         fft1d(image[j]);
  201.  
  202.     //swap(input);
  203.  
  204.     return image;
  205. }
  206. void IFFT(CImg<unsigned char> &image, vector<vector<complex<double>>>&factors)
  207. {
  208.     const int width = image.width();
  209.     const int height = image.height();
  210.  
  211.     //unsigned int *reverse1 = new unsigned int[width];
  212.     //for (int i = 0; i<width; i++) reverse1[i] = i;
  213.     //reverseBinary(width, reverse1);
  214.  
  215.     //const double PI = acos(-1.0);
  216.  
  217.     //bool *swapped = new bool[width];
  218.     //for (int i = 0; i<width; i++)
  219.     //  swapped[i] = false;
  220.     //for (int i = 0; i<width; i++)
  221.     //{
  222.     //  for (int j = 0; j<height; j++)
  223.     //      if (!(swapped[i] || swapped[reverse1[i]]))
  224.     //          swap(factors[i][j], factors[reverse1[i]][j]);
  225.     //  if (i != reverse1[i])
  226.     //  {
  227.     //      swapped[i] = true;
  228.     //      swapped[reverse1[i]] = true;
  229.     //  }
  230.     //}
  231.     //delete[] swapped;
  232.  
  233.     //swapped = new bool[height];
  234.     //for (int i = 0; i<width; i++)
  235.     //  swapped[i] = false;
  236.     //for (int j = 0; j<height; j++)
  237.     //{
  238.     //  for (int i = 0; i<width; i++)
  239.     //      if (!(swapped[j] || swapped[reverse1[j]]))
  240.     //          swap(factors[i][j], factors[i][reverse1[j]]);
  241.     //  if (j != reverse1[j])
  242.     //  {
  243.     //      swapped[j] = true;
  244.     //      swapped[reverse1[j]] = true;
  245.     //  }
  246.     //}
  247.     //delete[] swapped;
  248.  
  249.     //for (int j = 0; j<height; j++)
  250.     //{
  251.     //  complex<double> *W = new complex<double>[width / 2];
  252.     //  for (int n = 0; n<width / 2; n++)
  253.     //  {
  254.     //      W[n].real(cos(2 * PI*n / width));
  255.     //      W[n].imag(sin(2 * PI*n / width));
  256.     //  }
  257.  
  258.     //  const int levelMax = log2(width);
  259.  
  260.     //  for (int level = 0; level < levelMax; level++)
  261.     //  {
  262.     //      const int blockMax = width / pow(2.0, level + 1);
  263.  
  264.     //      for (int block = 0; block < blockMax; block++)
  265.     //      {
  266.     //          const int operationMax = pow(2.0, level);
  267.  
  268.     //          for (int operation = 0; operation < operationMax; operation++)
  269.     //          {
  270.     //              const int i = operation + (block * 2 * pow(2.0, level));
  271.  
  272.     //              const complex<double> a = factors[i][j];
  273.     //              complex<double> b = factors[i + operationMax][j];
  274.     //              b *= W[operation*(int)pow(2.0, levelMax - level - 1)];
  275.     //              factors[i][j] = a + b;
  276.     //              factors[i + operationMax][j] = a - b;
  277.     //          }
  278.     //      }
  279.     //  }
  280.     //  delete[] W;
  281.     //}
  282.  
  283.     //for (int i = 0; i<width; i++)
  284.     //{ //tweedle factor
  285.     //  complex<double> *W = new complex<double>[height / 2];
  286.     //  for (int m = 0; m<height / 2; m++)
  287.     //  {
  288.     //      W[m].real(cos(2 * PI*m / height));
  289.     //      W[m].imag(sin(2 * PI*m / height));
  290.     //  }
  291.  
  292.     //  const int levelMax = log2(height);
  293.  
  294.     //  for (int level = 0; level < levelMax; level++)
  295.     //  {
  296.     //      const int blockMax = height / pow(2.0, level + 1);
  297.  
  298.     //      for (int block = 0; block < blockMax; block++)
  299.     //      {
  300.     //          const int operationMax = pow(2.0, level);
  301.  
  302.     //          for (int operation = 0; operation < operationMax; operation++)
  303.     //          {
  304.     //              const int j = operation + (block * 2 * pow(2.0, level));
  305.  
  306.     //              const complex<double> a = factors[i][j];
  307.     //              complex<double> b = factors[i][j + operationMax];
  308.     //              b *= W[operation*(int)pow(2.0, levelMax - level - 1)];
  309.     //              factors[i][j] = a + b;
  310.     //              factors[i][j + operationMax] = a - b;
  311.     //          }
  312.     //      }
  313.     //  }
  314.     //  delete[] W;
  315.     //}
  316.  
  317.     //for (int j = 0; j<height; j++)
  318.     //  for (int i = 0; i<width; i++)
  319.     //      for (int k = 0; k<3; k++)
  320.     //      {
  321.     //          factors[i][j] /= sqrt(width*height + 0.0);
  322.     //          image(i, j) = trunc(factors[i][j].real());
  323.     //      }
  324.     //delete[] reverse1;
  325. }
  326.  
  327. void lowpass(CImg<float> &image, int range)
  328. {
  329.     vector<vector<complex<double>>> factors;
  330.     int width = image.width();
  331.     int height = image.height();
  332.     CImg<unsigned char> magnitude(width, height, 1, 3, 0);
  333.     CImg<unsigned char> phase(width, height, 1, 3, 0);
  334.     double distance;
  335.     CImg<unsigned char> mask(width, height, 1, 3, 0);
  336.     CImg<unsigned char> outputimage(width, height, 1, 3, 0);
  337.  
  338.     factors.resize(width);
  339.     for (int i = 0; i<width; i++)
  340.     {
  341.         factors[i].resize(height);
  342.         for (int j = 0; j<height; j++)
  343.         {
  344.                 factors[i][j]= 0;
  345.         }
  346.     }
  347.  
  348.     for (int i = 0; i < width; i++)
  349.     {
  350.         for (int j = 0; j < height; j++)
  351.         {
  352.             distance = sqrt((float)((i - width / 2)*(i - width / 2)) + ((j - height / 2)*(j - height / 2)));
  353.             if (distance < range)
  354.             {
  355.                 mask(i, j) = 1;
  356.             }
  357.             else
  358.             {
  359.                 mask(i, j) = 0;
  360.             }
  361.         }
  362.     }
  363.  
  364.     changePosition(mask);
  365.     //FFT(image);
  366.  
  367.     for (int i = 0; i < width; i++)
  368.     {
  369.         for (int j = 0; j < height; j++)
  370.         {
  371.                 factors[i][j].real(factors[i][j].real() * mask(i, j));
  372.                 factors[i][j].imag(factors[i][j].imag() * mask(i, j));
  373.         }
  374.     }
  375.  
  376.     magnitude = getMagnitude(factors);
  377.     phase = getPhase(factors);
  378.  
  379.     IFFT(outputimage, factors);
  380.  
  381.     outputimage.save("lowpass.bmp");
  382.     magnitude.save("magnitude_low.bmp");
  383.     phase.save("phase_low.bmp");
  384. }
  385.  
  386. void highpass(CImg<float> &image, double range)
  387. {
  388.  
  389.     vector<vector<complex<double>>>factors;
  390.  
  391.     CImg<unsigned char> magnitude(image.width(), image.height(), 1, 3, 0);
  392.     CImg<unsigned char> phase(image.width(), image.height(), 1, 3, 0);
  393.  
  394.     factors.resize(image.width());
  395.     for (int x = 0; x<image.width(); x++)
  396.     {
  397.         factors[x].resize(image.height());
  398.         for (int y = 0; y<image.height(); y++)
  399.         {
  400.                 factors[x][y] = 0;
  401.         }
  402.     }
  403.  
  404.     int width = image.width();
  405.     int height = image.height();
  406.     double distance;
  407.     CImg<unsigned char> mask(width, height, 1, 3, 0);
  408.     CImg<unsigned char> image_output(width, height, 1, 3, 0);
  409.  
  410.  
  411.     for (int x = 0; x < width; x++)
  412.     {
  413.         for (int y = 0; y < height; y++)
  414.         {
  415.                 distance = sqrt((float)((x - width / 2)*(x - width / 2)) + ((y - height / 2)*(y - height / 2)));
  416.                 if (distance > range)
  417.                     mask(x, y) = 1;
  418.                 else
  419.                     mask(x, y) = 0;
  420.         }
  421.     }
  422.  
  423.     changePosition(mask);
  424.     mask(0, 0, 0) = 1;
  425.     mask(0, 0, 1) = 1;
  426.     mask(0, 0, 2) = 1;
  427.  
  428.     //FFT(image, factors);
  429.  
  430.     for (int x = 0; x < width; x++)
  431.     {
  432.         for (int y = 0; y < height; y++)
  433.         {
  434.                 factors[x][y].real(factors[x][y].real() * mask(x, y));
  435.                 factors[x][y].imag(factors[x][y].imag() * mask(x, y));
  436.         }
  437.     }
  438.  
  439.     magnitude = getMagnitude(factors);
  440.     phase = getPhase(factors);
  441.  
  442.     IFFT(image_output, factors);
  443.  
  444.     image_output.save("highpass.bmp");
  445.     magnitude.save("magnitude_high.bmp");
  446.     phase.save("phase_high.bmp");
  447.  
  448. }
  449.  
  450. void bandpass(CImg<float> &image, int range, int thickness)
  451. {
  452.     vector<vector<complex<double>>> factors;
  453.     int width = image.width();
  454.     int height = image.height();
  455.     CImg<unsigned char> magnitude(width, height, 1, 3, 0);
  456.     CImg<unsigned char> phase(width, height, 1, 3, 0);
  457.     double distance;
  458.     CImg<unsigned char> mask(width, height, 1, 3, 0);
  459.     CImg<unsigned char> outputimage(width, height, 1, 3, 0);
  460.  
  461.     factors.resize(width);
  462.     for (int i = 0; i<width; i++)
  463.     {
  464.         factors[i].resize(height);
  465.         for (int j = 0; j<height; j++)
  466.         {
  467.                 factors[i][j] = 0;
  468.         }
  469.     }
  470.  
  471.     for (int i = 0; i < width; i++)
  472.     {
  473.         for (int j = 0; j < height; j++)
  474.         {
  475.             distance = sqrt((float)((i - width / 2)*(i - width / 2)) + ((j - height / 2)*(j - height / 2)));
  476.             if ((distance > range - (thickness / 2)) && (distance < range + (thickness / 2)))
  477.             {
  478.                 mask(i, j) = 1;
  479.             }
  480.             else
  481.             {
  482.                 mask(i, j) = 0;
  483.             }
  484.         }
  485.     }
  486.  
  487.     changePosition(mask);
  488.     mask(0, 0, 0) = 1;
  489.     mask(0, 0, 1) = 1;
  490.     mask(0, 0, 2) = 1;
  491.  
  492.     //FFT(image, factors);
  493.  
  494.     for (int i = 0; i < width; i++)
  495.     {
  496.         for (int j = 0; j < height; j++)
  497.         {
  498.                 factors[i][j].real(factors[i][j].real() * mask(i, j));
  499.                 factors[i][j].imag(factors[i][j].imag() * mask(i, j));
  500.         }
  501.     }
  502.  
  503.     magnitude = getMagnitude(factors);
  504.     phase = getPhase(factors);
  505.  
  506.     IFFT(outputimage, factors);
  507.  
  508.     outputimage.save("bandpass.bmp");
  509.     magnitude.save("magnitude_bpass.bmp");
  510.     phase.save("phase_bpass.bmp");
  511. }
  512.  
  513. void bandcut(CImg<float> &image, int range, int thickness)
  514. {
  515.     vector<vector<complex<double>>>factors;
  516.     int width = image.width();
  517.     int height = image.height();
  518.     CImg<unsigned char> magnitude(width, height, 1, 3, 0);
  519.     CImg<unsigned char> phase(width, height, 1, 3, 0);
  520.     double distance;
  521.     CImg<unsigned char> mask(width, height, 1, 3, 0);
  522.     CImg<unsigned char> outputimage(width, height, 1, 3, 0);
  523.  
  524.     factors.resize(width);
  525.     for (int i = 0; i<width; i++)
  526.     {
  527.         factors[i].resize(height);
  528.         for (int j = 0; j<height; j++)
  529.         {
  530.             factors[i][j] = 0;
  531.         }
  532.     }
  533.  
  534.     for (int i = 0; i < width; i++)
  535.     {
  536.         for (int j = 0; j < height; j++)
  537.         {
  538.             distance = sqrt((float)((i - width / 2)*(i - width / 2)) + ((j - height / 2)*(j - height / 2)));
  539.             if ((distance < range - (thickness / 2)) || (distance > range + (thickness / 2)))
  540.             {
  541.                 mask(i, j) = 1;
  542.             }
  543.             else
  544.             {
  545.                 mask(i, j) = 0;
  546.             }
  547.         }
  548.     }
  549.  
  550.     changePosition(mask);
  551.  
  552.     //FFT(image, factors);
  553.  
  554.     for (int i = 0; i < width; i++)
  555.     {
  556.         for (int j = 0; j < height; j++)
  557.         {
  558.                 factors[i][j].real(factors[i][j].real() * mask(i, j));
  559.                 factors[i][j].imag(factors[i][j].imag() * mask(i, j));
  560.         }
  561.     }
  562.  
  563.     magnitude = getMagnitude(factors);
  564.     phase = getPhase(factors);
  565.  
  566.     IFFT(outputimage, factors);
  567.  
  568.     outputimage.save("bandcut.bmp");
  569.     magnitude.save("magnitude_bcut.bmp");
  570.     phase.save("phase_bcut.bmp");
  571. }
  572.  
  573. CImg<unsigned char>& resize2(CImg<unsigned char>& image, double value)
  574. {
  575.     int width_after = static_cast<int>(image.width() * value);
  576.     int height_after = static_cast<int>(image.height() * value);
  577.     CImg<float> temp(width_after, height_after, 1, 3, 0); // (width, height, depth, channels, default vale)
  578.  
  579.     for (int i = 0; i < width_after; i++)
  580.     {
  581.         for (int j = 0; j < height_after; j++)
  582.         {
  583.                 temp(i, j) = image(i / value, j / value);
  584.         }
  585.     }
  586.     image = temp;
  587.     return image;
  588. }
  589.  
  590. vector<vector<complex<double>>> Vec2Img(CImg<double>image)
  591. {
  592.     vector<vector<complex<double>>> result;
  593.     for (int x = 0; x < image.width(); x++)
  594.     {
  595.         vector<complex<double>>temp;
  596.         for (int y = 0; y < image.height(); y++)
  597.         {
  598.             temp.push_back(complex<double>(image(x, y), 0.0));
  599.         }
  600.         result.push_back(temp);
  601.     }
  602.     return result;
  603. }
  604.  
  605.  
  606. void high_edge(CImg<float> &image, int number, int rotate_indx)
  607. {
  608.     int width = image.width();
  609.     int height = image.height();
  610.     CImg<unsigned char> outputimage(image.width(), image.height(), 1, 3, 0);
  611.     CImg<unsigned char> mask;
  612.     vector<vector<complex<double>>> factors;
  613.     CImg<unsigned char> magnitude(width, height, 1, 3, 0);
  614.     CImg<unsigned char> phase(width, height, 1, 3, 0);
  615.  
  616.     if (number == 0)
  617.     {
  618.         mask.load("F5mask1.bmp");
  619.     }
  620.     else if (number == 1)
  621.     {
  622.         mask.load("F5mask2.bmp");
  623.     }
  624.  
  625.     double scale = width / mask.width();
  626.     resize2(mask, scale);
  627.  
  628.     mask.save("maskF5.bmp");
  629.  
  630.     changePosition(mask);
  631.     mask(0, 0) = 255;
  632.  
  633.     factors.resize(width);
  634.     for (int i = 0; i<width; i++)
  635.     {
  636.         factors[i].resize(height);
  637.         for (int j = 0; j<height; j++)
  638.         {
  639.                 factors[i][j] = 0;
  640.         }
  641.     }
  642.  
  643.     //FFT(image, factors);
  644.  
  645.     for (int i = 0; i < width; i++)
  646.     {
  647.         for (int j = 0; j < height; j++)
  648.         {
  649.                 factors[i][j].real(factors[i][j].real() * (mask(i, j) / 255));
  650.                 factors[i][j].imag(factors[i][j].imag() * (mask(i, j) / 255));
  651.         }
  652.     }
  653.  
  654.     magnitude = getMagnitude(factors);
  655.     phase = getPhase(factors);
  656.  
  657.     IFFT(outputimage, factors);
  658.  
  659.     outputimage.save("high_edge.bmp");
  660.     magnitude.save("magnitude_highedge.bmp");
  661.     phase.save("phase_highedge.bmp");
  662.  
  663. }
  664.  
  665. void phase(CImg<float> &image, int q, int w)
  666. {
  667.     int width = image.width();
  668.     int height = image.height();
  669.     CImg<unsigned char> outputimage(width, height, 1, 3, 0);
  670.     vector<vector<complex<double>>> mask;
  671.     vector<vector<complex<double>>> factors;
  672.     CImg<unsigned char> magnitude(width, height, 1, 3, 0);
  673.     CImg<unsigned char> phase(width, height, 1, 3, 0);
  674.     const double PI = acos(-1.0);
  675.  
  676.     mask.resize(width);
  677.     for (int i = 0; i<width; i++)
  678.     {
  679.         mask[i].resize(height);
  680.         for (int j = 0; j<height; j++)
  681.         {
  682.             mask[i][j].real(cos(((-i*q * 2 * PI) / width) + ((-j*w * 2 * PI) / height) + (q + w)*PI));
  683.             mask[i][j].imag(sin(((-i*q * 2 * PI) / width) + ((-j*w * 2 * PI) / height) + (q + w)*PI));
  684.         }
  685.     }
  686.  
  687.     factors.resize(width);
  688.     for (int i = 0; i<width; i++)
  689.     {
  690.         factors[i].resize(height);
  691.         for (int j = 0; j<height; j++)
  692.         {
  693.             factors[i][j] = 0;
  694.         }
  695.     }
  696.  
  697.     //FFT(image, factors);
  698.  
  699.     for (int i = 0; i < width; i++)
  700.     {
  701.         for (int j = 0; j < height; j++)
  702.         {
  703.             factors[i][j] *= mask[i][j];
  704.         }
  705.     }
  706.  
  707.     magnitude = getMagnitude(factors);
  708.     phase = getPhase(factors);
  709.  
  710.     IFFT(outputimage, factors);
  711.  
  712.     outputimage.save("phase_filter.bmp");
  713.     magnitude.save("magnitude_phasefilter.bmp");
  714.     phase.save("phase_phasefilter.bmp");
  715.  
  716. }
  717.  
  718. int main(int argc, char *argv[])
  719. {
  720.     const double PI = acos(-1.0);
  721.     vector<vector<complex<double>>> factors;
  722.  
  723.     try
  724.     {
  725.         if (argc < 2) throw "Too few arguments! See --help";
  726.         /////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
  727.         if ((string)argv[1] == "--DFT")
  728.         {
  729.             if (argc != 3)
  730.                 throw "Incorrect/missing arguments! See --help";
  731.             CImg<float> img;
  732.             img.assign(argv[2]);
  733.             //--------------------------------------------------------------------
  734.             // making 3-dimensional table to store DFT values
  735.             factors.resize(img.width());
  736.             for (int x = 0; x<img.width(); x++)
  737.             {
  738.                 factors[x].resize(img.height());
  739.                 for (int y = 0; y<img.height(); y++)
  740.                 {
  741.                     factors[x][y] = 0;
  742.                 }
  743.             }
  744.             //--------------------------------------------------------------------
  745.             DFT(img, factors); // calculate factors
  746.                                //--------------------------------------------------------------------
  747.                                // save calculated magnitudes and phases to corresponding images
  748.             CImg<unsigned char> magnitude(img.width(), img.height(), 1, 3, 0);
  749.             CImg<unsigned char> phase(img.width(), img.height(), 1, 3, 0);
  750.             magnitude = getMagnitude(factors);
  751.             phase = getPhase(factors);
  752.             //--------------------------------------------------------------------
  753.             IDFT(img, factors);
  754.             //--------------------------------------------------------------------
  755.             magnitude.save("magnitude.bmp");
  756.             phase.save("phase.bmp");
  757.             img.save("inverse.bmp");
  758.             return 0;
  759.         }
  760.         /////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
  761.         if ((string)argv[1] == "--FFT")
  762.         {
  763.             if (argc != 3)
  764.                 throw "Incorrect/missing arguments! See --help";
  765.             CImg<unsigned char> img;
  766.             vector<vector<complex<double>>>res;
  767.             img.assign(argv[2]);
  768.             //--------------------------------------------------------------------
  769.             // making 3-dimensional table to store DFT values
  770.             factors.resize(img.width());
  771.             for (int x = 0; x<img.width(); x++)
  772.             {
  773.                 factors[x].resize(img.height());
  774.                 for (int y = 0; y<img.height(); y++)
  775.                 {
  776.                     factors[x][y] = 0;
  777.                 }
  778.             }
  779.             res = Vec2Img(img);
  780.             //--------------------------------------------------------------------
  781.             FFT(res); // calculate factors
  782.                                //--------------------------------------------------------------------
  783.                                // save calculated magnitudes and phases to corresponding images
  784.             CImg<unsigned char> magnitude(img.width(), img.height(), 1, 3, 0);
  785.             CImg<unsigned char> phase(img.width(), img.height(), 1, 3, 0);
  786.             magnitude = getMagnitude(res);
  787.             phase = getPhase(res);
  788.             //--------------------------------------------------------------------
  789.             //IFFT(img, factors);
  790.             //--------------------------------------------------------------------
  791.             magnitude.save("magnitude.bmp");
  792.             phase.save("phase.bmp");
  793.             //img.save("inverse.bmp");
  794.             return 0;
  795.         }
  796.  
  797.         if (argc == 4)
  798.         {
  799.             CImg<float> foto;
  800.             std::string par = reinterpret_cast<char*>(argv[1]);
  801.             foto.load(argv[2]);
  802.             int par2 = atoi(argv[3]);
  803.  
  804.             if (par == "--highpass")
  805.             {
  806.                 highpass(foto, par2);
  807.                 return 0;
  808.             }
  809.  
  810.             if (par == "--lowpass")
  811.             {
  812.                 lowpass(foto, par2);
  813.                 return 0;
  814.             }
  815.         }
  816.  
  817.         if (argc == 5)
  818.         {
  819.             CImg<float> foto;
  820.             std::string par = reinterpret_cast<char*>(argv[1]);
  821.             foto.load(argv[2]);
  822.             int par2 = atoi(argv[3]);
  823.             int par3 = atoi(argv[4]);
  824.  
  825.             if (par == "--bandpass")
  826.             {
  827.                 bandpass(foto, par2, par3);
  828.                 return 0;
  829.             }
  830.  
  831.             if (par == "--bandcut")
  832.             {
  833.                 bandcut(foto, par2, par3);
  834.                 return 0;
  835.             }
  836.  
  837.             if (par == "--highedge")
  838.             {
  839.                 high_edge(foto, par2, par3);
  840.                 return 0;
  841.             }
  842.  
  843.             if (par == "--phase")
  844.             {
  845.                 phase(foto, par2, par3);
  846.                 return 0;
  847.             }
  848.         } throw "Incorrect argument!";
  849.  
  850.     }
  851.     catch (char const *e)
  852.     {
  853.         cout << e;
  854.         return 1;
  855.     }
  856.     catch (CImgException &e)
  857.     {
  858.         cout << e.what();
  859.         return 2;
  860.     }
  861.  
  862.     return 0;
  863. }
Advertisement
Add Comment
Please, Sign In to add comment
Advertisement
Public Pastes
Advertisement
We use cookies for various purposes including analytics. By continuing to use Pastebin, you agree to our use of cookies as described in the Cookies Policy.  OK, I Understand
Not a member of Pastebin yet?
Sign Up, it unlocks many cool features!