graft_cv/utils.h

#pragma once

#include <opencv2\imgproc\imgproc.hpp>
#include <vector>
#include <time.h>
#include <string>
#include <numeric>

#include "data_def_api.h"
#include "data_def.h"

using namespace cv;
using namespace std;

namespace graft_cv{

	// convert ImgInfo to cv mat 
	Mat imginfo2mat(ImgInfo*);
	ImgInfo* mat2imginfo(const Mat&);
	void imginfo_release(ImgInfo**);
	bool isvalid(const ImgInfo*);
	string currTime();
	string getImgId(int im_type);
	string getImgIdOa(string iid, int im_idx);
	inline void throw_msg(string& msg);

	// set bottom image rows_cnt rows to "value", for xiaoxu camera
	void image_set_bottom(Mat& img,unsigned char value,int rows_cnt);
	void image_set_top(	Mat& img,	unsigned char value,	int rows_cnt);
	void image_draw_line(Mat& img,int x0,int y0,int x1,int y1);
	// histogram
	void mat_histogram(
		Mat&, 
		std::vector<int>&hist,
		int axis=0, 
		int r0=-1, 
		int r1=-1, 
		int c0=-1, 
		int c1=-1
		);
	// weighted histogram
	void mat_histogram_w(
		Mat& img,	
		std::vector<int>&hist,
		int axis=0,
		int r0=-1, 
		int r1=-1, 
		int c0=-1, 
		int c1=-1,
		bool asc_w=true //true---ascending weight with x/y coordinate, [0-1.0], false --- descending weight [1.0--0.0]
	);
	void mat_histogram_yfork(
		Mat& img,	
		std::vector<int>&hist,	
		int r0, 
		int r1	
	);

	// histogram to image	
	void hist2image(
		vector<int>& hist,
		Mat&img,
		int aix=1,
		int grid_col=50,
		int grid_row=50
		);
	
	template<class T>
	void hist2image_line(
	vector<T>& hist, 
	Mat&img, 	
	int grid_col/*=50*/,
	int grid_row/*=50*/
	)
{
	if(hist.size()<2){return;}
	vector<T>::iterator maxit = max_element(begin(hist),end(hist)); 
	if( grid_col<10){grid_col=10;}
	if( grid_row<10){grid_row=10;}	
	
	img = Mat::zeros(*maxit,hist.size(), CV_8UC1);
	drawgrid(img,grid_col, grid_row,1);
	for(size_t i=1; i<hist.size();++i){
		line(img,Point(i-1,(int)*maxit - (int)hist[i-1]),Point(i,(int)*maxit - (int)hist[i]),Scalar(255));
	}
	
};

	template<class T>
	void hist2image_scale_line(
	vector<T>& hist, 
	Mat&img, 	
	T ymin,
	T ymax,
	T yratio,//数值像素比，一个像素代表的数值量
	T th,
	int grid_col/*=50*/,
	int grid_row/*=50*/
	)
{
	if(hist.size()<2){return;}	
	if( grid_col<10){grid_col=10;}
	if( grid_row<10){grid_row=10;}	
	int height =(int)((ymax-ymin)/yratio);
	
	img = Mat::zeros(height,hist.size(), CV_8UC1);
	drawgrid(img,grid_col, grid_row,1);
	for(size_t i=1; i<hist.size();++i){
		int y_1 = height-(int)((hist[i-1]-ymin)/yratio);
		int y_0 = height-(int)((hist[i]-ymin)/yratio);
		line(img,Point(i-1,y_1),Point(i,y_0),Scalar(255));
	}
	int h1 = height-(int)((th-ymin)/yratio);
	int h2 = height-(int)((1.0/th-ymin)/yratio);
	line(img,Point(0,h1),Point(hist.size()-1,h1),Scalar(200));
	line(img,Point(0,h2),Point(hist.size()-1,h2),Scalar(200));
	
};

