[GRASS-SVN] r57180 - in sandbox/turek/scatter_plot: . include include/defs lib lib/imagery

svn_grass at osgeo.org svn_grass at osgeo.org
Tue Jul 16 07:51:38 PDT 2013


Author: turek
Date: 2013-07-16 07:51:37 -0700 (Tue, 16 Jul 2013)
New Revision: 57180

Added:
   sandbox/turek/scatter_plot/include/
   sandbox/turek/scatter_plot/include/defs/
   sandbox/turek/scatter_plot/include/defs/imagery.h
   sandbox/turek/scatter_plot/include/imagery.h
   sandbox/turek/scatter_plot/lib/
   sandbox/turek/scatter_plot/lib/imagery/
   sandbox/turek/scatter_plot/lib/imagery/scatt_plot.c
Log:
scatter plot backend added

Added: sandbox/turek/scatter_plot/include/defs/imagery.h
===================================================================
--- sandbox/turek/scatter_plot/include/defs/imagery.h	                        (rev 0)
+++ sandbox/turek/scatter_plot/include/defs/imagery.h	2013-07-16 14:51:37 UTC (rev 57180)
@@ -0,0 +1,174 @@
+#ifndef GRASS_IMAGEDEFS_H
+#define GRASS_IMAGEDEFS_H
+
+/* alloc.c */
+void *I_malloc(size_t);
+void *I_realloc(void *, size_t);
+int I_free(void *);
+double **I_alloc_double2(int, int);
+int *I_alloc_int(int);
+int **I_alloc_int2(int, int);
+int I_free_int2(int **);
+int I_free_double2(double **);
+double ***I_alloc_double3(int, int, int);
+int I_free_double3(double ***);
+
+/* eol.c */
+int I_get_to_eol(char *, int, FILE *);
+
+/* find.c */
+int I_find_group(const char *);
+int I_find_group_file(const char *, const char *);
+int I_find_subgroup(const char *, const char *);
+int I_find_subgroup_file(const char *, const char *, const char *);
+
+/* fopen.c */
+FILE *I_fopen_group_file_new(const char *, const char *);
+FILE *I_fopen_group_file_append(const char *, const char *);
+FILE *I_fopen_group_file_old(const char *, const char *);
+FILE *I_fopen_subgroup_file_new(const char *, const char *, const char *);
+FILE *I_fopen_subgroup_file_append(const char *, const char *, const char *);
+FILE *I_fopen_subgroup_file_old(const char *, const char *, const char *);
+
+/* georef.c */
+int I_compute_georef_equations(struct Control_Points *, double *, double *,
+			       double *, double *, int);
+int I_georef(double, double, double *, double *, double *, double *, int);
+
+/* georef_tps.c */
+int I_compute_georef_equations_tps(struct Control_Points *, double **, double **,
+			       double **, double **);
+int I_georef_tps(double, double, double *, double *, double *, double *, 
+                 struct Control_Points *, int);
+
+/* group.c */
+int I_get_group(char *);
+int I_put_group(const char *);
+int I_get_subgroup(const char *, char *);
+int I_put_subgroup(const char *, const char *);
+int I_get_group_ref(const char *, struct Ref *);
+int I_get_subgroup_ref(const char *, const char *, struct Ref *);
+int I_init_ref_color_nums(struct Ref *);
+int I_put_group_ref(const char *, const struct Ref *);
+int I_put_subgroup_ref(const char *, const char *, const struct Ref *);
+int I_add_file_to_group_ref(const char *, const char *, struct Ref *);
+int I_transfer_group_ref_file(const struct Ref *, int, struct Ref *);
+int I_init_group_ref(struct Ref *);
+int I_free_group_ref(struct Ref *);
+
+/* iclass.c */
+struct Map_info;
+int I_iclass_analysis(IClass_statistics *, struct Ref *, struct Map_info *, const char *, const char *, const char *);
+int I_iclass_init_group(const char *, const char *, struct Ref *);
+void I_iclass_create_raster(IClass_statistics *, struct Ref *, const char *);
+
+/* iclass_statistics.c */
+void I_iclass_statistics_get_nbands(IClass_statistics *, int *);
+void I_iclass_statistics_get_cat(IClass_statistics *, int *);
+void I_iclass_statistics_get_name(IClass_statistics *, const char **);
+void I_iclass_statistics_get_color(IClass_statistics *, const char **);
+void I_iclass_statistics_get_ncells(IClass_statistics *, int *);
+int I_iclass_statistics_get_max(IClass_statistics *, int, int *);
+int I_iclass_statistics_get_range_max(IClass_statistics *, int, int *);
+int I_iclass_statistics_get_min(IClass_statistics *, int, int *);
+int I_iclass_statistics_get_range_min(IClass_statistics *, int, int *);
+int I_iclass_statistics_get_sum(IClass_statistics *, int, float *);
+int I_iclass_statistics_get_mean(IClass_statistics *, int, float *);
+int I_iclass_statistics_get_stddev(IClass_statistics *, int, float *);
+void I_iclass_statistics_get_nstd(IClass_statistics *, float *);
+void I_iclass_statistics_set_nstd(IClass_statistics *, float);
+int I_iclass_statistics_get_histo(IClass_statistics *, int, int, int *);
+int I_iclass_statistics_get_product(IClass_statistics *, int, int, float *);
+void I_iclass_init_statistics(IClass_statistics *, int, const char *, const char *, float);
+void I_iclass_free_statistics(IClass_statistics *);
+
+/* iclass_signatures.c */
+int I_iclass_init_signatures(struct Signature *, struct Ref *);
+void I_iclass_add_signature(struct Signature *, IClass_statistics *);
+int I_iclass_write_signatures(struct Signature *, const char *, const char *, const char *);
+
+/* list_gp.c */
+int I_list_group(const char *, const struct Ref *, FILE *);
+int I_list_group_simple(const struct Ref *, FILE *);
+
+/* list_subgp.c */
+int I_list_subgroup(const char *, const char *, const struct Ref *, FILE *);
+int I_list_subgroup_simple(const struct Ref *, FILE *);
+
+/* loc_info.c */
+char *I_location_info(const char *);
+
+/* points.c */
+int I_new_control_point(struct Control_Points *, double, double, double,
+			double, int);
+int I_get_control_points(const char *, struct Control_Points *);
+int I_put_control_points(const char *, const struct Control_Points *);
+
+/* ref.c */
+FILE *I_fopen_group_ref_new(const char *);
+FILE *I_fopen_group_ref_old(const char *);
+FILE *I_fopen_subgroup_ref_new(const char *, const char *);
+FILE *I_fopen_subgroup_ref_old(const char *, const char *);
+
+/* scatt_plt.c */
+void I_sc_init_cats(struct scCats *, int, int);
+void I_sc_free_cats(struct scCats *);
+void I_sc_get_active_categories(int *, int *, struct scCats *);
+int I_sc_add_cat(struct scCats *, int);
+int I_sc_delete_cat(struct scCats *, int);
+int I_sc_insert_scatt_data(struct scCats *, struct scdScattData *, int, int);
+int I_sc_remove_scatt_data(struct scCats *, struct scdScattData *, int, int);
+int I_sc_set_value(struct scCats *, int, int, int, int);
+
+void I_scd_init_scatt_data(struct scdScattData *, int, int);
+void I_scd_free_scatt_data(struct scdScattData *);
+int I_scd_get_data_size(struct scdScattData *);
+void * I_scd_get_data_ptr(struct scdScattData *);
+int I_scd_set_value(struct scdScattData *, unsigned int, unsigned int);
+
+int I_compute_scatts(struct Cell_head *, struct scCats *, const char **, 
+                     int, struct scCats *, const char **);
+
+/* sig.c */
+int I_init_signatures(struct Signature *, int);
+int I_new_signature(struct Signature *);
+int I_free_signatures(struct Signature *);
+int I_read_one_signature(FILE *, struct Signature *);
+int I_read_signatures(FILE *, struct Signature *);
+int I_write_signatures(FILE *, struct Signature *);
+
+/* sigfile.c */
+FILE *I_fopen_signature_file_new(const char *, const char *, const char *);
+FILE *I_fopen_signature_file_old(const char *, const char *, const char *);
+
+/* sigset.c */
+int I_SigSetNClasses(struct SigSet *);
+struct ClassData *I_AllocClassData(struct SigSet *, struct ClassSig *, int);
+int I_InitSigSet(struct SigSet *);
+int I_SigSetNBands(struct SigSet *, int);
+struct ClassSig *I_NewClassSig(struct SigSet *);
+struct SubSig *I_NewSubSig(struct SigSet *, struct ClassSig *);
+int I_ReadSigSet(FILE *, struct SigSet *);
+int I_SetSigTitle(struct SigSet *, const char *);
+const char *I_GetSigTitle(const struct SigSet *);
+int I_SetClassTitle(struct ClassSig *, const char *);
+const char *I_GetClassTitle(const struct ClassSig *);
+int I_WriteSigSet(FILE *, const struct SigSet *);
+
+/* sigsetfile.c */
+FILE *I_fopen_sigset_file_new(const char *, const char *, const char *);
+FILE *I_fopen_sigset_file_old(const char *, const char *, const char *);
+
+/* target.c */
+int I_get_target(const char *, char *, char *);
+int I_put_target(const char *, const char *, const char *);
+
+/* title.c */
+int I_get_group_title(const char *, char *, int);
+int I_put_group_title(const char *, const char *);
+
+/* var.c */
+double I_variance(double, double, int);
+double I_stddev(double, double, int);
+
+#endif

