5 #ifndef FML_CPU_DIMOPS_H 6 #define FML_CPU_DIMOPS_H 21 template <
typename REAL>
22 static inline void col_sum(
const len_t j,
const len_t m,
const REAL *x, REAL &sum)
26 #pragma omp simd reduction(+:sum) 27 for (len_t i=0; i<m; i++)
31 template <
typename REAL>
32 static inline void col_mean(
const len_t j,
const len_t m,
const REAL *x, REAL &mean)
35 col_sum(j, m, x, mean);
39 template <
typename REAL>
40 static inline void col_mean_and_var(
const len_t j,
const len_t m,
const REAL *x, REAL &mean, REAL &var)
45 for (len_t i=0; i<m; i++)
47 REAL dt = x[i + m*j] - mean;
48 mean += dt/((REAL) i+1);
49 var += dt * (x[i + m*j] - mean);
52 var = sqrt(var / ((REAL) m-1));
66 template <
typename REAL>
70 const len_t m = x.
nrows();
71 const len_t n = x.
ncols();
77 for (len_t j=0; j<n; j++)
80 for (len_t i=0; i<m; i++)
81 s_d[i] += x_d[i + m*j];
95 template <
typename REAL>
100 const len_t m = x.
nrows();
101 const len_t n = x.
ncols();
105 for (len_t i=0; i<m; i++)
119 template <
typename REAL>
123 const len_t m = x.
nrows();
124 const len_t n = x.
ncols();
130 for (len_t j=0; j<n; j++)
131 internals::col_sum(j, m, x_d, s_d[j]);
144 template <
typename REAL>
148 const len_t m = x.
nrows();
149 const len_t n = x.
ncols();
155 for (len_t j=0; j<n; j++)
156 internals::col_mean(j, m, x_d, s_d[j]);
163 template <
typename REAL>
167 const len_t m = x.
nrows();
168 const len_t n = x.
ncols();
170 for (len_t j=0; j<n; j++)
173 internals::col_mean(j, m, x_d, mean);
176 for (len_t i=0; i<m; i++)
177 x_d[i + m*j] -= mean;
181 template <
typename REAL>
185 const len_t m = x.
nrows();
186 const len_t n = x.
ncols();
188 for (len_t j=0; j<n; j++)
192 internals::col_mean_and_var(j, m, x_d, mean, var);
195 for (len_t i=0; i<m; i++)
200 template <
typename REAL>
204 const len_t m = x.
nrows();
205 const len_t n = x.
ncols();
207 for (len_t j=0; j<n; j++)
211 internals::col_mean_and_var(j, m, x_d, mean, var);
214 for (len_t i=0; i<m; i++)
215 x_d[i + m*j] = (x_d[i + m*j] - mean) / var;
231 template <
typename REAL>
234 if (rm_mean && rm_sd)
235 internals::center_and_scale(x);
237 internals::center(x);
void rowsums(const cpumat< REAL > &x, cpuvec< REAL > &s)
Compute the row sums.
Definition: dimops.hh:67
Vector class for data held on a single CPU.
Definition: cpuvec.hh:27
void colmeans(const cpumat< REAL > &x, cpuvec< REAL > &s)
Compute the column means.
Definition: dimops.hh:145
void resize(len_t size)
Resize the internal object storage.
Definition: cpuvec.hh:177
REAL * data_ptr()
Pointer to the internal array.
Definition: unimat.hh:35
T * data_ptr()
Pointer to the internal array.
Definition: univec.hh:26
void rowmeans(const cpumat< REAL > &x, cpuvec< REAL > &s)
Compute the row means.
Definition: dimops.hh:96
len_t nrows() const
Number of rows.
Definition: unimat.hh:31
void scale(const bool rm_mean, const bool rm_sd, cpumat< REAL > &x)
Remove the mean and/or the sd from a matrix.
Definition: dimops.hh:232
void colsums(const cpumat< REAL > &x, cpuvec< REAL > &s)
Compute the column sums.
Definition: dimops.hh:120
void fill_zero()
Set all values to zero.
Definition: cpuvec.hh:283
Row/column operations.
Definition: dimops.hh:17
Matrix class for data held on a single CPU.
Definition: cpumat.hh:32
len_t ncols() const
Number of columns.
Definition: unimat.hh:33