	void drawgrid(		
		Mat&img,		
		int grid_col=50,
		int grid_row=50,
		int line_thick=1
		);

	//r2 index calculation for line regression
	//         y = beta0 + beta1 * x
	template<class T>
	double r2index(
		const vector<T>& x, 
		const vector<T>& y, 
		size_t i0, 
		size_t i1,
		double& beta0,
		double& beta1)
	{
		double r2 = 0.0;
		assert(x.size() == y.size());
		assert(i0<i1);
		assert(i1-i0>1);
		assert(i0>=0 && i1<y.size());
    
		// y = beta0 + beta1 * x
		//double beta0=0, beta1=0;
		double t1=0, t2=0,t3=0,t4=0;
		double n = (double)(i1-i0+1);
		for(size_t i=i0; i<=i1; ++i)
		{
			t1 += (double)x[i] * (double)x[i];
			t2 += (double)x[i];
			t3 += (double)x[i] * (double)y[i];
			t4 += (double)y[i];
		}
		beta0 = (t1*t4 - t2*t3) / (t1*n - t2*t2);
		beta1 = (t3*n - t2*t4) / (t1*n -t2*t2);
    
		double sst=0, ssr=0, sse=0;
		double y_mu = t4/n;
		double y_hat;
		for(size_t i=i0; i<=i1; ++i)
		{
			y_hat = beta0 + beta1*(double)x[i];
			ssr += (y_hat - y_mu)*(y_hat - y_mu);
			sse += ((double)y[i] - y_hat)* ((double)y[i] - y_hat);
		}
		sst = ssr + sse;
		r2 = 1.0 - sse/sst;   
    
		return r2;
	};

	void get_stem_x_range_oa(
		const std::vector<int>& hist_h, 
		double h_ratio, 
		int padding, 
		int cent_x,
		int width_x,
		int&x0, 
		int&x1, 
		int&stem_x0,
		int&stem_x1
		);
	void get_stem_x_range_rscut(
		const std::vector<int>& hist_h, 
		double h_ratio, 
		int padding, 
		int cent_x,
		int width_x,
		int&x0, 
		int&x1, 
		int&stem_x0,
		int&stem_x1
		);
	
	void get_stem_x_range_scion(
		const std::vector<int>& hist_h, 
		double h_ratio, 
		int padding, 
		int&x0, 
		int&x1, 
		int&stem_x0,
		int&stem_x1
		);

	////////////////////////////////////////////////////////////////////////
	// y-fork detect
	void get_stem_y_fork(
		const std::vector<int>& hist, 
		double ratio,
		int stem_dia_min, 
		int stem_fork_y_min,
		double stem_dia_mp,
		int& fork_y,     //output, 茎分叉y像素位置
		int& end_y,      //output，茎根部y像素位置
		int& stem_dia    //output， end_y和fork_y间茎粗
	);

	void get_stem_y_fork_rs(
	const std::vector<int>& hist, 
	double ratio,
	int stem_dia_min, 
	int stem_fork_y_min,
	double stem_dia_mp,
	int& fork_y, //output
	int& end_y,  //output
	int& stem_dia,    //output
	int& roi_max_y   //output
	);
	void get_stem_y_fork_rs_update(
	const Mat& bin_img,
	int stem_x_padding,
	int x0,
	int x1,
	int max_y,
	int stem_root_y,
	int stem_dia,
	bool image_show,
	double cut_point_offset_ratio,
	Point& fork_cent,
	double& max_radius,
	double& stem_angle
	);
	double yfork_validity_position(
		int max_y, 
		int end_y, 
		int fork_y);
	double yfork_validity_stemdiaratio(
		int stem_dia, 
		int dia_y,
		double opt_dia_ratio);
	void yfork_validity_stemdiachange(
	    const std::vector<int>& hist,
		int stem_dia_min,
		std::vector<double>& var_ratio		
		);