Added: sandbox/turek/scatter_plot/include/imagery.h
===================================================================
--- sandbox/turek/scatter_plot/include/imagery.h	                        (rev 0)
+++ sandbox/turek/scatter_plot/include/imagery.h	2013-07-16 14:51:37 UTC (rev 57180)
@@ -0,0 +1,180 @@
+#ifndef GRASS_IMAGERY_H
+#define GRASS_IMAGERY_H
+
+#include <grass/gis.h>
+#include <grass/raster.h>
+
+/* File/directory name lengths */
+#define INAME_LEN GNAME_MAX	/* coupled to raster map name length */
+
+struct Ref_Color
+{
+    unsigned char *table;	/* color table for min-max values */
+    unsigned char *index;	/* data translation index */
+    unsigned char *buf;	/* data buffer for reading color file */
+    int fd;			/* for image i/o */
+    CELL min, max;		/* min,max CELL values */
+    int n;			/* index into Ref_Files */
+};
+
+struct Ref_Files
+{
+    char name[INAME_LEN];	/* length is not in sync with other definitions */
+    char mapset[INAME_LEN];
+};
+
+struct Ref
+{
+    int nfiles;
+    struct Ref_Files *file;
+    struct Ref_Color red, grn, blu;
+};
+
+struct Tape_Info
+{
+    char title[75];
+    char id[2][75];
+    char desc[5][75];
+};
+
+struct Control_Points
+{
+    int count;
+    double *e1;
+    double *n1;
+    double *e2;
+    double *n2;
+    int *status;
+};
+
+struct One_Sig
+{
+    char desc[100];
+    int npoints;
+    double *mean;		/* one mean for each band */
+    double **var;		/* covariance band-band   */
+    int status;		/* may be used to 'delete' a signature */
+    float r, g, b;		/* color */
+    int have_color;
+};
+
+struct Signature
+{
+    int nbands;
+    int nsigs;
+    char title[100];
+    struct One_Sig *sig;
+};
+
+struct SubSig
+{
+    double N;
+    double pi;
+    double *means;
+    double **R;
+    double **Rinv;
+    double cnst;
+    int used;
+};
+
+struct ClassData
+{
+    int npixels;
+    int count;
+    double **x;		/* pixel list: x[npixels][nbands] */
+    double **p;		/* prob        p[npixels][subclasses] */
+};
+
+struct ClassSig
+{
+    long classnum;
+    char *title;
+    int used;
+    int type;
+    int nsubclasses;
+    struct SubSig *SubSig;
+    struct ClassData ClassData;
+};
+
+struct SigSet
+{
+    int nbands;
+    int nclasses;
+    char *title;
+    struct ClassSig *ClassSig;
+};
+
+/* IClass */
+
+/*! Holds statistical values for creating histograms and raster maps for one class.
+
+  One class is represented by one category (cat).
+*/
+typedef struct
+{
+    int cat;                /*!< class */
+    const char *name;       /*!< signature description (class name) */
+    const char *color;      /*!< class color (RRR:GGG:BBB)*/
+    int nbands;             /*!< number of bands */
+    
+    int ncells;             /*!< number of cells in training areas */
+
+    int *band_min;          /*!< minimum value for each band */
+    int *band_max;          /*!< maximum value for each band */
+    float *band_sum;        /*!< sum of values for each band */
+    float *band_mean;       /*!< mean of values for each band */
+    float *band_stddev;     /*!< standard deviation for each band */
+
+    float **band_product;   /*!< sum of products of cell category values of 2 bands */
+    int **band_histo;       /*!< number of cells for cell category value (0-256) for each band */
+
+    int *band_range_min;    /*!< min range of values to create raster map */
+    int *band_range_max;    /*!< max range of values to create raster map */
+    float nstd;             /*!< multiplier of standard deviation */
+
+    
+} IClass_statistics;
+
+/* Scatter Plot backend */
+struct scCats 
+{
+    int type; //TODO do not identify it with number
+
+    int n_cats;
+    
+    int n_bands;
+    int n_scatts;
+
+    int   n_a_cats;
+    int * cats_ids;
+    int * cats_idxs;
+
+    struct scScatts ** cats_arr;
+};
+
+struct scScatts
+{
+    int n_a_scatts;
+
+    int * scatts_bands;
+    int * scatt_idxs;
+
+    struct scdScattData ** scatts_arr;
+};
+
+struct scdScattData
+{
+    int n_vals;
+
+    unsigned int  * b_conds_arr;
+    unsigned int  * scatt_vals_arr;
+};
+
+#define SIGNATURE_TYPE_MIXED 1
+
+#define GROUPFILE "CURGROUP"
+#define SUBGROUPFILE "CURSUBGROUP"
+
+#include <grass/defs/imagery.h>
+
+#endif

