Patterns in static

Apophenia

Outline | Index |
 
Loading...
Searching...
No Matches
apop_db_sqlite.c
1
6#include <sqlite3.h>
7#include <string.h>
8
9sqlite3 *db=NULL; //There's only one SQLite database handle. Here it is.
10
11
12
14typedef struct StdDevCtx StdDevCtx;
15struct StdDevCtx {
16 double avg; /* avg of terms */
17 double avg2; /* avg of the squares of terms */
18 double avg3; /* avg of the cube of terms */
19 double avg4; /* avg of the fourth-power of terms */
20 int cnt; /* Number of terms counted */
21};
24static void twoStep(sqlite3_context *context, int argc, sqlite3_value **argv){
25 if (argc<1) return;
26 StdDevCtx *p = sqlite3_aggregate_context(context, sizeof(*p));
27 if (p && argv[0]){
28 double x = sqlite3_value_double(argv[0]);
29 double ratio = p->cnt/(p->cnt+1.0);
30 p->cnt++;
31 p->avg *= ratio;
32 p->avg2 *= ratio;
33 p->avg += x/(p->cnt +0.0);
34 p->avg2 += gsl_pow_2(x)/(p->cnt +0.0);
35 }
36}
37
38static void threeStep(sqlite3_context *context, int argc, sqlite3_value **argv){
39 if (argc<1) return;
40 StdDevCtx *p = sqlite3_aggregate_context(context, sizeof(*p));
41 if (p && argv[0]){
42 double x = sqlite3_value_double(argv[0]);
43 double ratio = p->cnt/(p->cnt+1.0);
44 p->cnt++;
45 p->avg *= ratio;
46 p->avg2 *= ratio;
47 p->avg3 *= ratio;
48 p->avg += x/p->cnt;
49 p->avg2 += gsl_pow_2(x)/p->cnt;
50 p->avg3 += gsl_pow_3(x)/p->cnt;
51 }
52}
53
54static void fourStep(sqlite3_context *context, int argc, sqlite3_value **argv){
55 if( argc<1 ) return;
56 StdDevCtx *p = sqlite3_aggregate_context(context, sizeof(*p));
57 if (p && argv[0]){
58 double x = sqlite3_value_double(argv[0]);
59 p->cnt++;
60 p->avg = (x + p->avg * (p->cnt-1.))/p->cnt;
61 p->avg2 = (gsl_pow_2(x)+ p->avg2 * (p->cnt-1.))/p->cnt;
62 p->avg3 = (gsl_pow_3(x)+ p->avg3 * (p->cnt-1.))/p->cnt;
63 p->avg4 = (gsl_pow_4(x)+ p->avg4 * (p->cnt-1.))/p->cnt;
64 }
65}
66
67static void stdDevFinalizePop(sqlite3_context *context){
68 StdDevCtx *p = sqlite3_aggregate_context(context, sizeof(*p));
69 if (p && p->cnt>1)
70 sqlite3_result_double(context, sqrt((p->avg2 - gsl_pow_2(p->avg))));
71 else if (p->cnt == 1)
72 sqlite3_result_double(context, 0);
73}
74
75static void varFinalizePop(sqlite3_context *context){
76 StdDevCtx *p = sqlite3_aggregate_context(context, sizeof(*p));
77 if( p && p->cnt>1 )
78 sqlite3_result_double(context, (p->avg2 - gsl_pow_2(p->avg)));
79 else if (p->cnt == 1)
80 sqlite3_result_double(context, 0);
81}
82
83static void stdDevFinalize(sqlite3_context *context){
84 StdDevCtx *p = sqlite3_aggregate_context(context, sizeof(*p));
85 if( p && p->cnt>1 ){
86 double rCnt = p->cnt;
87 sqlite3_result_double(context,
88 sqrt((p->avg2 - gsl_pow_2(p->avg))*rCnt/(rCnt-1.0)));
89 } else if (p->cnt == 1)
90 sqlite3_result_double(context, 0);
91}
92
93static void varFinalize(sqlite3_context *context){
94 StdDevCtx *p = sqlite3_aggregate_context(context, sizeof(*p));
95 if( p && p->cnt>1 ){
96 double rCnt = p->cnt;
97 sqlite3_result_double(context,
98 (p->avg2 - gsl_pow_2(p->avg))*rCnt/(rCnt-1.0));
99 } else if (p->cnt == 1)
100 sqlite3_result_double(context, 0);
101}
102
103static void skewFinalize(sqlite3_context *context){
104 StdDevCtx *p = sqlite3_aggregate_context(context, sizeof(*p));
105 if( p && p->cnt>1 ){
106 double rCnt = p->cnt;
107 sqlite3_result_double(context,
108 (p->avg3*rCnt - 3*p->avg2*p->avg*rCnt
109 + 2*rCnt * gsl_pow_3(p->avg)) * rCnt/((rCnt-1.0)*(rCnt-2.0)));
110 } else if (p->cnt == 1)
111 sqlite3_result_double(context, 0);
112}
113
114static void kurtFinalize(sqlite3_context *context){
115 StdDevCtx *p = sqlite3_aggregate_context(context, sizeof(*p));
116 if( p && p->cnt>1 ){
117 double n = p->cnt;
118 double kurtovern = p->avg4 - 4*p->avg3*p->avg
119 + 6 * p->avg2*gsl_pow_2(p->avg)
120 - 3* gsl_pow_4(p->avg);
121 double var = p->avg2 - gsl_pow_2(p->avg);
122 long double coeff0= n*n/(gsl_pow_3(n)*(gsl_pow_2(n)-3*n+3));
123 long double coeff1= n*gsl_pow_2(n-1)+ (6*n-9);
124 long double coeff2= n*(6*n-9);
125 sqlite3_result_double(context, coeff0*(coeff1 * kurtovern + coeff2 * gsl_pow_2(var)));
126 } else if (p->cnt == 1)
127 sqlite3_result_double(context, 0);
128}
129
130static void powFn(sqlite3_context *context, int argc, sqlite3_value **argv){
131 double base = sqlite3_value_double(argv[0]);
132 double exp = sqlite3_value_double(argv[1]);
133 sqlite3_result_double(context, pow(base, exp));
134}
135
136static void rngFn(sqlite3_context *context, int argc, sqlite3_value **argv){
137 Staticdef(gsl_rng *, rng, apop_rng_alloc(apop_opts.rng_seed++));
138 //sqlite3_result_double(context, gsl_rng_uniform(rng));
139 sqlite3_result_double(context, gsl_rng_uniform(apop_rng_get_thread(-1)));
140}
141
142#define sqfn(name) static void name##Fn(sqlite3_context *context, int argc, sqlite3_value **argv){ \
143 sqlite3_result_double(context, name(sqlite3_value_double(argv[0]))); }
144
145sqfn(sqrt) sqfn(exp) sqfn(log) sqfn(log10) sqfn(sin)
146sqfn(cos) sqfn(tan) sqfn(asin) sqfn(acos) sqfn(atan)
147
148
149static int apop_sqlite_db_open(char const *filename){
150 int status = sqlite3_open(filename ? filename : ":memory:", &db);
151 Apop_stopif(status, db=NULL; return status,
152 0, "The database %s didn't open.", filename ? filename : "in memory");
153 sqlite3_create_function(db, "stddev", 1, SQLITE_ANY, NULL, NULL, &twoStep, &stdDevFinalize);
154 sqlite3_create_function(db, "std", 1, SQLITE_ANY, NULL, NULL, &twoStep, &stdDevFinalizePop);
155 sqlite3_create_function(db, "stddev_samp", 1, SQLITE_ANY, NULL, NULL, &twoStep, &stdDevFinalize);
156 sqlite3_create_function(db, "stddev_pop", 1, SQLITE_ANY, NULL, NULL, &twoStep, &stdDevFinalizePop);
157 sqlite3_create_function(db, "var", 1, SQLITE_ANY, NULL, NULL, &twoStep, &varFinalize);
158 sqlite3_create_function(db, "var_samp", 1, SQLITE_ANY, NULL, NULL, &twoStep, &varFinalize);
159 sqlite3_create_function(db, "var_pop", 1, SQLITE_ANY, NULL, NULL, &twoStep, &varFinalizePop);
160 sqlite3_create_function(db, "variance", 1, SQLITE_ANY, NULL, NULL, &twoStep, &varFinalizePop);
161 sqlite3_create_function(db, "skew", 1, SQLITE_ANY, NULL, NULL, &threeStep, &skewFinalize);
162 sqlite3_create_function(db, "kurt", 1, SQLITE_ANY, NULL, NULL, &fourStep, &kurtFinalize);
163 sqlite3_create_function(db, "kurtosis", 1, SQLITE_ANY, NULL, NULL, &fourStep, &kurtFinalize);
164 sqlite3_create_function(db, "ln", 1, SQLITE_ANY, NULL, &logFn, NULL, NULL);
165 sqlite3_create_function(db, "ran", 0, SQLITE_ANY, NULL, &rngFn, NULL, NULL);
166 sqlite3_create_function(db, "pow", 2, SQLITE_ANY, NULL, &powFn, NULL, NULL);
167
168#define sqlink(name) sqlite3_create_function(db, #name , 1, SQLITE_ANY, NULL, &name##Fn, NULL, NULL);
169 sqlink(sqrt) sqlink(exp) sqlink(sin) sqlink(cos)
170 sqlink(tan) sqlink(asin) sqlink(acos) sqlink(atan) sqlink(log) sqlink(log10)
171 apop_query("pragma short_column_names");
172 return 0;
173}
174
176typedef struct { //for the apop_query_to_... functions.
177 int firstcall, namecol;
178 size_t currentrow;
179 apop_data *outdata;
180} callback_t;
183//This is the callback for apop_query_to_text.
184static int db_to_chars(void *qinfo,int argc, char **argv, char **column){
185 callback_t *qi= qinfo;
186 apop_data* d = qi->outdata; //alias. Allocated in calling fn.
187 int addnames = 0, ncshift=0;
188 if (!d->names->textct) addnames++;
189 if (qi->firstcall){
190 qi->firstcall = 0;
191 for(int i=0; i<argc; i++)
192 if (apop_opts.db_name_column && !strcasecmp(column[i], apop_opts.db_name_column)){
193 qi->namecol = i;
194 break;
195 }
196 }
197 int rows = d->textsize[0];
198 int cols = argc - (qi->namecol >= 0);
199 apop_text_alloc(d, rows+1, cols);//doesn't move d.
200 for (size_t jj=0; jj<argc; jj++)
201 if (jj == qi->namecol){
202 apop_name_add(d->names, argv[jj], 'r');
203 ncshift ++;
204 } else {
205 apop_text_set(d, rows, jj-ncshift, (argv[jj]==NULL)? apop_opts.nan_string: argv[jj]);
206 //Asprintf(&(d->text[rows][jj-ncshift]), "%s", (argv[jj]==NULL)? "NaN": argv[jj]);
207 if(addnames)
208 apop_name_add(d->names, column[jj], 't');
209 }
210 return 0;
211}
212
213apop_data * apop_sqlite_query_to_text(char *query){
214 char *err = NULL;
215 callback_t qinfo = {.outdata=apop_data_alloc(), .namecol=-1, .firstcall=1};
216 if (db==NULL) apop_db_open(NULL);
217 sqlite3_exec(db, query, db_to_chars, &qinfo, &err); ERRCHECK_SET_ERROR(qinfo.outdata)
218 if (qinfo.outdata->textsize[0]==0){
219 apop_data_free(qinfo.outdata);
220 return NULL;
221 }
222 return qinfo.outdata;
223}
224
226typedef struct {
227 apop_data *d;
228 int intypes[5];//names, vectors, mcols, textcols, weights.
229 int current, thisrow, error_thrown;
230 const char *instring;
231} apop_qt;
234static void count_types(apop_qt *in, const char *intypes){
235 int i = 0;
236 char c;
237 in->instring = intypes;
238 while ((c=intypes[i++]))
239 if (c=='n'||c=='N') in->intypes[0]++;
240 else if (c=='v'||c=='V') in->intypes[1]++;
241 else if (c=='m'||c=='M') in->intypes[2]++;
242 else if (c=='t'||c=='T') in->intypes[3]++;
243 else if (c=='w'||c=='W') in->intypes[4]++;
244 if (in->intypes[0]>1)
245 Apop_notify(1, "You asked apop_query_to_mixed data for multiple row names. I'll ignore all but the last one.");
246 if (in->intypes[1]>1)
247 Apop_notify(1, "You asked apop_query_to_mixed for multiple vectors. I'll ignore all but the last one.");
248 if (in->intypes[4]>1)
249 Apop_notify(1, "You asked apop_query_to_mixed for multiple weighting vectors. I'll ignore all but the last one.");
250}
251
252static int multiquery_callback(void *instruct, int argc, char **argv, char **column){
253 apop_qt *in = instruct;
254 char c;
255 int thistcol = 0,
256 thismcol = 0,
257 colct = 0,
258 i, addnames = 0;
259 in->thisrow ++;
260 if (!in->d) {
261 in->d = in->intypes[2]
262 ? apop_data_alloc(in->intypes[1], 1, in->intypes[2])
263 : apop_data_alloc(in->intypes[1]);
264 if (in->intypes[4])
265 in->d->weights = gsl_vector_alloc(1);
266 if (in->intypes[3]){
267 in->d->textsize[0] = 1;
268 in->d->textsize[1] = in->intypes[3];
269 in->d->text = malloc(sizeof(char***));
270 }
271 }
272 if (!(in->d->names->colct + in->d->names->textct + (in->d->names->vector!=NULL)))
273 addnames++;
274 if (in->d->textsize[1]){
275 in->d->textsize[0] = in->thisrow;
276 in->d->text = realloc(in->d->text, sizeof(char ***)*in->thisrow);
277 in->d->text[in->thisrow-1] = malloc(sizeof(char**) * in->d->textsize[1]);
278 }
279 if (in->intypes[2])
280 apop_matrix_realloc(in->d->matrix, in->thisrow, in->intypes[2]);
281 for (i=in->current=0; i< argc; i++){
282 c = in->instring[in->current++];
283 if (c=='n'||c=='N'){
284 apop_name_add(in->d->names, (argv[i]? argv[i] : "NaN") , 'r');
285 if(addnames)
286 apop_name_add(in->d->names, column[i], 'h');
287 } else if (c=='v'||c=='V'){
288 apop_vector_realloc(in->d->vector, in->thisrow);
289 apop_data_set(in->d, in->thisrow-1, -1,
290 argv[i] ? atof(argv[i]) : GSL_NAN);
291 if(addnames)
292 apop_name_add(in->d->names, column[i], 'v');
293 } else if (c=='m'||c=='M'){
294 apop_data_set(in->d, in->thisrow-1, thismcol++,
295 argv[i] ? atof(argv[i]) : GSL_NAN);
296 if(addnames)
297 apop_name_add(in->d->names, column[i], 'c');
298 } else if (c=='t'||c=='T'){
299 Asprintf(&(in->d->text[in->thisrow-1][thistcol++]), "%s",
300 argv[i] ? argv[i] : "NaN");
301 if(addnames)
302 apop_name_add(in->d->names, column[i], 't');
303 } else if (c=='w'||c=='W'){
304 apop_vector_realloc(in->d->weights, in->thisrow);
305 gsl_vector_set(in->d->weights, in->thisrow-1,
306 argv[i] ? atof(argv[i]) : GSL_NAN);
307 }
308 colct++;
309 }
310 int requested = in->intypes[0]+in->intypes[1]+in->intypes[2]+in->intypes[3]+in->intypes[4];
311 Apop_stopif(colct != requested, in->error_thrown='d'; return 1, 1,
312 "you asked for %i columns in your list of types(%s), but your query produced %u columns. "
313 "The remainder will be placed in the text section. Output data set's ->error element set to 'd'." , requested, in->instring, colct);
314 return 0;
315}
316
317apop_data *apop_sqlite_multiquery(const char *intypes, char *query){
318 Apop_stopif(!intypes, apop_return_data_error('t'), 0, "You gave me NULL for the list of input types. I can't work with that.");
319 Apop_stopif(!query, apop_return_data_error('q'), 0, "You gave me a NULL query. I can't work with that.");
320 char *err = NULL;
321 apop_qt info = { };
322 count_types(&info, intypes);
323 if (!db) apop_db_open(NULL);
324 sqlite3_exec(db, query, multiquery_callback, &info, &err);
325 Apop_stopif(info.error_thrown, if (!info.d) apop_data_alloc(); info.d->error='d'; return info.d,
326 0, "dimension error");
327 ERRCHECK_SET_ERROR(info.d)
328 return info.d;
329}
#define apop_data_free(freeme)
Definition docs/include/apop.h:192
int apop_db_open(char const *filename)
Definition apop_db.c:86
int apop_data_set(apop_data *data, size_t row, int col, const double val, const char *rowname, const char *colname, const char *page)
Definition apop_data.c:915
#define Apop_notify(verbosity,...)
Definition docs/include/apop.h:993
apop_data * apop_text_alloc(apop_data *in, const size_t row, const size_t col)
Definition apop_data.c:1063
gsl_matrix * apop_matrix_realloc(gsl_matrix *m, size_t newheight, size_t newwidth)
Definition apop_data.c:1259
#define Apop_stopif(test, onfail, level,...)
Definition docs/include/apop.h:1030
gsl_vector * apop_vector_realloc(gsl_vector *v, size_t newheight)
Definition apop_data.c:1309
apop_data * apop_data_alloc(const size_t size1, const size_t size2, const int size3)
Definition apop_data.c:34
gsl_rng * apop_rng_alloc(int seed)
Definition apop_bootstrap.c:19
int apop_name_add(apop_name *n, char const *add_me, char type)
Definition apop_name.c:42
int apop_text_set(apop_data *in, const size_t row, const size_t col, const char *fmt,...)
Definition apop_data.c:1035
Definition docs/include/apop.h:72
char * nan_string
Definition docs/include/apop.h:125
char * db_name_column
Definition docs/include/apop.h:124
Autogenerated by doxygen (Debian ).