![]() |
|
Macros | |
#define | Set_gsl_handler gsl_error_handler_t *prior_handler = gsl_set_error_handler(apop_gsl_error); |
#define | Unset_gsl_handler gsl_set_error_handler(prior_handler); |
Typedefs | |
typedef int(* | apop_fn_ir )(apop_data *, void *) |
Functions | |
apop_data * | apop_data_alloc (const size_t size1, const size_t size2, const int size3) |
apop_data * | apop_data_calloc (const size_t size1, const size_t size2, const int size3) |
apop_data * | apop_matrix_to_data (gsl_matrix *m) |
apop_data * | apop_vector_to_data (gsl_vector *v) |
void | apop_text_free (char ***freeme, int rows, int cols) |
char | apop_data_free_base (apop_data *freeme) |
void | apop_data_memcpy (apop_data *out, const apop_data *in) |
apop_data * | apop_data_copy (const apop_data *in) |
apop_data * | apop_data_stack (apop_data *m1, apop_data *m2, char posn, char inplace) |
apop_data ** | apop_data_split (apop_data *in, int splitpoint, char r_or_c) |
void | apop_data_rm_columns (apop_data *d, int *drop) |
apop_data * | apop_data_prune_columns_base (apop_data *d, char **colnames) |
double * | apop_data_ptr (apop_data *data, int row, int col, const char *rowname, const char *colname, const char *page) |
double | apop_data_get (const apop_data *data, size_t row, int col, const char *rowname, const char *colname, const char *page) |
void | apop_gsl_error_for_set (const char *reason, const char *file, int line, int gsl_errno) |
int | apop_data_set (apop_data *data, size_t row, int col, const double val, const char *colname, const char *rowname, const char *page) |
int | apop_data_set_row (apop_data *d, apop_data *row, int row_number) |
void | apop_data_add_named_elmt (apop_data *d, char *name, double val) |
void | apop_data_add_names_base (apop_data *d, const char type, char const **names) |
int | apop_text_add (apop_data *in, const size_t row, const size_t col, const char *fmt,...) |
apop_data * | apop_text_alloc (apop_data *in, const size_t row, const size_t col) |
apop_data * | apop_data_transpose (apop_data *in, char transpose_text, char inplace) |
gsl_matrix * | apop_matrix_realloc (gsl_matrix *m, size_t newheight, size_t newwidth) |
gsl_vector * | apop_vector_realloc (gsl_vector *v, size_t newheight) |
apop_data * | apop_data_get_page (const apop_data *data, const char *title, const char match) |
apop_data * | apop_data_add_page (apop_data *dataset, apop_data *newpage, const char *title) |
apop_data * | apop_data_rm_page (apop_data *data, const char *title, const char free_p) |
apop_data * | apop_data_rm_rows (apop_data *in, int *drop, apop_fn_ir do_drop, void *drop_parameter) |
Variables | |
char * | apop_nul_string = "" |
void apop_data_add_named_elmt | ( | apop_data * | d, |
char * | name, | ||
double | val | ||
) |
A convenience function to add a named element to a data set. Many of Apophenia's testing procedures use this to easily produce a column of named parameters. It is public as a convenience.
d | The apop_data structure. Must not be NULL , but may be blank (as per allocation via apop_data_alloc ( ) ). |
name | The name to add |
val | the value to add to the set. |
NULL
), I will call apop_vector_realloc internally to make space. Add a page to a apop_data set. It gets a name so you can find it later.
dataset | The input data set, to which a page will be added. |
newpage | The page to append |
title | The name of the new page. Remember, this is truncated at 100 characters. |
NULL
data set and apop_opts.verbose >=1
.<Covariance>
, will be considered informational and ignored by search routines, missing data routines, et cetera. This is achieved by a rule in apop_data_pack and apop_data_unpack.Here is a toy example that establishes a baseline data set, adds a page, modifies it, and then later retrieves it.
apop_data* apop_data_alloc | ( | const size_t | size1, |
const size_t | size2, | ||
const int | size3 | ||
) |
Allocate a apop_data structure, to be filled with data.
apop_data_alloc(2,3,4)
: vector size, matrix rows, matrix cols. If the first argument is zero, you get a NULL
vector. apop_data_alloc(2,3)
, would allocate just a matrix, leaving the vector NULL
. apop_data_alloc(2)
, would allocate just a vector, leaving the matrix NULL
. apop_data_alloc()
, will produce a basically blank set, with out->matrix==out->vector==NULL
.For allocating the text part, see apop_text_alloc.
The weights
vector is set to NULL
. If you need it, allocate it via
out->error=='a' | Allocation error. The matrix, vector, or names couldn't be malloc ed, which probably means that you requested a very large data set. |
NULL
. But if even this much fails, your computer may be on fire and you should go put it out.apop_data* apop_data_calloc | ( | const size_t | size1, |
const size_t | size2, | ||
const int | size3 | ||
) |
Allocate a apop_data structure, to be filled with data; set everything in the allocated portion to zero. See apop_data_alloc for details.
out->error=='m' | malloc error; probably out of memory. |
Copy one apop_data structure to another. That is, all data is duplicated.
Basically a front-end for apop_data_memcpy for those who prefer this sort of syntax.
Unlike apop_data_memcpy, I do follow the more
pointer.
in | the input data |
out.error='a' | Allocation error. |
out.error='c' | Cyclic link: D->more == D (may be later in the chain, e.g., D->more->more = D->more ) You'll have only a partial copy. |
out.error='d' | Dimension error; should never happen. |
out.error='p' | Missing part error; should never happen.
|
char apop_data_free_base | ( | apop_data * | freeme | ) |
Free the elements of the given apop_data set and then the apop_data set itself. Intended to be used by apop_data_free, a macro that calls this to free elements, then sets the value to NULL
.
NULL
. For typical cases, that's slightly more useful than this function.freeme.error='c' | Circular linking is against the rules. If freeme->more == freeme , then I set freeme.error='c' and return. If you send in a structure like A -> B -> B, then both data sets A and B will be marked. |
0
on OK, 'c'
on error. It's good form to get a page from your data set by name, because you may not know the order for the pages, and the stepping through makes for dull code anyway (apop_data *page = dataset; while (page->more) page= page->more;
).
data | The apop_data set to use. No default; if NULL , gives a warning if apop_opts.verbose >=1 and returns NULL . |
title | The name of the page to retrieve. Default="Info" , which is the name of the page of additional estimation information returned by estimation routines (log likelihood, status, AIC, BIC, confidence intervals, ...). |
match | If 'c' , case-insensitive match (via strcasecmp ); if 'e' , exact match, if 'r' regular expression substring search (via apop_regex). Default='c' . |
NULL
.Copy one apop_data structure to another.
This function does not allocate the output structure or the vector, matrix, text, or weights elements—I assume you have already done this and got the dimensions right. I will assert that there is at least enough room in the destination for your data, and fail if the copy would write more elements than there are bins.
in
and out
have a more
pointer, also copy subsequent page(s). out | A structure that this function will fill. Must be preallocated with the appropriate sizes. |
in | The input data. |
out.error='d' | Dimension error; couldn't copy. |
out.error='p' | Part missing; e.g., in->matrix exists but out->matrix doesn't; couldn't copy. |
Keep only the columns of a data set that you name. This is the function called internally by the apop_data_prune_columns macro. In most cases, you'll want to use that macro. An example of the two uses demonstrating the difference:
d | The data set to prune. |
colnames | A null-terminated list of names to retain (i.e. the columns that shouldn't be pruned out). |
void apop_data_rm_columns | ( | apop_data * | d, |
int * | drop | ||
) |
Remove the columns set to one in the drop
vector. The returned data structure looks like it was modified in place, but the data matrix and the names are duplicated before being pared down, so if your data is taking up more than half of your memory, this may not work.
d | the apop_data structure to be pared down. |
drop | an array of ints. If use[7]==1, then column seven will be cut from the output. A reminder: calloc(in->size2 , sizeof(int)) will fill your array with zeros on allocation, and memset(use, 1, in->size2 * sizeof(int)) will quickly fill an array of ints with nonzero values. |
Remove the first page from an apop_data set that matches a given name.
data | The input data set, to which a page will be added. No default. If NULL , I return silently if apop_opts.verbose < 1 ; print an error otherwise. |
title | The case-insensitive name of the page to remove. Default: "Info" |
free_p | If 'y' , then apop_data_free the page. Default: 'y' . |
apop_data
page that I just pulled out. Thus, you can use this to pull a single page from a data set. I set that page's more
pointer to NULL
, to minimize any confusion about more-than-linear linked list topologies. If free_p=='y'
(the default) or the page is not found, return NULL
.->more
pointer in the apop_data set is not to fully implement a linked list, but primarily to allow you to staple auxiliary information to a main data set.apop_opts.verbose >= 1
. apop_data* apop_data_rm_rows | ( | apop_data * | in, |
int * | drop, | ||
apop_fn_ir | do_drop, | ||
void * | drop_parameter | ||
) |
Remove the rows set to one in the drop
vector or for which the do_drop
function returns one.
in | the apop_data structure to be pared down |
drop | a vector with as many elements as the max of the vector, matrix, or text parts of in , with a one marking those columns to be removed. |
do_drop | A function that returns one for rows to drop and zero for rows to not drop. A sample function: (onerow, .row=0, .col=0) , which can help to keep things readable). |
drop_parameter | If your do_drop function requires additional input, put it here and it will be passed through. |
NULL
vector
, matrix
, weight
, and text. Therefore, you may wish to check for NULL
elements after use. I remove rownames, but leave the other names, in case you want to add new data rows.NULL
, I return without doing anything, and print a warning if apop_opts.verbose >=1
. If you provide both, I will drop the row if either the vector has a one in that row's position, or if the function returns a nonzero value. Now that you've used Apop_r to pull a row from an apop_data set, this function lets you write that row to another position in the same data set or a different data set entirely.
The set written to must have the same form as the original:
d->names->rowct == row_number
(all rows up to row_number
have row names), then extend the list of row names by one to add the new name. Else, don't add the row name. apop_data_copy
or apop_name_copy
if you need to copy these names.If any of the source elements are NULL
, I won't bother to check that element in the destination.
_Thread_local
keyword) or a version of GCC with the __thread
extension enabled. Split one input apop_data structure into two.
For the opposite operation, see apop_data_stack.
in | The apop_data structure to split |
splitpoint | The index of what will be the first row/column of the second data set. E.g., if this is -1 and r_or_c=='c' , then the whole data set will be in the second data set; if this is the length of the matrix then the whole data set will be in the first data set. Another way to put it is that splitpoint will equal the number of rows/columns in the first matrix (unless it is -1, in which case the first matrix will have zero rows, or it is greater than the matrix's size, in which case it will have as many rows as the original). |
r_or_c | If this is 'r' or 'R', then put some rows in the first data set and some in the second; of 'c' or 'C', split columns into first and second data sets. |
NULL
pointer will be returned in that position. For example, for a data set of 50 rows, apop_data **out = apop_data_split(data, 100, 'r')
sets out[0] = apop_data_copy(data)
and out[1] = NULL
.more
pointer is ignored. apop_data->vector
is taken to be the -1st element of the matrix. apop_data_free(in)
after this. Put the first data set either on top of or to the left of the second data set.
The fn returns a new data set, meaning that at the end of this function, until you apop_data_free() the original data sets, you will be taking up twice as much memory. Plan accordingly.
For the opposite operation, see apop_data_split.
m1 | the upper/rightmost data set (default = NULL ) |
m2 | the second data set (default = NULL ) |
posn | If 'r', stack rows of m1's matrix above rows of m2's if 'c', stack columns of m1's matrix to left of m2's (default = 'r') |
inplace | If 'y' , use apop_matrix_realloc and apop_vector_realloc to modify m1 in place; see the caveats on those function. Otherwise, allocate a new vector, leaving m1 unmolested. (default='n') |
m1
out->error=='a' | Allocation error. |
out->error=='d' | Dimension error; couldn't make a complete copy. |
m2
is NULL
and inplace
is 'y'
, where you'll get the original m1
pointer back) more
is ignored. m1
, with the names for m2
appended to the row or column list, as appropriate. Transpose the matrix and text elements of the input data set, including the row/column names.
The vector and weights elements of the input data set are completely ignored (but see also apop_vector_to_matrix, which can convert a vector to a 1 X N matrix.) If copying, these other elements won't be present; if .inplace='y'
, it is up to you to handle these not-transposed elements correctly.
in | The input apop_data set. If NULL , I return NULL . Default is NULL . |
transpose_text | If 'y' , then also transpose the text element. Default is 'y' . |
inplace | If 'y' , transpose the input in place; if 'n' , produce a transposed copy, leaving the original untouched. Due to how gsl_matrix_transpose_memcpy works, a copy will still be made, then copied to the original location. Default is 'y' . |
inplace=='n'
, a newly alloced apop_data set, with the appropriately transposed matrix and/or text. The vector and weights elements will be NULL
. If transpose_text='n'
, then the text element of the output set will also be NULL
.inplace=='y'
, a pointer to the original data set, with matrix and (if transpose_text='y'
) text transposed and vector and weights left in place untouched.gsl_matrix
with no names or text, you may prefer to use gsl_matrix_transpose_memcpy
.gsl_matrix* apop_matrix_realloc | ( | gsl_matrix * | m, |
size_t | newheight, | ||
size_t | newwidth | ||
) |
This function will resize a gsl_matrix
to a new height or width.
Data in the matrix will be retained. If the new height or width is smaller than the old, then data in the later rows/columns will be cropped away (in a non–memory-leaking manner). If the new height or width is larger than the old, then new cells will be filled with garbage; it is your responsibility to zero out or otherwise fill new rows/columns before use.
Warning I: Using this function is basically bad form—especially when used in a for
loop that adds a column each time. A large number of realloc
s can take a noticeable amount of time. You are thus encouraged to make an effort to determine the size of your data beforehand.
Warning II: The gsl_matrix
is a versatile struct that can represent submatrices and other cuts from parent data. I can't deal with those, and check for such situations beforehand. [Besides, resizing a portion of a parent matrix makes no sense.]
m | The already-allocated matrix to resize. If you give me NULL , this becomes equivalent to gsl_matrix_alloc |
newheight,newwidth | The height and width you'd like the matrix to be. |
apop_data* apop_matrix_to_data | ( | gsl_matrix * | m | ) |
Deprecated; please do not use. Just use a compound literal:
int apop_text_add | ( | apop_data * | in, |
const size_t | row, | ||
const size_t | col, | ||
const char * | fmt, | ||
... | |||
) |
Add a string to the text element of an apop_data set. If you send me a NULL
string, I will write the value of apop_opts.nan_string
in the given slot. If there is already something in that slot, that string is freed, preventing memory leaks.
in | The apop_data set, that already has an allocated text element. |
row | The row |
col | The col |
fmt | The text to write. |
... | You can use a printf-style fmt and follow it with the usual variables to fill in. |
asprintf
), not a pointer to the input(s). asprintf(&(your_dataset->text[row][col]), "your string")
.This allocates an array of strings and puts it in the text
element of an apop_data set.
If the text
element already exists, then this is effectively a realloc
function, reshaping to the size you specify.
in | An apop_data set. It's OK to send in NULL , in which case an apop_data set with NULL matrix and vector elements is returned. |
row | the number of rows of text. |
col | the number of columns of text. |
NULL
, then this is a repeat of the input pointer. out->error=='a' | Allocation error. |
void apop_text_free | ( | char *** | freeme, |
int | rows, | ||
int | cols | ||
) |
Free a matrix of chars* (i.e., a char***). This is the form of the text element of the apop_data set, so you can use this for:
This is what apop_data_free
uses internally.
gsl_vector* apop_vector_realloc | ( | gsl_vector * | v, |
size_t | newheight | ||
) |
This function will resize a gsl_vector
to a new length.
Data in the vector will be retained. If the new height is smaller than the old, then data at the end of the vector will be cropped away (in a non–memory-leaking manner). If the new height is larger than the old, then new cells will be filled with garbage; it is your responsibility to zero out or otherwise fill them before use.
Warning I: Using this function is basically bad form—especially when used in a for
loop that adds an element each time. A large number of realloc
s can take a noticeable amount of time. You are thus encouraged to make an effort to determine the size of your data beforehand.
Warning II: The gsl_vector
is a versatile struct that can represent subvectors, matrix columns and other cuts from parent data. I can't deal with those, and check for such situations beforehand. [Besides, resizing a portion of a parent matrix makes no sense.]
v | The already-allocated vector to resize. If you give me NULL , this is equivalent to gsl_vector_alloc |
newheight | The height you'd like the vector to be. |
apop_data* apop_vector_to_data | ( | gsl_vector * | v | ) |
Deprecated; please do not use. Just use a compound literal, as in the code sample in the documentation for apop_matrix_to_data.