Added: sandbox/turek/scatter_plot/lib/imagery/scatt_plot.c
===================================================================
--- sandbox/turek/scatter_plot/lib/imagery/scatt_plot.c	                        (rev 0)
+++ sandbox/turek/scatter_plot/lib/imagery/scatt_plot.c	2013-07-16 14:51:37 UTC (rev 57180)
@@ -0,0 +1,616 @@
+/*!
+   \file lib/imagery/scatt_plot.c
+
+   \brief Imagery library - functions for wx Scatter Plot Tool.
+
+   Low level functions used by wx Scatter Plot Tool.
+
+   Copyright (C) 2013 by the GRASS Development Team
+
+   This program is free software under the GNU General Public License
+   (>=v2).  Read the file COPYING that comes with GRASS for details.
+
+   \author Stepan Turek <stepan.turek at seznam.cz> (Mentor: Martin Landa)
+ */
+
+#include <grass/raster.h>
+#include <grass/imagery.h>
+#include <grass/gis.h>
+
+#include <stdio.h>
+#include <stdlib.h>
+#include <math.h>
+#include <string.h>
+#include <time.h>  
+
+
+//TODO split formulas to be more understandable
+static void id_scatt_to_bands(const int scatt_id, const int n_bands, int * band_1, int * band_2)
+{   
+    int n_b1 = n_bands - 1;
+
+    * band_1 = (int) ((2 * n_b1 + 1 - sqrt((double)((2 * n_b1 + 1) * (2 * n_b1 + 1) - 8 * scatt_id))) / 2);
+
+    * band_2 = scatt_id - ((* band_1) * (2 * n_b1 + 1) - (* band_1) * (* band_1)) / 2 + (* band_1) + 1;
+
+    return;
+}
+
+static void bands_to_id_scatt(const int band_1, const int band_2, const int n_bands, int * scatt_id)
+{   
+    int n_b1 = n_bands - 1;
+
+    * scatt_id = (band_1 * (2 * n_b1 + 1) - band_1 * band_1) / 2 + band_2 - band_1 - 1;
+
+    return;
+}
+
+void I_sc_init_cats(struct scCats * cats, int n_bands, int type)
+{
+    int i_cat;
+
+    cats->type = type;
+
+    cats->n_cats = 10;
+    cats->n_a_cats = 0;
+
+    cats->n_bands = n_bands;
+    cats->n_scatts = (n_bands - 1) * n_bands / 2;
+
+    cats->cats_arr = (struct scScatts **) malloc(cats->n_cats * sizeof(struct scScatts *));
+    memset(cats->cats_arr, 0, cats-> n_cats * sizeof(struct scScatts *));
+
+    cats->cats_ids = (int *)  malloc(cats->n_cats * sizeof(int));
+    cats->cats_idxs =(int *)  malloc(cats->n_cats * sizeof(int));
+
+    for(i_cat = 0; i_cat < cats->n_cats; i_cat++)
+        cats->cats_idxs[i_cat] = -1;
+
+    return;
+}
+
+void I_sc_free_cats(struct scCats * cats)
+{
+    int i_cat;
+
+    for(i_cat = 0; i_cat < cats->n_a_cats; i_cat++)
+    {        
+        if(cats->cats_arr[i_cat])
+        {   
+            free(cats->cats_arr[i_cat]->scatt_idxs);
+            free(cats->cats_arr[i_cat]->scatts_bands);
+            free(cats->cats_arr[i_cat]->scatts_arr);
+            free(cats->cats_arr[i_cat]);
+        }
+    }
+
+    free(cats->cats_ids);
+    free(cats->cats_idxs);
+    free(cats->cats_arr);
+
+    cats->n_cats = 0;
+    cats->n_a_cats = 0;
+    cats->n_bands = 0;
+    cats->n_scatts = 0;
+    cats->type = -1;
+
+    return;
+}
+
+void I_sc_get_active_categories(int * a_cats_ids, int * n_a_cats, struct scCats * cats)
+{
+    a_cats_ids = cats->cats_ids;
+    * n_a_cats = cats->n_a_cats;
+}
+
+int I_sc_add_cat(struct scCats * cats, int cat_id)
+{
+    int i_scatt;
+    int n_a_cats = cats->n_a_cats;
+
+    if(cat_id < 0 || cat_id >= cats->n_cats)
+        return -1;
+
+    if(cats->cats_idxs[cat_id] >= 0) 
+        return -1;
+    G_message("Adding category  %d", cat_id);
+
+    cats->cats_ids[n_a_cats] = cat_id;
+    cats->cats_idxs[cat_id] = n_a_cats;
+    
+    cats->cats_arr[n_a_cats] = (struct scScatts *) malloc(sizeof(struct scScatts));
+
+    cats->cats_arr[n_a_cats]->scatts_arr = (struct scdScattData **) malloc(cats->n_scatts * sizeof(struct scdScattData *));
+    memset((cats->cats_arr[n_a_cats]->scatts_arr), 0, cats->n_scatts * sizeof(struct scdScattData *));
+    
+    cats->cats_arr[n_a_cats]->n_a_scatts = 0;
+
+    cats->cats_arr[n_a_cats]->scatts_bands = (int *) malloc(cats->n_scatts * 2 * sizeof(int));
+     
+    cats->cats_arr[n_a_cats]->scatt_idxs = (int *) malloc(cats->n_scatts * sizeof(int));
+    for(i_scatt = 0; i_scatt < cats->n_scatts; i_scatt++)
+        cats->cats_arr[n_a_cats]->scatt_idxs[i_scatt] = -1;
+    
+    ++cats->n_a_cats;
+    G_message("nacats %d ", cats->n_a_cats);
+
+    return 0;
+}
+
+int I_sc_delete_cat(struct scCats * cats, int cat_id)
+{
+    int cat_idx, i_cat;
+
+    if(cat_id < 0 || cat_id >= cats->n_cats)
+        return -1;
+
+    cat_idx = cats->cats_idxs[cat_id];
+    if(cat_idx < 0)
+        return -1;
+
+    free(cats->cats_arr[cat_idx]->scatt_idxs);
+    free(cats->cats_arr[cat_idx]->scatts_bands);
+    free(cats->cats_arr[cat_idx]->scatts_arr);
+    free(cats->cats_arr[cat_idx]);
+
+    for(i_cat = cat_idx; i_cat < cats->n_a_cats - 1; i_cat++)
+    {
+        cats->cats_arr[i_cat] = cats->cats_arr[i_cat + 1];
+        cats->cats_ids[i_cat] = cats->cats_ids[i_cat + 1];
+    }
+    cats->cats_idxs[cat_id] = -1; 
+
+    --cats->n_a_cats;
+    
+    return 0;
+}
+    
+
+int I_sc_insert_scatt_data(struct scCats * cats, struct scdScattData * scatt_vals, int cat_id, int scatt_id)
+{
+    int band_1, band_2, cat_idx, n_a_scatts;
+    struct scScatts * scatts;
+
+    if(cat_id < 0 || cat_id >= cats->n_cats)
+        return -1;
+
+    cat_idx = cats->cats_idxs[cat_id];
+    if(cat_idx < 0)
+        return -1;
+
+    if(scatt_id < 0 && scatt_id >= cats->n_scatts)
+        return -1;
+
+    scatts = cats->cats_arr[cat_idx];
+    if(scatts->scatt_idxs[scatt_id] >= 0)
+        return -1;
+
+    if(!scatt_vals->b_conds_arr && cats->type == 1)
+        return -1;
+
+    if(!scatt_vals->scatt_vals_arr && cats->type == 0)
+        return -1;
+
+    n_a_scatts = scatts->n_a_scatts;
+
+    scatts->scatt_idxs[scatt_id] = n_a_scatts;
+
+    id_scatt_to_bands(scatt_id, cats->n_bands, &band_1, &band_2);
+
+    scatts->scatts_bands[n_a_scatts * 2] = band_1;
+    scatts->scatts_bands[n_a_scatts * 2 + 1] = band_2;
+
+    scatts->scatts_arr[n_a_scatts] = scatt_vals;
+    ++scatts->n_a_scatts;
+
+    G_message("I_sc_insert_scatt_data cat_id %d, scatt_id %d, scatt_idx", cat_id, scatt_id, n_a_scatts);
+
+    return 0;
+}
+
+int I_sc_remove_scatt_data(struct scCats * cats, struct scdScattData * scatt_vals, int cat_id, int scatt_id)
+{
+    int cat_idx, scatt_idx, n_init_scatts, i_scatt;
+    struct scScatts * scatts;
+
+    if(cat_id < 0 && cat_id >= cats->n_cats)
+        return -1;
+
+    cat_idx = cats->cats_idxs[cat_id];
+    if(cat_idx < 0)
+        return -1;
+
+    if(scatt_id < 0 || scatt_id >= cats->n_scatts)
+        return -1;
+
+    scatts = cats->cats_arr[cat_idx];
+    if(scatts->scatt_idxs[scatt_id] < 0)
+        return -1;
+
+    scatt_vals = scatts->scatts_arr[scatt_idx];
+
+    for(i_scatt = scatt_idx; i_scatt < scatts->n_a_scatts - 1; i_scatt++)
+    {
+        scatts->scatts_arr[i_scatt] = scatts->scatts_arr[i_scatt + 1];
+        scatts->scatts_bands[i_scatt * 2] = scatts->scatts_bands[(i_scatt + 1)* 2];
+        scatts->scatts_bands[i_scatt * 2 + 1] = scatts->scatts_bands[(i_scatt + 1) * 2 + 1];
+    }
+    scatts->scatts_arr[scatts->n_a_scatts] = NULL;
+
+    scatts->scatt_idxs[scatt_id] = -1;
+
+    scatt_vals = scatts->scatts_arr[scatt_id];
+    scatts->n_a_scatts--;
+
+    G_message("I_sc_remove_scatt_data %d, %d", cat_id, scatt_id);
+
+    return 0;
+}
+
+int I_sc_set_value(struct scCats * cats, int cat_id, int scatt_id, int value_idx, int value)
+{
+    int n_a_scatts = cats->cats_arr[cat_id]->n_a_scatts;
+    int cat_idx, scatt_idx, ret;
+
+    cat_idx = cats->cats_idxs[cat_id];
+    if(cat_idx < 0)
+        return -1;
+    
+    if(cats->cats_arr[cat_idx]->scatt_idxs[scatt_id] < 0)
+        return -1;
+
+    cat_idx = cats->cats_idxs[cat_id];
+    scatt_idx = cats->cats_arr[cat_idx]->scatt_idxs[scatt_id];
+
+    I_scd_set_value(cats->cats_arr[cat_idx]->scatts_arr[scatt_idx], value_idx, value);
+
+    return 0;
+}
+
+void I_scd_init_scatt_data(struct scdScattData * scatt_vals, int n_vals, int type)
+{
+    scatt_vals->n_vals = n_vals;
+
+    if(type == 0)
+    {   
+        scatt_vals->scatt_vals_arr = (unsigned int *) malloc(n_vals * sizeof(unsigned int));
+        memset(scatt_vals->scatt_vals_arr, 0, n_vals * sizeof(unsigned int));
+        scatt_vals->b_conds_arr = NULL;
+    }
+    else if(type == 1)
+    {   
+        scatt_vals->b_conds_arr = (unsigned short *) malloc(n_vals * sizeof(unsigned short));
+        memset(scatt_vals->b_conds_arr, 0, n_vals * sizeof(unsigned short));
+        scatt_vals->scatt_vals_arr = NULL;
+    }
+
+    return;
+}
+
+void I_scd_free_scatt_data(struct scdScattData * scatt_vals)
+{
+
+    free(scatt_vals->b_conds_arr);
+    free(scatt_vals->scatt_vals_arr);
+
+    scatt_vals = NULL;
+
+    return;
+}
+
+int I_scd_get_data_size(struct scdScattData * scatt_vals)
+{
+
+    return scatt_vals->n_vals;
+}
+
+void * I_scd_get_data_ptr(struct scdScattData * scatt_vals)
+{
+    if(!scatt_vals->b_conds_arr)
+        return scatt_vals->b_conds_arr;
+    else if(!scatt_vals->scatt_vals_arr)
+        return scatt_vals->scatt_vals_arr;
+
+    return NULL;
+}
+
+int I_scd_set_value(struct scdScattData * scatt_vals, unsigned int val_idx, unsigned int val)
+{
+    if(val_idx < 0 && val_idx >  scatt_vals->n_vals)
+        return -1;
+
+    if(scatt_vals->b_conds_arr)
+        scatt_vals->b_conds_arr[val_idx] = val;
+    else if(scatt_vals->scatt_vals_arr)
+        scatt_vals->scatt_vals_arr[val_idx] = val;
+    else
+        return -1;
+
+    return 0;
+}
+
+static inline void update_cat_scatt_plt(CELL ** chunks, int chunk_size, unsigned short * belongs_pix, struct scScatts * scatts)
+{
+    int band_axis_1, band_axis_2, i_scatt, array_idx, cat_idx, i_chunks_pix;
+    int r_bits = 256;
+
+    CELL * band_1_chunks;
+    CELL * band_2_chunks;
+
+    int * scatts_bands = scatts->scatts_bands;
+
+    for(i_scatt = 0; i_scatt < scatts->n_a_scatts; i_scatt++)
+    {   
+        band_1_chunks = chunks[scatts_bands[i_scatt * 2]];
+        band_2_chunks = chunks[scatts_bands[i_scatt * 2 + 1]];
+
+        for(i_chunks_pix = 0; i_chunks_pix < chunk_size; i_chunks_pix++)
+        {
+            if(!belongs_pix[i_chunks_pix])
+                continue;
+
+            array_idx = band_1_chunks[i_chunks_pix] + band_2_chunks[i_chunks_pix] * r_bits;
+            ++scatts->scatts_arr[i_scatt]->scatt_vals_arr[array_idx];
+        }
+    }
+}
+
+static inline int compute_scatts_from_chunk(struct scCats * scatt_plts, struct scCats * scatt_conds, 
+                                            CELL ** chunks, char ** null_chunks, int chunk_size, FILE ** f_cats_rasts, struct Cell_head *region, int i_chunk)
+{
+
+    int i_chunks_pix, i_cat, i_scatt, n_a_scatts, i_cond;
+    int cat_id, scatt_plts_cat_idx, array_idx;
+    char * band_1_null_chunks,* band_2_null_chunks;
+    
+    int r_bits = 256;
+
+    struct scScatts * scatts_conds;
+    struct scScatts * scatt_plts_scatts;
+    struct scdScattData * conds;
+
+    int * scatts_bands;
+    struct scdScattData ** scatts_arr;
+
+    CELL * band_1_chunks;
+    CELL * band_2_chunks;
+    unsigned int * i_scatt_conds;
+
+    unsigned short * belongs_pix = (unsigned short *) malloc(chunk_size * sizeof(unsigned short)); 
+    unsigned char * rast_pixs = (unsigned char *) malloc(chunk_size * sizeof(unsigned char));
+
+    for(i_cat = 0; i_cat < scatt_conds->n_a_cats; i_cat++)
+    {
+        scatts_conds = scatt_conds->cats_arr[i_cat];
+
+        cat_id = scatt_conds->cats_ids[i_cat];
+
+        scatt_plts_cat_idx = scatt_plts->cats_idxs[cat_id];
+        if(scatt_plts_cat_idx < 0)
+            continue;
+
+        scatts_bands = scatts_conds->scatts_bands;
+
+        scatt_plts_scatts = scatt_plts->cats_arr[scatt_plts_cat_idx];
+
+        // if no condition for cats_arr is defined, all pixels are taken
+        memset(belongs_pix, 0, chunk_size * sizeof(unsigned short));
+        if(!scatts_conds->n_a_scatts) {
+
+            for(i_scatt = 0; i_scatt < scatt_plts_scatts->n_a_scatts; i_scatt++)
+            {
+                band_1_null_chunks =  null_chunks[scatts_bands[i_scatt * 2]];
+                band_2_null_chunks =  null_chunks[scatts_bands[i_scatt * 2 + 1]];
+                    
+                for(i_chunks_pix = 0; i_chunks_pix < chunk_size; i_chunks_pix++){
+                    if(band_1_null_chunks[i_chunks_pix] != 1 && 
+                        band_2_null_chunks[i_chunks_pix] != 1)                   
+                        belongs_pix[i_chunks_pix] = 1;
+                }
+            }
+        }
+        else
+        {
+            // test every defined conditions in scatter plots
+            for(i_scatt = 0; i_scatt <  scatts_conds->n_a_scatts; i_scatt++)
+            {   
+                band_1_chunks =  chunks[scatts_bands[i_scatt * 2]];
+                band_2_chunks =  chunks[scatts_bands[i_scatt * 2 + 1]];
+
+                band_1_null_chunks =  null_chunks[scatts_bands[i_scatt * 2]];
+                band_2_null_chunks =  null_chunks[scatts_bands[i_scatt * 2 + 1]];
+
+                i_scatt_conds = scatts_conds->scatts_arr[i_scatt]->b_conds_arr;
+
+                for(i_chunks_pix = 0; i_chunks_pix < chunk_size; i_chunks_pix++)
+                {
+                    if(belongs_pix[i_chunks_pix] || 
+                       band_1_null_chunks[i_chunks_pix] == 1 || 
+                       band_2_null_chunks[i_chunks_pix] == 1)
+                        continue;
+                    if(i_scatt_conds[band_1_chunks[i_chunks_pix] + band_2_chunks[i_chunks_pix] * r_bits])
+                        belongs_pix[i_chunks_pix] = 1;
+                }                
+            }
+        }
+
+        for(i_chunks_pix = 0; i_chunks_pix < chunk_size; i_chunks_pix++)
+            rast_pixs[i_chunks_pix] = belongs_pix[i_chunks_pix]  & 255;
+
+        fwrite(rast_pixs, sizeof(unsigned char), (chunk_size)/sizeof(unsigned char), f_cats_rasts[i_cat]);
+        // TODO when code will be paralelized put it out of the loop
+        if (ferror(f_cats_rasts[i_cat]))
+        {
+            free(rast_pixs);
+            free(belongs_pix);
+            G_debug(3, "Unable to write into file.");
+            return -1;
+        }
+        update_cat_scatt_plt(chunks, chunk_size, belongs_pix, scatt_plts_scatts);
+    }
+
+    free(rast_pixs);
+    free(belongs_pix);
+
+    return 0;
+}
+
+static void get_needed_bands(struct scCats * cats, int * b_needed_bands)
+{
+    // results in b_needed_bands - array of bools - if item has value 1, band (defined by item index) is needed to be opened
+    int i_cat, i_scatt, cat_id;
+
+    for(i_cat = 0;  i_cat < cats->n_a_cats; i_cat++)
+    {   
+        for(i_scatt = 0;  i_scatt < cats->cats_arr[i_cat]->n_a_scatts; i_scatt++)
+        {   
+            G_debug(3, "Active scatt %d in catt %d", i_scatt, i_cat);
+
+            b_needed_bands[cats->cats_arr[i_cat]->scatts_bands[i_scatt * 2]] = 1;
+            b_needed_bands[cats->cats_arr[i_cat]->scatts_bands[i_scatt * 2 + 1]] = 1;
+        }
+    }
+    return;
+}
+
+static void free_compute_scatts_data(int * fd_bands, CELL ** chunks, char ** null_chunks, int * n_a_bands, int * bands_ids, 
+                               FILE ** f_cats_rasts, int n_a_cats)
+{
+    int i, band_id;
+
+    for(i = 0; i < (* n_a_bands); i++)
+    {   
+        band_id = bands_ids[i];
+        if(band_id >= 0) {
+            Rast_close(fd_bands[i]);
+            G_free(chunks[band_id]);
+            G_free(null_chunks[band_id]);
+        }
+    }
+
+    if(f_cats_rasts)
+        for(i = 0; i < n_a_cats; i++)
+            fclose(f_cats_rasts[i]);
+}
+
+
+//TODO change name: I_UpdateScattData???
+int I_compute_scatts(struct Cell_head *region, struct scCats * scatt_conds, const char ** bands, 
+                    int n_bands, struct scCats * scatt_plts, const char ** cats_rasts) 
+{/*scatt_plts must be zeros! */
+    if (n_bands != scatt_plts->n_bands ||
+        n_bands != scatt_conds->n_bands)
+        return -1;
+
+    const char *mapset;
+
+    FILE * f_cats_rasts[scatt_conds->n_a_cats];
+
+    char rast_cats_header[1024];//TODO magic number 
+
+    CELL * chunks[n_bands];
+    char * null_chunks[n_bands];
+
+    RASTER_MAP_TYPE data_type;
+
+    int max_chunk_size = 11234;
+    int nrows, i_band, id_band, n_a_bands, band_id, 
+        chunk_size, i_cat, id_cat, i, i_chunk, i_chars, i_row;
+    int fd_bands[n_bands];
+    for(i_band = 0; i_band < n_bands; i_band++)
+        fd_bands[i_band] = -1;
+    
+    int bands_ids[n_bands];
+    for(i_band = 0; i_band < n_bands; i_band++)
+        bands_ids[i_band] = -1;
+
+    int b_needed_bands[n_bands];  
+    memset(b_needed_bands, 0, (size_t)n_bands * sizeof(int));
+
+    get_needed_bands(scatt_conds, &b_needed_bands[0]);
+    get_needed_bands(scatt_plts, &b_needed_bands[0]);
+
+    Rast_set_window(region);
+
+    n_a_bands = 0;/*TODO realy needed?*/
+    for(id_band = 0; id_band < n_bands; id_band++)
+    {
+        if(b_needed_bands[id_band])
+        {
+            G_debug(3, "Opening raster no. %d with name: %s", id_band, bands[id_band]);
+
+            /* TODO solve  returns*/
+            if ((mapset = G_find_raster2(bands[id_band],"")) == NULL) {
+                free_compute_scatts_data(&fd_bands[0], &chunks[0], &null_chunks[0], &n_a_bands, 
+                                   &bands_ids[0], NULL, scatt_conds->n_a_cats);
+                return -1;
+            }
+            
+            if ((fd_bands[n_a_bands] = Rast_open_old(bands[id_band], mapset)) < 0) {
+                free_compute_scatts_data(&fd_bands[0], &chunks[0], &null_chunks[0], &n_a_bands, 
+                                   &bands_ids[0], NULL, scatt_conds->n_a_cats);
+                return -1;
+            }
+
+            //TODO check data type, minimum and maximum value
+            data_type = Rast_get_map_type(fd_bands[n_a_bands]);
+            //if(data_type != CELL)
+            //    return -1;
+
+            //What happens if it is not CELL tyoe
+            chunks[id_band] =  Rast_allocate_c_buf();
+            null_chunks[id_band] =  Rast_allocate_null_buf();
+
+            bands_ids[n_a_bands] = id_band;
+
+            ++n_a_bands;
+        }
+    }
+    
+    i_chars = sprintf(rast_cats_header, "P5\n%d\n%d\n1\n", region->cols, region->rows);
+    for(i_cat = 0; i_cat < scatt_conds->n_a_cats; i_cat++)
+    {
+
+        id_cat = scatt_conds->cats_ids[i_cat];
+        f_cats_rasts[i_cat] = fopen(cats_rasts[id_cat], "wb");
+        if(!f_cats_rasts[i_cat])
+        {
+            free_compute_scatts_data(&fd_bands[0], &chunks[0], &null_chunks[0], &n_a_bands, 
+                               &bands_ids[0], &f_cats_rasts[0], scatt_conds->n_a_cats);
+            return -1;
+        }
+
+        fwrite(rast_cats_header, sizeof(char), i_chars/sizeof(char), f_cats_rasts[i_cat]);
+        if (ferror(f_cats_rasts[i_cat])){
+            free_compute_scatts_data(&fd_bands[0], &chunks[0], &null_chunks[0], &n_a_bands, 
+                               &bands_ids[0], &f_cats_rasts[0], scatt_conds->n_a_cats);
+            return -1;
+        }
+    }
+
+    i_chunk = 0;
+
+    nrows = Rast_window_rows();
+    chunk_size = Rast_window_cols();
+
+    for (i_row = 0; i_row < nrows; i_row++)
+    {
+        for(i_band = 0; i_band < n_a_bands; i_band++)
+        {   
+            band_id = bands_ids[i_band];
+            G_debug(3, "Reading data for band %d", band_id);
+            Rast_get_c_row(fd_bands[i_band], chunks[band_id], i_row);
+            Rast_get_null_value_row (fd_bands[i_band], null_chunks[band_id], i_row);
+        }
+        if(compute_scatts_from_chunk(scatt_plts, scatt_conds, chunks, null_chunks, chunk_size, &f_cats_rasts[0], region, i_chunk) == -1) 
+        {
+            free_compute_scatts_data(&fd_bands[0], &chunks[0], &null_chunks[0], &n_a_bands, 
+                               &bands_ids[0], &f_cats_rasts[0], scatt_conds->n_a_cats);
+            return -1;
+        }
+ 
+    }
+    free_compute_scatts_data(&fd_bands[0], &chunks[0], &null_chunks[0], &n_a_bands, 
+                             &bands_ids[0], &f_cats_rasts[0], scatt_conds->n_a_cats);
+    return 0;    
+}
\ No newline at end of file



More information about the grass-commit mailing list