	//通过计算茎均值，运用统计过程控制方法检测显著点
	void get_stem_y_fork_3sigma(
		const std::vector<int>& hist, 
		double ratio,
		int stem_dia_min, 
		int stem_fork_y_min,
		double stem_dia_mp,
		int& fork_y,     //output, 茎分叉y像素位置
		int& end_y,      //output，茎根部y像素位置
		int& stem_dia    //output， end_y和fork_y间茎粗
	);

	double stem_y_fork_validity(
		const std::vector<int>& hist,
		int fork_y,
		int end_y);

	void get_hist_segment(
		const std::vector<int>& hist, 
		int threshold,
		std::vector<gcv_point<int>>& segments
	);
	void get_hist_mean_var_local(
		const std::vector<int>& hist, 
		int x0,
		int x1,
		double& mean,
		double& var
	);
	//void get_stem_y_index(
	//	Mat& bin_img,
	//	int x0,
	//	int x1,
	//	int stem_dia_min, 
	//	int stem_fork_y_min,
	//	double stem_dia_mp,
	//	int& fork_y, //output
	//	int& end_y,  //output
	//	int& stem_dia    //output
	//);
	void imshow_wait(
		const char* winname, 
		Mat&, 
		int waittype=-1
		);

	int get_stem_fork_right(
		Mat&stem_img, 
		int stem_fork_y, 
		int stem_x0, 
		int stem_x1,
		int detect_window);
	int get_stem_fork_left(
		Mat&stem_img, 
		int stem_fork_y, 
		int stem_x0, 
		int stem_x1,
		int detect_window);
	void get_stem_fork_xs(
		Mat&stem_img, 
		int stem_fork_y, 
		int stem_x0, 
		int stem_x1,
		int x0, 
		int x1,		
		int & fork_x0,
		int & fork_x1
	);

	void get_next_pt(
	    Mat& edge_img,
	    gcv_point<int>&start_pt, 
	    gcv_point<int>&pre_pt, 
		gcv_point<int>&nxt_pt, //output
	    bool is_up=true);

	// 二次曲线拟合，返回二次曲线对称轴x的坐标
	double qua_fitting(vector<double>&xx, vector<double>&yy);

	template <typename T>
    vector<size_t> sort_indexes_e(vector<T> &v, bool ascending=true)
    {
        vector<size_t> idx(v.size());
        iota(idx.begin(), idx.end(), 0);
		if(ascending){
        sort(idx.begin(), idx.end(),
            [&v](size_t i1, size_t i2) {return v[i1] < v[i2]; });
		}
		else{
			sort(idx.begin(), idx.end(),
            [&v](size_t i1, size_t i2) {return v[i1] > v[i2]; });		
		}
        return idx;
    }
	template<typename T>
	void hist_filter(vector<T>& hist, int method=0, int radius=2)
	{
		//mothod: 0---mid; 1--mean
		assert(radius>=1);
		if(hist.size()<3){return;}
		vector<T>tmp;
		
		if(method==0){
			T* buff = new T[2+2*radius];
			for(int i=0; i<hist.size();++i){
				int cnt=0;
				for(int r = -radius;r<=radius;++r){
					int idx=i+r;
					if(idx<0 || idx>=hist.size()){continue;}
					buff[cnt++]=hist[idx];
				}
				sort(buff,buff+cnt);
				int cent = cnt/2;
				tmp.push_back(buff[cent]);
			}
			delete []buff;
		}
		if(method==1){
			for(int i=0; i<hist.size();++i){
				int cnt=0;
				T mean = T(0);
				for(int r = -radius;r<=radius;++r){
					int idx=i+r;
					if(idx<0 || idx>=hist.size()){continue;}
					mean+=hist[idx];
				}
				mean =mean/cnt;
				tmp.push_back(mean);
			}
		}
		
		for(size_t j=0;j<hist.size();++j){hist[j]=tmp[j];}
	}

	


	
};