[GRASS-SVN] r67003 - grass-addons/grass7/imagery/i.spec.sam
svn_grass at osgeo.org
svn_grass at osgeo.org
Wed Dec 2 02:38:57 PST 2015
Author: ychemin
Date: 2015-12-02 02:38:57 -0800 (Wed, 02 Dec 2015)
New Revision: 67003
Added:
grass-addons/grass7/imagery/i.spec.sam/test1_spec_angle.py
Modified:
grass-addons/grass7/imagery/i.spec.sam/main.c
grass-addons/grass7/imagery/i.spec.sam/open.c
Log:
Added a manual test for specangle in Python
Modified: grass-addons/grass7/imagery/i.spec.sam/main.c
===================================================================
--- grass-addons/grass7/imagery/i.spec.sam/main.c 2015-12-02 09:49:26 UTC (rev 67002)
+++ grass-addons/grass7/imagery/i.spec.sam/main.c 2015-12-02 10:38:57 UTC (rev 67003)
@@ -57,7 +57,7 @@
char result_name[80];
char *result_prefix;
-mat_struct *open_files(char * matrixfile, char *img_grp);
+mat_struct *open_files(char * matrixfile, char *img_grp, int type);
float spectral_angle(vec_struct * Avector1, vec_struct * Avector2, int vtype);
DCELL myround(double x);
@@ -75,6 +75,7 @@
} parm;
mat_struct *A; /*first use in open.c G_matrix_set()*/
+ mat_struct *AT; /*first use in open.c G_matrix_set()*/
char *group;
float spectangle; /*numerical value of the spectral angle*/
@@ -88,7 +89,7 @@
parm.group = G_define_standard_option(G_OPT_I_GROUP);
- parm.group->description = "Imagery group to target for Spectral Mixture Analyis";
+ parm.group->description = "Imagery group to target for Spectral Mixture Analysis";
parm.matrixfile = G_define_standard_option(G_OPT_F_INPUT);
parm.matrixfile->description = "Matrix file containing spectral signatures";
@@ -107,7 +108,7 @@
G_message("%s",result_prefix);
/*Creating A, the spectral signature matrix here*/
- A = open_files(parm.matrixfile->answer, parm.group->answer);
+ A = open_files(parm.matrixfile->answer, parm.group->answer, 1);
/* Spectral Matrix is stored in A now */
/* Check matrix orthogonality
@@ -115,17 +116,21 @@
* using spectral unmixing and modeling spectral mixtrues with
* TM data. Photogrammetric Engineering & Remote Sensing, Vol.63, No6.
*/
-
+ AT = open_files(parm.matrixfile->answer, parm.group->answer, 0);
G_message("/* Check matrix orthogonality*/");
for (i = 0; i < Ref.nfiles; i++) /* Ref.nfiles = matrixsize*/
{
- Avector = G_matvect_get_row(A, i); /* go columnwise through matrix*/
+ Avector = G_matvect_get_column(AT, i); /* go columnwise through matrix*/
+ G_message("Avector rows:%d cols:%d",Avector->rows,Avector->cols);
for (j = 0; j < Ref.nfiles ; j++)
{
if (j !=i)
{
- b = G_matvect_get_row(A, j); /* compare with next col in A */
+ b = G_matvect_get_column(AT, j); /* compare with next col in A */
+ G_message("b rows:%d cols:%d",b->rows,b->cols);
+ G_message("process spectangle %d %d %d",i,j, Ref.nfiles);
spectangle = spectral_angle(Avector, b, RVEC);
+ G_message("processed spectangle");
anglefield[i][j]= spectangle;
G_vector_free(b);
}
Modified: grass-addons/grass7/imagery/i.spec.sam/open.c
===================================================================
--- grass-addons/grass7/imagery/i.spec.sam/open.c 2015-12-02 09:49:26 UTC (rev 67002)
+++ grass-addons/grass7/imagery/i.spec.sam/open.c 2015-12-02 10:38:57 UTC (rev 67003)
@@ -13,7 +13,7 @@
int G_matrix_read2(FILE * fp, mat_struct * out); /* Modified version of G_matrix_read(..). */
-mat_struct *open_files(char *matrixfile, char *img_grp)
+mat_struct *open_files(char *matrixfile, char *img_grp, int type)
{
char *name, *mapset;
FILE *fp;
@@ -35,11 +35,14 @@
A = G_matrix_init(A_input.cols, A_input.rows, A_input.cols);/*changed r/c*/
if (A == NULL)
G_fatal_error(_("Unable to allocate memory for matrix"));
-
- A = G_matrix_transpose(&A_input);/*transposed for spec angle process*/
-
- if(A->rows < A->cols)
- G_fatal_error("Need m (rows) >= n (cols) to obtain least squares fit\n");
+ if(type==1)
+ A = G_matrix_transpose(&A_input);/*transposed for spec angle process*/
+ if(type==0){
+ A = G_matrix_copy(&A_input);/*for matrix orthogonality check*/
+ return A;
+ }
+/* if(A->rows < A->cols)
+ G_fatal_error("Need m (rows) >= n (cols) to obtain least squares fit\n");*/
/*Only for debug, so temporary disabled*/
/*
G_verbose_message("Your spectral matrix = ");
@@ -48,7 +51,7 @@
m_output(A);
}
*/
- matrixsize=A->rows; /*changed cols to rows*/
+ matrixsize=A->rows;
G_message("/* open input files from group */");
/* open input files from group */
Added: grass-addons/grass7/imagery/i.spec.sam/test1_spec_angle.py
===================================================================
--- grass-addons/grass7/imagery/i.spec.sam/test1_spec_angle.py (rev 0)
+++ grass-addons/grass7/imagery/i.spec.sam/test1_spec_angle.py 2015-12-02 10:38:57 UTC (rev 67003)
@@ -0,0 +1,65 @@
+#!/usr/bin/env
+import numpy as np
+
+#Spectral data to load as maps
+# Band: r g b i1 i2 i3
+# Spektren zeilenweise eingeben!
+# 1. Sagebrush
+# 2. Saltbush
+# 3. Ground
+# 4. Dry Grass
+#row0: 8.87 13.14 11.71 35.85
+#Matrix: 4 by 4
+#row0: 8.87 13.14 11.71 35.85
+#row1: 13.59 20.12 19.61 50.66
+#row2: 28.26 34.82 38.27 40.1
+#row3: 10.54 16.35 23.7 38.98
+
+#Define the spectral signatures for each land use class
+cls1=[8.87, 13.14, 11.71, 35.85]
+cls2=[13.59, 20.12, 19.61, 50.66]
+cls3=[28.26, 34.82, 38.27, 40.10]
+cls4=[10.54, 16.35, 23.70, 38.98]
+
+#Generate a disturbed class manually
+cls1_wannabe=[9, 14, 12, 36]
+cls2_wannabe=[14, 21, 20, 51]
+cls3_wannabe=[29, 35, 39, 41]
+cls4_wannabe=[11, 17, 24, 39]
+
+
+#Define Norm of a Vector product
+def prod(v1, v2):
+ return np.sum(np.multiply(v1,v2))
+
+#Define Euclidian Norm of a vector ("L2-norm")
+def l2(v):
+ return np.sqrt((np.multiply(v,v)).sum(axis=0))
+
+
+#Define Spectral Angle of a vector pair
+def specangle(v1,v2):
+ return np.arccos(prod(v1,v2)/(l2(v1)*l2(v2)))*180/np.pi
+
+
+#Run for all 4 classes
+print("cls1")
+print(cls1)
+print(cls1_wannabe)
+print (specangle(cls1,cls1_wannabe))
+
+print("cls2")
+print(cls2)
+print(cls2_wannabe)
+print (specangle(cls2,cls2_wannabe))
+
+print("cls3")
+print(cls3)
+print(cls3_wannabe)
+print (specangle(cls3,cls3_wannabe))
+
+print("cls4")
+print(cls4)
+print(cls4_wannabe)
+print (specangle(cls4,cls4_wannabe))
+
More information about the grass-commit
mailing list