1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61
|
#define USING_R 1
#include "S.h"
#ifdef USING_R
# include <R.h>
# include <Rdefines.h>
/* # include <Rinternals.h> */
# define S_EVALUATOR
#endif
#include "sp.h"
SEXP sp_zerodist(SEXP pp, SEXP pncol, SEXP zero) {
int i, j, k, ncol, nrow, nzero = 0, *which = NULL;
double **x, *xi, *xj, d, dist, zerodist2;
SEXP ret = NULL;
S_EVALUATOR
ncol = INTEGER_POINTER(pncol)[0];
nrow = LENGTH(pp)/ncol;
zerodist2 = NUMERIC_POINTER(zero)[0] * NUMERIC_POINTER(zero)[0];
x = (double **) malloc(nrow * sizeof(double *));
if (x == NULL)
error("could not allocate memory in zerodist");
for (i = 0; i < nrow; i++)
x[i] = &(NUMERIC_POINTER(pp)[i*ncol]);
for (i = 0; i < nrow; i++) {
xi = x[i];
for (j = 0; j < i; j++) {
xj = x[j];
for (k = 0, dist = 0.0; k < ncol; k++) {
d = (xi[k] - xj[k]);
dist += d * d;
}
if (dist <= zerodist2) {
which = (int *) realloc(which, (nzero+2) * sizeof(int));
if (which == NULL)
error("could not allocate memory in zerodist");
which[nzero] = j; /* lowest */
which[nzero + 1] = i;
nzero += 2;
}
}
R_CheckUserInterrupt();
}
free(x);
PROTECT(ret = NEW_INTEGER(nzero));
for (i = 0; i < nzero; i++)
INTEGER_POINTER(ret)[i] = which[i];
UNPROTECT(1);
if (which != NULL)
free(which);
return(ret);
}
/*
static SEXP sp_push_zerodist(SEXP ret, int i, int j) {
}
*/
|