class_days_pc_clusters.c File Reference

Classification subroutine find the closest cluster of the principal components for each given day in EOF space. More...

#include <classif.h>
Include dependency graph for class_days_pc_clusters.c:

Go to the source code of this file.

Functions

void class_days_pc_clusters (int *days_class_cluster, double *pc_eof_days, double *eof_days_cluster, char *type, int neof, int ncluster, int ndays)
 Output a vector (dimension days) containing the closer (Euclidian distance) cluster number.

Detailed Description

Classification subroutine find the closest cluster of the principal components for each given day in EOF space.

Definition in file class_days_pc_clusters.c.


Function Documentation

void class_days_pc_clusters ( int *  days_class_cluster,
double *  pc_eof_days,
double *  eof_days_cluster,
char *  type,
int  neof,
int  ncluster,
int  ndays 
)

Output a vector (dimension days) containing the closer (Euclidian distance) cluster number.

The distance is computed as the distance between a day's principal components for each EOF and the cluster principal components for each EOF. To evaluate the closest one, all the square of the distances for all EOFs are summed for each cluster, before the square root is applied.

Parameters:
[out] days_class_cluster Cluster number associated for each day.
[in] pc_eof_days Principal Components of EOF (daily data).
[in] eof_days_cluster Clusters' centroid positions for each eof.
[in] type Type of distance used. Possible values: euclidian.
[in] neof Number of EOFs.
[in] ncluster Number of clusters.
[in] ndays Number of days in the pc_eof_days vector.

Definition at line 60 of file class_days_pc_clusters.c.

Referenced by generate_clusters(), wt_downscaling(), and wt_learning().

00061                                                           {
00072   double dist_min; /* Minimum distance found between a given day PC (summed over all EOF) and each cluster centroid. */
00073   int clust_dist_min; /* Cluster number which has the minimum distance dist_min */
00074   double dist_sum; /* Sum of distances (partial computation) over all EOFs */
00075   double val; /* Distance between a given day PC (for a particular EOF) and one cluster centroid. */
00076   double dist_clust; /* Distance (full computation of dist_sum). */
00077   
00078   int day; /* Loop counter for days */
00079   int clust; /* Loop counter for cluster */
00080   int eof; /* Loop counter for eofs */
00081 
00082   if ( !strcmp(type, "euclidian") ) {
00083     /* Euclidian distance type */
00084 
00085     /* Parse each day */
00086     for (day=0; day<ndays; day++) {
00087 
00088       /* Initialize */
00089       dist_min = 9999999999.0;
00090       clust_dist_min = 999;
00091 
00092 #if DEBUG >= 7
00093       (void) fprintf(stderr, "day=%d\n", day);
00094 #endif
00095 
00096       /* Parse each cluster */
00097       for (clust=0; clust<ncluster; clust++) {
00098 
00099 #if DEBUG >= 7
00100         (void) fprintf(stderr, "clust=%d\n", clust);
00101 #endif
00102 
00103         dist_sum = 0.0;
00104         /* Sum all distances (over EOF) between the PC of the day and the PC of the cluster centroid for each EOF respectively */
00105         for (eof=0; eof<neof; eof++) {
00106           val = pc_eof_days[day+eof*ndays] - eof_days_cluster[eof+clust*neof];
00107 #if DEBUG >= 9
00108           printf("%d %d %lf %lf\n",clust,eof,pc_eof_days[day+eof*ndays],eof_days_cluster[eof+clust*neof]);
00109 #endif
00110           /* Euclidian distance: square */
00111           dist_sum += (val * val);
00112         }
00113         /* Euclidian distance: square root of squares */
00114         dist_clust = sqrt(dist_sum);
00115 
00116 #if DEBUG >= 7
00117         (void) fprintf(stderr, "dist_clust=%lf\n", dist_clust);
00118 #endif
00119 
00120         /* Is it a cluster which has less distance as the minimum found yet ? */
00121         if (dist_clust < dist_min) {
00122           /* Save cluster number */
00123           clust_dist_min = clust;
00124           dist_min = dist_clust;
00125         }
00126       }
00127       if (clust_dist_min == 999) {
00128         /* Failing algorithm */
00129         (void) fprintf(stderr, "%s: ABORT: Impossible: no cluster was selected!! Problem in algorithm...\n", __FILE__);
00130         (void) abort();
00131       }
00132       /* Assign cluster with minimum distance to all EOFs for this day */
00133       days_class_cluster[day] = clust_dist_min;
00134 #if DEBUG >= 9
00135       (void) fprintf(stderr, "%s: day %d cluster %d\n", __FILE__, day, clust_dist_min);
00136 #endif
00137     }
00138   }
00139   else {
00140     /* Unknown distance type */
00141     (void) fprintf(stderr, "%s: ABORT: Unknown distance type=%s!!\n", __FILE__, type);
00142     (void) abort();
00143   }
00144 }


Generated on 12 May 2016 for DSCLIM by  doxygen 1.6.1