]> sjero.net Git - qpsnr/commitdiff
Add siti computation support. Note that this is extremely hacky right now.
authorSamuel Jero <sj323707@ohio.edu>
Tue, 21 May 2013 23:34:03 +0000 (19:34 -0400)
committerSamuel Jero <sj323707@ohio.edu>
Tue, 21 May 2013 23:34:03 +0000 (19:34 -0400)
src/stats.cpp

index c537caf44c3a010c3f72596b718b784573608d74..fc0c968e601e5a429c2908b5fb1cc4a8365891d1 100644 (file)
@@ -21,6 +21,9 @@
 #include "mt.h"
 #include "shared_ptr.h"
 #include "settings.h"
+#include <stdio.h>
+#include <climits>
+#include <string.h>
 #include <cmath>
 #include <string>
 #include <stdexcept>
@@ -93,6 +96,89 @@ namespace stats {
                return avg/ssim_accum.size();
        }
 
+       double std_dev(int *buffer, int numelem)
+       {
+               // find sum and mean
+               long sum = 0;
+               for (int i=0; i<numelem; i++) sum+=buffer[i];
+               int mean = sum/numelem;
+
+               // sum of the squares
+               long long sum_squares = 0;
+               for (int i=0; i<numelem; i++)
+                       sum_squares += (long long)pow((float)(buffer[i]-mean),2);
+
+               // return the square root
+               return sqrt((double)(sum_squares/numelem));
+       }
+
+       void conv_sobel(int *buffer, int width, int height)
+       {
+               int mask1[9] = {-1,-2,-1,0,0,0,1,2,1}; // TODO right values
+               int mask2[9] = {-1,0,1,-2,0,2,-1,0,1};
+
+               // convolution
+               int buffsize = width*height;
+               int *buffer_tmp = (int *)malloc(buffsize*sizeof(int));
+               memcpy(buffer_tmp, buffer, buffsize*sizeof(int));
+               for (int l=1; l<height-1; l++) {
+               for (int c=1; c<width-1; c++) {
+
+                       int p1 = (l-1)*width+(c-1);     // pos first pixel
+                       int p2 = l*width+(c-1);         // pos 2nd line first pixel
+                       int p3 = (l+1)*width+(c-1);     // pos 3nd line first pixel
+
+                       // apply the mask and set the result
+                       long gv = 0, gh = 0;
+                       gv += buffer[p1]*mask1[0] + buffer[p1+1]*mask1[1] + buffer[p1+2]*mask1[2];
+                       gv += buffer[p2]*mask1[3] + buffer[p2+1]*mask1[4] + buffer[p2+2]*mask1[5];
+                       gv += buffer[p3]*mask1[6] + buffer[p3+1]*mask1[7] + buffer[p3+2]*mask1[8];
+                       gh += buffer[p1]*mask2[0] + buffer[p1+1]*mask2[1] + buffer[p1+2]*mask2[2];
+                       gh += buffer[p2]*mask2[3] + buffer[p2+1]*mask2[4] + buffer[p2+2]*mask2[5];
+                       gh += buffer[p3]*mask2[6] + buffer[p3+1]*mask2[7] + buffer[p3+2]*mask2[8];
+
+                       buffer_tmp[l*width+c] = sqrt(pow((double)gv,2)+pow((double)gh,2));
+               }
+               }
+               memcpy(buffer, buffer_tmp, buffsize*sizeof(int));
+               free(buffer_tmp);
+       }
+
+       static double compute_si(const unsigned char *ref, const unsigned char *cmp, const unsigned int& x, const unsigned int& y){
+               int sizebuff = x*y;
+               int *buff_sob1 = (int *)malloc(sizebuff*sizeof(int));
+               char c;
+
+               // create new buffer and apply both sobel filters
+               for (int i=0; i<sizebuff; i++)
+                       buff_sob1[i] = cmp[3*i];
+
+               conv_sobel(buff_sob1, x, y);
+
+               // get the standard deviation
+               double sdev = std_dev(buff_sob1, sizebuff);
+               return sdev;
+       }
+
+       static double compute_ti(const unsigned char *ref, const unsigned char *cmp, const unsigned int& x, const unsigned int& y){
+               int numelem = x*y;
+               int *buff_tmp = (int *)malloc(numelem*sizeof(int));
+
+               if(ref==NULL){
+                       return 0;
+               }
+
+               // difference between 2nd frame and 1st one to buff_tmp
+               for (int i=0; i<numelem; i++) {
+                       buff_tmp[i] = (int)cmp[3*i]-(int)ref[3*i];
+               }
+
+               // find and return the standard deviation of buff_tmp
+               double sdev = std_dev(buff_tmp, numelem);
+               free(buff_tmp);
+               return sdev;
+       }
+
        static inline double r_0_1(const double& d) {
                return std::max(0.0, std::min(1.0, d));
        }
