[GRASS-SVN] r68553 - sandbox/bo/i.segment.gsoc2016/i.segment
svn_grass at osgeo.org
svn_grass at osgeo.org
Tue May 31 21:39:01 PDT 2016
Author: hao2309
Date: 2016-05-31 21:39:01 -0700 (Tue, 31 May 2016)
New Revision: 68553
Added:
sandbox/bo/i.segment.gsoc2016/i.segment/mean_shift_2016_05_30_yb.c
Log:
basic mean shift codes 2016_05_30_yb
Added: sandbox/bo/i.segment.gsoc2016/i.segment/mean_shift_2016_05_30_yb.c
===================================================================
--- sandbox/bo/i.segment.gsoc2016/i.segment/mean_shift_2016_05_30_yb.c (rev 0)
+++ sandbox/bo/i.segment.gsoc2016/i.segment/mean_shift_2016_05_30_yb.c 2016-06-01 04:39:01 UTC (rev 68553)
@@ -0,0 +1,224 @@
+#include <stdio.h>
+#include <stdlib.h>
+#include <string.h>
+#include <math.h>
+#define FILEPATH "testing/test_map.ascii"
+
+/* =============input parameters=============== */
+//input spatial bandwidth = moving window radius
+#define RADIUS 2
+
+//input range bandwidth in normalized difference
+#define RANGEWIDTH 0.1
+
+//EPSILON to control when iterations should stop
+#define EPSILON 0.01
+
+//max number of iterations
+#define MAXITER 20
+
+//Threhold to cluster image
+#define CLUSTER_BANDWIDTH 0.5
+
+// row number:
+#define ROWNUM 8
+//column number:
+#define COLNUM 8
+
+
+
+void read_matrix(float **data_matrix){
+ //open the matrix file
+ FILE *data_file;
+ int i,j,m;
+ if ((data_file = fopen(FILEPATH, "r"))==NULL){
+ printf("ERROR: can not open file");
+ exit(0);
+ };
+ //loading the head information
+ char head[8][100];
+ for(i=0;i<8;i++){
+ fgets(head[i], 100, data_file);
+ printf (head[i]);
+ }
+ //loading the image matrix
+
+ for (i=0;i<ROWNUM;i++){
+ for (j=0;j<COLNUM;j++){
+ if(!fscanf(data_file,"%f",&data_matrix[i][j]))
+ break;
+ }
+ }
+ fclose(data_file);
+}
+
+//print out the matrix
+void print_matrix(float **mtr, int mtr_rows, int mtr_cols){
+ int i,j;
+ for (i=0;i<mtr_rows;i++)
+ {
+ for (j=0;j<mtr_cols;j++)
+ {
+ printf("\t%.2f", mtr[i][j]);
+ }
+ printf("\n");
+ }
+}
+
+//copy matrix
+void copy_matrix(float **newmtr, float **srcmtr)
+{
+ int i,j;
+ for (i=0;i<ROWNUM;i++){
+ for (j=0;j<COLNUM;j++){
+ newmtr[i][j] = srcmtr[i][j];
+ }
+ }
+}
+
+//cluster bandwidth simplified method to clustering image matrix
+void clump_matrix(float **mtr){
+ //the value of class
+ int class_value = 0;
+ //how many pixels clumped
+ int px_count = 0;
+ //remember the first stop location
+ int first_stop_row = 0;
+ int first_stop_col = 0;
+ int pixel_indicator_matrix[ROWNUM][COLNUM] = {{0}};//the indicator matrix to mark if the pixel has been clumped
+ //start iterations
+ while(px_count<ROWNUM*COLNUM){
+ class_value+=1;// change class value
+ int first_change = 0;//mark if the pixel is the first non-clump pixel
+ //go through all rows
+ float first_stop_value = mtr[first_stop_row][first_stop_col];//based pixel to be compared with other pixels
+ for (int i=0;i<ROWNUM;i++){
+ //go through all cols
+ for (int j=0;j<COLNUM;j++){
+ //condition: use cluster-bandwidth threhold to cluster the image
+ if (fabs(mtr[i][j]-first_stop_value) <= CLUSTER_BANDWIDTH && pixel_indicator_matrix[i][j]==0){
+ pixel_indicator_matrix[i][j] = 1;//mark pixel in indicator matrix 0 as non-clump, 1 as clump
+ mtr[i][j] = class_value;// assign the cluster value
+ px_count+=1;
+ }
+ //mark the first non-clumped pixel
+ else if(first_change == 0 && pixel_indicator_matrix[i][j]==0){
+ first_change =1;// mark the first time non-clump
+ first_stop_row = i;//record the pixel location of first time non-clump
+ first_stop_col = j;
+ }
+ }//end of cols
+ }//end of rows
+ }//end iteration
+}
+
+int main(){
+ float **in_matrix;//initialte input matrix
+ float **temp_matrix;//initialte last iteration matrix
+ float **out_matrix;//initialte output matrix
+ //window size in pixel number
+ int window_size = RADIUS*2+1;
+ // iter counter
+ int iter = 0;
+ // number of changed pixel > epsilon
+ /* ================assign mem to matrix==================== */
+ int n_changes = 1;
+ in_matrix = malloc(sizeof(*in_matrix) *ROWNUM);
+ temp_matrix= malloc(sizeof(*in_matrix) *ROWNUM);
+ out_matrix= malloc(sizeof(*in_matrix) *ROWNUM);
+ for (int i=0;i<8;i++){
+ in_matrix[i] = malloc(sizeof(*in_matrix[i]) * COLNUM);
+ temp_matrix[i]= malloc(sizeof(*in_matrix) *ROWNUM);
+ out_matrix[i]= malloc(sizeof(*in_matrix) *ROWNUM);
+ }
+ //end assign mem
+ printf ("The head file info of the image is: \n");
+ read_matrix(in_matrix);//read head info and data
+ copy_matrix(out_matrix, in_matrix);//copy input raster to output raster
+ while (iter<MAXITER && n_changes>0){
+ //initialte the in raster, outraster and last iteration raster.
+ iter += 1;
+ //make output raster of last iteration input raster of this iteration
+ copy_matrix(temp_matrix, in_matrix);
+ copy_matrix(in_matrix, out_matrix);
+ copy_matrix(out_matrix, temp_matrix);
+ //go through all rows
+ for (int current_row=0;current_row<ROWNUM;current_row++){
+ //# go through all columns
+ for (int current_col=0;current_col<COLNUM;current_col++){
+ //create moving window, resize in the image edge
+ int window_row_west = current_row - RADIUS;//left bound of the window
+ int window_row_east = window_row_west + window_size;//right bound of the window
+ //clip window in the edge
+ if (window_row_west<0){
+ window_row_west = 0;
+ }
+ if (window_row_east>ROWNUM){
+ window_row_east = ROWNUM;
+ }
+ int window_col_north = current_col - RADIUS;//upper bound of the window
+ int window_col_south = current_col + window_size;//lower bound of the window
+ if (window_col_north<0){
+ window_col_north = 0;
+ }
+ if (window_col_south>ROWNUM){
+ window_col_south = ROWNUM;
+ }
+ float sum_of_weights = 0;
+ float new_raster_value = 0;
+ //start of the moving window
+ for (int current_window_row = window_row_west;current_window_row<window_row_east;current_window_row++){
+ for (int current_window_col = window_col_north;current_window_col<window_col_south;current_window_col++){
+ //check spatial distance (bandwidth)
+ float weight = 1;
+ if (sqrt(pow((current_window_row - current_row),2)+pow((current_window_col - current_col),2))<=RADIUS){
+ //check range bandwidth
+ if (fabs((in_matrix[current_window_row][current_window_col]-in_matrix[current_row][current_col])/in_matrix[current_row][current_col])<=RANGEWIDTH){
+ weight *=1;
+ sum_of_weights += weight;
+ new_raster_value += in_matrix[current_window_row][current_window_col]*weight;
+ }
+ }
+ }
+ }//end of moving window
+ new_raster_value /= sum_of_weights;
+ out_matrix[current_row][current_col] = new_raster_value;
+ //condition: if the pixel converged(less than epsilon)
+ if (fabs(new_raster_value - temp_matrix[current_row][current_col])>EPSILON){
+ n_changes++;
+ }
+
+ }//end of all all columns
+ }//end of all rows
+ }//end of iteration
+ printf ("The image after mean filtering is: \n");
+ print_matrix(out_matrix,ROWNUM,COLNUM);
+/* =================clustering image ============================ */
+ printf ("The image after clustering is: \n");
+ clump_matrix(out_matrix);
+ print_matrix(out_matrix,ROWNUM,COLNUM);
+ //===================free matrix=============================
+ for (int i=0;i<8;i++){
+ free(in_matrix[i]);
+ free(temp_matrix[i]);
+ free(out_matrix[i]);
+ }
+ free(in_matrix);
+ free(temp_matrix);
+ free(out_matrix);
+ printf ("Program finished running\n");
+ return 0;
+}
+
+
+
+
+
+
+
+
+
+
+
+
+
More information about the grass-commit
mailing list