![]() |
|
Functions | |
apop_data * | apop_map (apop_data *in, apop_fn_d *fn_d, apop_fn_v *fn_v, apop_fn_r *fn_r, apop_fn_dp *fn_dp, apop_fn_vp *fn_vp, apop_fn_rp *fn_rp, apop_fn_dpi *fn_dpi, apop_fn_vpi *fn_vpi, apop_fn_rpi *fn_rpi, apop_fn_di *fn_di, apop_fn_vi *fn_vi, apop_fn_ri *fn_ri, void *param, int inplace, char part, int all_pages) |
gsl_vector * | apop_matrix_map (const gsl_matrix *m, double(*fn)(gsl_vector *)) |
void | apop_matrix_apply (gsl_matrix *m, void(*fn)(gsl_vector *)) |
gsl_vector * | apop_vector_map (const gsl_vector *v, double(*fn)(double)) |
void | apop_vector_apply (gsl_vector *v, void(*fn)(double *)) |
gsl_matrix * | apop_matrix_map_all (const gsl_matrix *in, double(*fn)(double)) |
void | apop_matrix_apply_all (gsl_matrix *in, void(*fn)(double *)) |
double | apop_vector_map_sum (const gsl_vector *in, double(*fn)(double)) |
double | apop_matrix_map_all_sum (const gsl_matrix *in, double(*fn)(double)) |
double | apop_matrix_map_sum (const gsl_matrix *in, double(*fn)(gsl_vector *)) |
double | apop_map_sum (apop_data *in, apop_fn_d *fn_d, apop_fn_v *fn_v, apop_fn_r *fn_r, apop_fn_dp *fn_dp, apop_fn_vp *fn_vp, apop_fn_rp *fn_rp, apop_fn_dpi *fn_dpi, apop_fn_vpi *fn_vpi, apop_fn_rpi *fn_rpi, apop_fn_di *fn_di, apop_fn_vi *fn_vi, apop_fn_ri *fn_ri, void *param, char part, int all_pages) |
These functions will pull each element of a vector or matrix, or each row of a matrix, and apply a function to the given element. See the data->map/apply section of the outline for many examples.
There are two types, which were developed at different times. The apop_map and apop_map_sum functions use variadic function inputs to cover a lot of different types of process depending on the inputs. Other functions with types in their names, like apop_matrix_map and apop_vector_apply, may be easier to use in some cases. With one exception, they use the same guts, so use whichever type is convenient.
Here are a few technical details of usage:
apop_opts.thread_count
is greater than one, then the matrix will be broken into chunks and each sent to a different thread. Notice that the GSL is generally threadsafe, and SQLite is threadsafe conditional on several commonsense caveats that you'll find in the SQLite documentation.
...sum functions are convenience functions that just call
...map and then add up the contents. Thus, you will need to have adequate memory for the allocation of the temp matrix/vector. apop_data* apop_map | ( | apop_data * | in, |
apop_fn_d * | fn_d, | ||
apop_fn_v * | fn_v, | ||
apop_fn_r * | fn_r, | ||
apop_fn_dp * | fn_dp, | ||
apop_fn_vp * | fn_vp, | ||
apop_fn_rp * | fn_rp, | ||
apop_fn_dpi * | fn_dpi, | ||
apop_fn_vpi * | fn_vpi, | ||
apop_fn_rpi * | fn_rpi, | ||
apop_fn_di * | fn_di, | ||
apop_fn_vi * | fn_vi, | ||
apop_fn_ri * | fn_ri, | ||
void * | param, | ||
int | inplace, | ||
char | part, | ||
int | all_pages | ||
) |
Apply a function to every element of a data set, matrix or vector; or, apply a vector-taking function to every row or column of a matrix.
There are a lot of options: your function could take any combination of a gsl_vector
, a double
, an apop_data, a parameter set, and the position of the element in the vector or matrix. As such, the function takes twelve function inputs, one for each combination of vector/matrix, params/no params, index/no index. Fortunately, because this function uses the Designated initializers syntax for inputs, you will specify only one.
For example, here is a function that will cut off each element of the input data to between .
fn_v | A function of the form double your_fn(gsl_vector *in) |
fn_d | A function of the form double your_fn(double in) |
fn_r | A function of the form double your_fn(apop_data *in) |
fn_vp | A function of the form double your_fn(gsl_vector *in, void *param) |
fn_dp | A function of the form double your_fn(double in, void *param) |
fn_rp | A function of the form double your_fn(apop_data *in, void *param) |
fn_vpi | A function of the form double your_fn(gsl_vector *in, void *param, int index) |
fn_dpi | A function of the form double your_fn(double in, void *param, int index) |
fn_rpi | A function of the form double your_fn(apop_data *in, void *param, int index) |
fn_vi | A function of the form double your_fn(gsl_vector *in, int index) |
fn_di | A function of the form double your_fn(double in, int index) |
fn_ri | A function of the form double your_fn(apop_data *in, int index) |
in | The input data set. If NULL , I'll return NULL immediately. |
param | A pointer to the parameters to be passed to those function forms taking a *param . |
part | Which part of the apop_data struct should I use?'v'==Just the vector 'm'==Every element of the matrix, in turn 'a'==Both 'v' and 'm' 'r'==Apply a function gsl_vector ![]() double to each row of the matrix'c'==Apply a function gsl_vector ![]() double to each column of the matrixDefault is 'a', but notice that I'll ignore a NULL vector or matrix, so if your data set has only a vector (for example), that's what I'll use. |
all_pages | If 'y' , then I follow the more pointer to subsequent pages, else I handle only the first page of data. [I abuse this for an internal semaphore, by the way, so your input must always be nonnegative and less than 1,000. Of course, 'y' and 'n' fit these rules fine.] Default: 'n' . |
inplace | If 'n' (the default), generate a new apop_data set for output, which will contain the mapped values (and the names from the original set). If 'y', modify in place. The double ![]() double versions, 'v' , 'm' , and 'a' , write to exactly the same location as before. The gsl_vector ![]() double versions, 'r' , and 'c' , will write to the vector. Be careful: if you are writing in place and there is already a vector there, then the original vector is lost.If 'v' (as in void), return NULL . (Default = 'n') |
r
in them, like fn_ri
, are row-by-row. I'll use Apop_r to get each row in turn, and send it to the function. The first implication is that your function should be expecting a apop_data set with exactly one row in it. The second is that part
is ignored: it only makes sense to go row-by-row.inplace='y'
, then you will be modifying your input data set, row by row; if you set inplace='n'
, then I will return an apop_data set whose vector
element is as long as your data set (i.e., as long as the longest of your text, vector, or matrix parts).apop_opts.thread_count
to a value greater than one, I will split the data set into as many chunks as you specify, and process them simultaneously. You need to watch out for the usual hang-ups about multithreaded programming, but if your data is iid, and each row's processing is independent of the others, you should have no problems. Bear in mind that generating threads takes some small overhead, so simple cases like adding a few hundred numbers will actually be slower when threading.out->error='p' | missing or mismatched parts error, such as NULL matrix when you sent a function acting on the matrix element. |
double apop_map_sum | ( | apop_data * | in, |
apop_fn_d * | fn_d, | ||
apop_fn_v * | fn_v, | ||
apop_fn_r * | fn_r, | ||
apop_fn_dp * | fn_dp, | ||
apop_fn_vp * | fn_vp, | ||
apop_fn_rp * | fn_rp, | ||
apop_fn_dpi * | fn_dpi, | ||
apop_fn_vpi * | fn_vpi, | ||
apop_fn_rpi * | fn_rpi, | ||
apop_fn_di * | fn_di, | ||
apop_fn_vi * | fn_vi, | ||
apop_fn_ri * | fn_ri, | ||
void * | param, | ||
char | part, | ||
int | all_pages | ||
) |
A function that effectively calls apop_map and returns the sum of the resulting elements. Thus, this function returns a single double
. See the apop_map page for details of the inputs, which are the same here, except that inplace
doesn't make sense—this function will always just add up the input function outputs.
See also the map/apply page for details.
NULL
. If apop_opts.verbose >= 2
I print a warning. void apop_matrix_apply | ( | gsl_matrix * | m, |
void(*)(gsl_vector *) | fn | ||
) |
Apply a function to every row of a matrix. The function that you input takes in a gsl_vector and returns nothing. apop_matrix_apply
will produce a vector view of each row, and send each row to your function.
m | The matrix |
fn | A function of the form void fn(gsl_vector* in) |
NULL
, this is a no-op and returns immediately. void apop_matrix_apply_all | ( | gsl_matrix * | in, |
void(*)(double *) | fn | ||
) |
Applies a function to every element in a matrix (as opposed to every row)
in | The matrix whose elements will be inputs to the function |
fn | A function with a form like void f(double *in) which will modify the data at the in-pointer in place. |
NULL
, this is a no-op and returns immediately. gsl_vector* apop_matrix_map | ( | const gsl_matrix * | m, |
double(*)(gsl_vector *) | fn | ||
) |
Map a function onto every row of a matrix. The function that you input takes in a gsl_vector and returns a double
. apop_matrix_map
will produce a vector view of each row, and send each row to your function. It will output a gsl_vector
holding your function's output for each row.
m | The matrix |
fn | A function of the form double fn(gsl_vector* in) |
gsl_vector
with the corresponding value for each row.NULL
matrix, I return NULL
. gsl_matrix* apop_matrix_map_all | ( | const gsl_matrix * | in, |
double(*)(double) | fn | ||
) |
Maps a function to every element in a matrix (as opposed to every row)
in | The matrix whose elements will be inputs to the function |
fn | A function with a form like double f(double in) . |
NULL
matrix, I return NULL
. double apop_matrix_map_all_sum | ( | const gsl_matrix * | in, |
double(*)(double) | fn | ||
) |
Like apop_matrix_map_all
, but returns the sum of the resulting mapped function. For example, apop_matrix_map_all_sum(v, isnan)
returns the number of elements of m
that are NaN
.
NULL
matrix, I return the sum of zero items: zero. double apop_matrix_map_sum | ( | const gsl_matrix * | in, |
double(*)(gsl_vector *) | fn | ||
) |
Like apop_matrix_map
, but returns the sum of the resulting mapped vector. For example, let log_like
be a function that returns the log likelihood of an input vector; then apop_matrix_map_sum(m, log_like)
returns the total log likelihood of the rows of m
.
NULL
matrix, I return the sum of zero items: zero. void apop_vector_apply | ( | gsl_vector * | v, |
void(*)(double *) | fn | ||
) |
Apply a function to every row of a matrix. The function that you input takes in a gsl_vector and returns nothing. apop_apply
will send a pointer to each element of your vector to your function.
v | The input vector |
fn | A function of the form void fn(double in) |
NULL
, this is a no-op and returns immediately. gsl_vector* apop_vector_map | ( | const gsl_vector * | v, |
double(*)(double) | fn | ||
) |
Map a function onto every element of a vector. The function that you input takes in a double
and returns a double
. apop_apply
will send each element to your function, and will output a gsl_vector
holding your function's output for each row.
v | The input vector |
fn | A function of the form double fn(double in) |
gsl_vector
(allocated by this function) with the corresponding value for each row.NULL
vector, I return NULL
. double apop_vector_map_sum | ( | const gsl_vector * | in, |
double(*)(double) | fn | ||
) |
Like apop_vector_map
, but returns the sum of the resulting mapped function. For example, apop_vector_map_sum(v, isnan)
returns the number of elements of v
that are NaN
.
NULL
vector, I return the sum of zero items: zero.