@@ -224,6 +310,42 @@ namespace stats {
                }
        }
 
+       class siti_job : public mt::ThreadPool::Job {
+                       const unsigned char     *_ref,
+                                               *_cmp;
+                       const unsigned int      _x,
+                                               _y;
+                       std::pair<double,double>                        &_res;
+               public:
+                       siti_job(const unsigned char *ref, const unsigned char *cmp, const unsigned int& x, const unsigned int& y, std::pair<double,double>& res) :
+                       _ref(ref), _cmp(cmp), _x(x), _y(y), _res(res) {
+                       }
+
+                       virtual void run(void) {
+                               _res.first = compute_si(_ref, _cmp, _x, _y);
+                               _res.second= compute_ti(_ref, _cmp, _x, _y);
+                       }
+               };
+
+       static void get_siti_tp(const VUCHAR* ref, const std::vector<bool>& v_ok, const std::vector<VUCHAR>& streams, std::vector<std::pair<double,double> >& res, const unsigned int& x, const unsigned int& y) {
+                       const unsigned int                      sz = v_ok.size();
+                       std::vector<shared_ptr<siti_job> >      v_jobs;
+                       for(unsigned int i =0; i < sz; ++i) {
+                               if (v_ok[i]) {
+                                       v_jobs.push_back(new siti_job(&((*ref)[0]), &(streams[i][0]), x, y, res[i]));
+                                       __stats_tp.add(v_jobs.rbegin()->get());
+                               } else{
+                                       res[i].first = 0.0;
+                                       res[i].second=0.0;
+                               }
+                       }
+                       //wait for all
+                       for(std::vector<shared_ptr<siti_job> >::iterator it = v_jobs.begin(); it != v_jobs.end(); ++it) {
+                               (*it)->wait();
+                               (*it) = 0;
+                       }
+               }
+
        class hsi_job : public mt::ThreadPool::Job {
                unsigned char           *_buf;
                const unsigned int      _sz;
@@ -527,6 +649,38 @@ namespace stats {
                        }
                }
        };
+
+       class siti : public s_base {
+               protected:
+                       VUCHAR* _prev;
+                       void print(const int& ref_frame, const std::vector<std::pair<double,double> >& v_res) {
+                               _ostr << ref_frame << ',';
+                               for(int i = 0; i < _n_streams; ++i)
+                                       _ostr << v_res[i].first << ',' << v_res[i].second << ",";
+                               _ostr << std::endl;
+                       }
+               public:
+                       siti(const int& n_streams, const int& i_width, const int& i_height, std::ostream& ostr) :
+                       s_base(n_streams, i_width, i_height, ostr){
+                               _prev=NULL;
+                       }
+
+                       virtual void process(const int& ref_frame, VUCHAR& ref, const std::vector<bool>& v_ok, std::vector<VUCHAR>& streams) {
+                               if (v_ok.size() != streams.size() || v_ok.size() != _n_streams) throw std::runtime_error("Invalid data size passed to analyzer");
+                               // convert to Y colorspace
+                               rgb_2_Y_tp(ref, v_ok, streams);
+                               //
+                               if(_prev==NULL){
+                                       _prev=new std::vector<unsigned char> (streams[0]);
+                               }
+                               std::vector<std::pair<double,double> >  v_res(_n_streams);
+                               get_siti_tp(_prev, v_ok, streams, v_res, _i_width, _i_height);
+                               delete(_prev);
+                               _prev=new std::vector<unsigned char> (streams[0]);
+                               //
+                               print(ref_frame, v_res);
+                       }
+               };
 }
 
 stats::s_base* stats::get_analyzer(const char* id, const int& n_streams, const int& i_width, const int& i_height, std::ostream& ostr) {
@@ -535,5 +689,6 @@ stats::s_base* stats::get_analyzer(const char* id, const int& n_streams, const i
        else if (s_id == "avg_psnr") return new avg_psnr(n_streams, i_width, i_height, ostr);
        else if (s_id == "ssim") return new ssim(n_streams, i_width, i_height, ostr);
        else if (s_id == "avg_ssim") return new avg_ssim(n_streams, i_width, i_height, ostr);
+       else if (s_id == "siti") return new siti(n_streams, i_width, i_height, ostr);
        throw std::runtime_error("Invalid analyzer id");
 }