annotate query.cpp @ 770:c54bc2ffbf92 tip

update tags
author convert-repo
date Fri, 16 Dec 2011 11:34:01 +0000
parents b9dbe4611dde
children
rev   line source
mas01cr@509 1 extern "C" {
mas01cr@509 2 #include "audioDB_API.h"
mas01cr@509 3 }
mas01cr@498 4 #include "audioDB-internals.h"
mas01cr@498 5 #include "accumulators.h"
mas01cr@239 6
mas01cr@498 7 bool audiodb_powers_acceptable(const adb_query_refine_t *r, double p1, double p2) {
mas01cr@498 8 if (r->flags & ADB_REFINE_ABSOLUTE_THRESHOLD) {
mas01cr@498 9 if ((p1 < r->absolute_threshold) || (p2 < r->absolute_threshold)) {
mas01cr@239 10 return false;
mas01cr@239 11 }
mas01cr@239 12 }
mas01cr@498 13 if (r->flags & ADB_REFINE_RELATIVE_THRESHOLD) {
mas01cr@498 14 if (fabs(p1-p2) > fabs(r->relative_threshold)) {
mas01cr@239 15 return false;
mas01cr@239 16 }
mas01cr@239 17 }
mas01cr@239 18 return true;
mas01cr@239 19 }
mas01cr@239 20
mas01cr@498 21 adb_query_results_t *audiodb_query_spec(adb_t *adb, const adb_query_spec_t *qspec) {
mas01cr@498 22 adb_qstate_internal_t qstate = {0};
mas01cr@498 23 qstate.allowed_keys = new std::set<std::string>;
mas01cr@498 24 adb_query_results_t *results;
mas01cr@498 25 if(qspec->refine.flags & ADB_REFINE_INCLUDE_KEYLIST) {
mas01cr@498 26 for(unsigned int k = 0; k < qspec->refine.include.nkeys; k++) {
mas01cr@498 27 qstate.allowed_keys->insert(qspec->refine.include.keys[k]);
mas01cr@498 28 }
mas01cr@498 29 } else {
mas01cr@498 30 for(unsigned int k = 0; k < adb->header->numFiles; k++) {
mas01cr@498 31 qstate.allowed_keys->insert((*adb->keys)[k]);
mas01cr@498 32 }
mas01cr@498 33 }
mas01cr@498 34 if(qspec->refine.flags & ADB_REFINE_EXCLUDE_KEYLIST) {
mas01cr@498 35 for(unsigned int k = 0; k < qspec->refine.exclude.nkeys; k++) {
mas01cr@498 36 qstate.allowed_keys->erase(qspec->refine.exclude.keys[k]);
mas01cr@498 37 }
mas01cr@498 38 }
mas01mc@292 39
mas01cr@498 40 switch(qspec->params.distance) {
mas01cr@498 41 case ADB_DISTANCE_DOT_PRODUCT:
mas01cr@498 42 switch(qspec->params.accumulation) {
mas01cr@498 43 case ADB_ACCUMULATION_DB:
mas01cr@498 44 qstate.accumulator = new DBAccumulator<adb_result_dist_gt>(qspec->params.npoints);
mas01cr@498 45 break;
mas01cr@498 46 case ADB_ACCUMULATION_PER_TRACK:
mas01cr@498 47 qstate.accumulator = new PerTrackAccumulator<adb_result_dist_gt>(qspec->params.npoints, qspec->params.ntracks);
mas01cr@498 48 break;
mas01cr@498 49 case ADB_ACCUMULATION_ONE_TO_ONE:
mas01cr@498 50 qstate.accumulator = new NearestAccumulator<adb_result_dist_gt>();
mas01cr@498 51 break;
mas01cr@498 52 default:
mas01cr@498 53 goto error;
mas01cr@239 54 }
mas01cr@239 55 break;
mas01cr@498 56 case ADB_DISTANCE_EUCLIDEAN_NORMED:
mas01cr@498 57 case ADB_DISTANCE_EUCLIDEAN:
mas01mc@768 58 case ADB_DISTANCE_KULLBACK_LEIBLER_DIVERGENCE:
mas01cr@498 59 switch(qspec->params.accumulation) {
mas01cr@498 60 case ADB_ACCUMULATION_DB:
mas01cr@498 61 qstate.accumulator = new DBAccumulator<adb_result_dist_lt>(qspec->params.npoints);
mas01cr@498 62 break;
mas01cr@498 63 case ADB_ACCUMULATION_PER_TRACK:
mas01cr@498 64 qstate.accumulator = new PerTrackAccumulator<adb_result_dist_lt>(qspec->params.npoints, qspec->params.ntracks);
mas01cr@498 65 break;
mas01cr@498 66 case ADB_ACCUMULATION_ONE_TO_ONE:
mas01cr@498 67 qstate.accumulator = new NearestAccumulator<adb_result_dist_lt>();
mas01cr@498 68 break;
mas01cr@498 69 default:
mas01cr@498 70 goto error;
mas01mc@263 71 }
mas01mc@263 72 break;
mas01cr@239 73 default:
mas01cr@498 74 goto error;
mas01mc@329 75 }
mas01cr@498 76
mas01cr@498 77 if((qspec->refine.flags & ADB_REFINE_RADIUS) && audiodb_index_exists(adb->path, qspec->refine.radius, qspec->qid.sequence_length)) {
mas01cr@498 78 if(audiodb_index_query_loop(adb, qspec, &qstate) < 0) {
mas01cr@498 79 goto error;
mas01cr@498 80 }
mas01cr@498 81 } else {
mas01cr@498 82 if(audiodb_query_loop(adb, qspec, &qstate)) {
mas01cr@498 83 goto error;
mas01cr@498 84 }
mas01mc@329 85 }
mas01mc@292 86
mas01cr@498 87 results = qstate.accumulator->get_points();
mas01cr@498 88
mas01cr@498 89 delete qstate.accumulator;
mas01cr@498 90 delete qstate.allowed_keys;
mas01cr@498 91
mas01cr@498 92 return results;
mas01cr@498 93
mas01cr@498 94 error:
mas01cr@498 95 if(qstate.accumulator)
mas01cr@498 96 delete qstate.accumulator;
mas01cr@498 97 if(qstate.allowed_keys)
mas01cr@498 98 delete qstate.allowed_keys;
mas01cr@498 99 return NULL;
mas01cr@239 100 }
mas01cr@239 101
mas01cr@498 102 int audiodb_query_free_results(adb_t *adb, const adb_query_spec_t *spec, adb_query_results_t *rs) {
mas01cr@498 103 free(rs->results);
mas01cr@498 104 free(rs);
mas01cr@498 105 return 0;
mas01cr@239 106 }
mas01cr@239 107
mas01cr@589 108 /* FIXME: we should check the return values from allocation */
mas01cr@498 109 static void audiodb_initialize_arrays(adb_t *adb, const adb_query_spec_t *spec, int track, unsigned int numVectors, double *query, double *data_buffer, double **D, double **DD) {
mas01cr@239 110 unsigned int j, k, l, w;
mas01cr@239 111 double *dp, *qp, *sp;
mas01mc@768 112 double a,b, tmp1;
mas01mc@768 113 #ifdef SYMMETRIC_KL
mas01mc@768 114 double tmp2;
mas01mc@768 115 #endif
mas01cr@239 116
mas01cr@498 117 const unsigned wL = spec->qid.sequence_length;
mas01cr@239 118
mas01cr@239 119 for(j = 0; j < numVectors; j++) {
mas01cr@239 120 // Sum products matrix
mas01cr@498 121 D[j] = new double[(*adb->track_lengths)[track]];
mas01cr@239 122 // Matched filter matrix
mas01cr@498 123 DD[j]=new double[(*adb->track_lengths)[track]];
mas01cr@239 124 }
mas01cr@239 125
mas01cr@239 126 // Dot product
mas01cr@239 127 for(j = 0; j < numVectors; j++)
mas01cr@498 128 for(k = 0; k < (*adb->track_lengths)[track]; k++){
mas01cr@498 129 qp = query + j * adb->header->dim;
mas01cr@498 130 sp = data_buffer + k * adb->header->dim;
mas01cr@239 131 DD[j][k] = 0.0; // Initialize matched filter array
mas01cr@239 132 dp = &D[j][k]; // point to correlation cell j,k
mas01cr@239 133 *dp = 0.0; // initialize correlation cell
mas01cr@498 134 l = adb->header->dim; // size of vectors
mas01mc@768 135 if (spec->params.distance!=ADB_DISTANCE_KULLBACK_LEIBLER_DIVERGENCE){
mas01mc@768 136 while(l--)
mas01mc@768 137 *dp += *qp++ * *sp++;
mas01mc@768 138 }
mas01mc@768 139 else{ // KL
mas01mc@768 140 while(l--){
mas01mc@768 141 a = *qp++;
mas01mc@768 142 b = *sp++;
mas01mc@768 143 tmp1 = a * log( a / b );
mas01mc@768 144 if(isnan(tmp1))
mas01mc@768 145 tmp1=0.0;
mas01mc@768 146 #ifdef SYMMETRIC_KL
mas01mc@768 147 tmp2 = b * log( b / a );
mas01mc@768 148 if(isnan(tmp2))
mas01mc@768 149 tmp2=0.0;
mas01mc@768 150 *dp += ( tmp1 + tmp2 ) / 2.0;
mas01mc@768 151 #else
mas01mc@768 152 *dp += tmp1;
mas01mc@768 153 #endif
mas01mc@768 154 }
mas01mc@768 155 }
mas01cr@239 156 }
mas01cr@239 157
mas01cr@239 158 double* spd;
mas01cr@663 159 if(!(spec->refine.flags & ADB_REFINE_HOP_SIZE)) {
mas01cr@239 160 for(w = 0; w < wL; w++) {
mas01cr@239 161 for(j = 0; j < numVectors - w; j++) {
mas01cr@239 162 sp = DD[j];
mas01cr@239 163 spd = D[j+w] + w;
mas01cr@498 164 k = (*adb->track_lengths)[track] - w;
mas01mc@292 165 while(k--)
mas01mc@292 166 *sp++ += *spd++;
mas01cr@239 167 }
mas01cr@239 168 }
mas01cr@663 169 } else {
mas01cr@675 170 uint32_t qhop = spec->refine.qhopsize;
mas01cr@675 171 qhop = qhop ? qhop : 1;
mas01cr@675 172 uint32_t ihop = spec->refine.ihopsize;
mas01cr@675 173 ihop = ihop ? ihop : 1;
mas01cr@239 174 for(w = 0; w < wL; w++) {
mas01cr@675 175 for(j = 0; j < numVectors - w; j += qhop) {
mas01cr@239 176 sp = DD[j];
mas01cr@239 177 spd = D[j+w]+w;
mas01cr@675 178 for(k = 0; k < (*adb->track_lengths)[track] - w; k += ihop) {
mas01cr@239 179 *sp += *spd;
mas01cr@675 180 sp += ihop;
mas01cr@675 181 spd += ihop;
mas01cr@239 182 }
mas01cr@239 183 }
mas01cr@239 184 }
mas01cr@239 185 }
mas01cr@239 186 }
mas01cr@239 187
mas01cr@498 188 static void audiodb_delete_arrays(int track, unsigned int numVectors, double **D, double **DD) {
mas01cr@239 189 if(D != NULL) {
mas01cr@239 190 for(unsigned int j = 0; j < numVectors; j++) {
mas01cr@239 191 delete[] D[j];
mas01cr@239 192 }
mas01cr@239 193 }
mas01cr@239 194 if(DD != NULL) {
mas01cr@239 195 for(unsigned int j = 0; j < numVectors; j++) {
mas01cr@239 196 delete[] DD[j];
mas01cr@239 197 }
mas01cr@239 198 }
mas01cr@239 199 }
mas01cr@239 200
mas01cr@498 201 int audiodb_read_data(adb_t *adb, int trkfid, int track, double **data_buffer_p, size_t *data_buffer_size_p) {
mas01cr@498 202 uint32_t track_length = (*adb->track_lengths)[track];
mas01cr@498 203 size_t track_size = track_length * sizeof(double) * adb->header->dim;
mas01cr@498 204 if (track_size > *data_buffer_size_p) {
mas01cr@239 205 if(*data_buffer_p) {
mas01cr@239 206 free(*data_buffer_p);
mas01cr@239 207 }
mas01cr@239 208 {
mas01cr@498 209 *data_buffer_size_p = track_size;
mas01cr@498 210 void *tmp = malloc(track_size);
mas01cr@239 211 if (tmp == NULL) {
mas01cr@498 212 goto error;
mas01cr@239 213 }
mas01cr@239 214 *data_buffer_p = (double *) tmp;
mas01cr@239 215 }
mas01cr@239 216 }
mas01cr@239 217
mas01cr@498 218 read_or_goto_error(trkfid, *data_buffer_p, track_size);
mas01cr@498 219 return 0;
mas01cr@498 220
mas01cr@498 221 error:
mas01cr@498 222 return 1;
mas01cr@239 223 }
mas01cr@239 224
mas01cr@498 225 int audiodb_track_id_datum(adb_t *adb, uint32_t track_id, adb_datum_t *d) {
mas01cr@498 226 off_t track_offset = (*adb->track_offsets)[track_id];
mas01cr@509 227 if(adb->header->flags & ADB_HEADER_FLAG_REFERENCES) {
mas01cr@498 228 /* create a reference/insert, then use adb_insert_create_datum() */
mas01cr@498 229 adb_reference_t reference = {0};
mas01cr@509 230 char features[ADB_MAXSTR], power[ADB_MAXSTR], times[ADB_MAXSTR];
mas01cr@596 231 lseek_set_or_goto_error(adb->fd, adb->header->dataOffset + track_id * ADB_FILETABLE_ENTRY_SIZE);
mas01cr@509 232 read_or_goto_error(adb->fd, features, ADB_MAXSTR);
mas01cr@498 233 reference.features = features;
mas01cr@509 234 if(adb->header->flags & ADB_HEADER_FLAG_POWER) {
mas01cr@596 235 lseek_set_or_goto_error(adb->fd, adb->header->powerTableOffset + track_id * ADB_FILETABLE_ENTRY_SIZE);
mas01cr@509 236 read_or_goto_error(adb->fd, power, ADB_MAXSTR);
mas01cr@498 237 reference.power = power;
mas01cr@498 238 }
mas01cr@509 239 if(adb->header->flags & ADB_HEADER_FLAG_TIMES) {
mas01cr@596 240 lseek_set_or_goto_error(adb->fd, adb->header->timesTableOffset + track_id * ADB_FILETABLE_ENTRY_SIZE);
mas01cr@509 241 read_or_goto_error(adb->fd, times, ADB_MAXSTR);
mas01cr@498 242 reference.times = times;
mas01cr@498 243 }
mas01cr@498 244 return audiodb_insert_create_datum(&reference, d);
mas01cr@498 245 } else {
mas01cr@498 246 /* initialize from sources of data that we already have */
mas01cr@498 247 d->nvectors = (*adb->track_lengths)[track_id];
mas01cr@498 248 d->dim = adb->header->dim;
mas01cr@498 249 d->key = (*adb->keys)[track_id].c_str();
mas01cr@498 250 /* read out stuff from the database tables */
mas01cr@596 251 malloc_and_fill_or_goto_error(double *, d->data, adb->header->dataOffset + track_offset, d->nvectors * d->dim * sizeof(double));
mas01cr@509 252 if(adb->header->flags & ADB_HEADER_FLAG_POWER) {
mas01cr@596 253 malloc_and_fill_or_goto_error(double *, d->power, adb->header->powerTableOffset + track_offset / d->dim, d->nvectors * sizeof(double));
mas01cr@598 254 } else {
mas01cr@598 255 d->power = NULL;
mas01cr@498 256 }
mas01cr@509 257 if(adb->header->flags & ADB_HEADER_FLAG_TIMES) {
mas01cr@596 258 malloc_and_fill_or_goto_error(double *, d->times, adb->header->timesTableOffset + 2 * track_offset / d->dim, 2 * d->nvectors * sizeof(double));
mas01cr@598 259 } else {
mas01cr@598 260 d->times = NULL;
mas01cr@498 261 }
mas01cr@498 262 return 0;
mas01cr@498 263 }
mas01cr@498 264 error:
mas01cr@580 265 audiodb_really_free_datum(d);
mas01cr@498 266 return 1;
mas01cr@498 267 }
mas01mc@292 268
mas01cr@498 269 int audiodb_datum_qpointers(adb_datum_t *d, uint32_t sequence_length, double **vector_data, double **vector, adb_qpointers_internal_t *qpointers) {
mas01cr@498 270 uint32_t nvectors = d->nvectors;
mas01cr@498 271
mas01cr@498 272 qpointers->nvectors = nvectors;
mas01cr@498 273
mas01cr@498 274 size_t vector_size = nvectors * sizeof(double) * d->dim;
mas01mc@762 275 *vector_data = new double[vector_size / sizeof(double) ];
mas01cr@498 276 memcpy(*vector_data, d->data, vector_size);
mas01cr@498 277
mas01mc@762 278 qpointers->l2norm_data = new double[vector_size / (sizeof(double)*d->dim)];
mas01cr@498 279 audiodb_l2norm_buffer(*vector_data, d->dim, nvectors, qpointers->l2norm_data);
mas01cr@498 280 audiodb_sequence_sum(qpointers->l2norm_data, nvectors, sequence_length);
mas01cr@498 281 audiodb_sequence_sqrt(qpointers->l2norm_data, nvectors, sequence_length);
mas01cr@498 282
mas01cr@498 283 if(d->power) {
mas01mc@762 284 qpointers->power_data = new double[vector_size / (sizeof(double)*d->dim)];
mas01cr@498 285 memcpy(qpointers->power_data, d->power, vector_size / d->dim);
mas01cr@498 286 audiodb_sequence_sum(qpointers->power_data, nvectors, sequence_length);
mas01cr@498 287 audiodb_sequence_average(qpointers->power_data, nvectors, sequence_length);
mas01cr@239 288 }
mas01cr@239 289
mas01cr@498 290 if(d->times) {
mas01cr@498 291 qpointers->mean_duration = new double[1];
mas01cr@498 292 *qpointers->mean_duration = 0;
mas01cr@498 293 for(unsigned int k = 0; k < nvectors; k++) {
mas01cr@498 294 *qpointers->mean_duration += d->times[2*k+1] - d->times[2*k];
mas01cr@239 295 }
mas01cr@498 296 *qpointers->mean_duration /= nvectors;
mas01cr@239 297 }
mas01cr@239 298
mas01cr@498 299 *vector = *vector_data;
mas01cr@498 300 qpointers->l2norm = qpointers->l2norm_data;
mas01cr@498 301 qpointers->power = qpointers->power_data;
mas01cr@498 302 return 0;
mas01cr@498 303 }
mas01cr@498 304
mas01cr@498 305 int audiodb_query_spec_qpointers(adb_t *adb, const adb_query_spec_t *spec, double **vector_data, double **vector, adb_qpointers_internal_t *qpointers) {
mas01cr@498 306 adb_datum_t *datum;
mas01cr@498 307 adb_datum_t d = {0};
mas01cr@498 308 uint32_t sequence_length;
mas01cr@498 309 uint32_t sequence_start;
mas01cr@498 310
mas01cr@498 311 datum = spec->qid.datum;
mas01cr@498 312 sequence_length = spec->qid.sequence_length;
mas01cr@498 313 sequence_start = spec->qid.sequence_start;
mas01cr@498 314
mas01cr@498 315 if(datum->data) {
mas01cr@498 316 if(datum->dim != adb->header->dim) {
mas01cr@498 317 return 1;
mas01cr@239 318 }
mas01cr@498 319 /* initialize d, and mark that nothing needs freeing later. */
mas01cr@498 320 d = *datum;
mas01cr@672 321 d.key = "";
mas01cr@498 322 datum = &d;
mas01cr@498 323 } else if (datum->key) {
mas01cr@498 324 uint32_t track_id;
mas01cr@498 325 if((track_id = audiodb_key_index(adb, datum->key)) == (uint32_t) -1) {
mas01cr@498 326 return 1;
mas01cr@498 327 }
mas01cr@498 328 audiodb_track_id_datum(adb, track_id, &d);
mas01cr@498 329 } else {
mas01cr@498 330 return 1;
mas01cr@239 331 }
mas01cr@239 332
mas01cr@498 333 /* FIXME: check the overflow logic here */
mas01cr@498 334 if(sequence_start + sequence_length > d.nvectors) {
mas01cr@498 335 if(datum != &d) {
mas01cr@580 336 audiodb_really_free_datum(&d);
mas01cr@498 337 }
mas01cr@498 338 return 1;
mas01cr@498 339 }
mas01cr@239 340
mas01cr@498 341 audiodb_datum_qpointers(&d, sequence_length, vector_data, vector, qpointers);
mas01cr@498 342
mas01cr@498 343 /* Finally, if applicable, set up the moving qpointers. */
mas01cr@498 344 if(spec->qid.flags & ADB_QID_FLAG_EXHAUSTIVE) {
mas01cr@498 345 /* the qpointers are already at the start, and so correct. */
mas01cr@498 346 } else {
mas01cr@498 347 /* adjust the qpointers to point to the correct place in the sequence */
mas01cr@498 348 *vector = *vector_data + spec->qid.sequence_start * d.dim;
mas01cr@498 349 qpointers->l2norm = qpointers->l2norm_data + spec->qid.sequence_start;
mas01cr@498 350 if(d.power) {
mas01cr@498 351 qpointers->power = qpointers->power_data + spec->qid.sequence_start;
mas01cr@239 352 }
mas01cr@498 353 qpointers->nvectors = sequence_length;
mas01cr@239 354 }
mas01cr@498 355
mas01cr@498 356 /* Clean up: free any bits of datum that we have ourselves
mas01cr@498 357 * allocated. */
mas01cr@498 358 if(datum != &d) {
mas01cr@580 359 audiodb_really_free_datum(&d);
mas01cr@498 360 }
mas01cr@498 361
mas01cr@498 362 return 0;
mas01cr@239 363 }
mas01cr@239 364
mas01cr@498 365 static int audiodb_set_up_dbpointers(adb_t *adb, const adb_query_spec_t *spec, adb_qpointers_internal_t *dbpointers) {
mas01cr@498 366 uint32_t nvectors = adb->header->length / (adb->header->dim * sizeof(double));
mas01cr@498 367 uint32_t sequence_length = spec->qid.sequence_length;
mas01mc@292 368
mas01cr@498 369 bool using_power = spec->refine.flags & (ADB_REFINE_ABSOLUTE_THRESHOLD|ADB_REFINE_RELATIVE_THRESHOLD);
mas01cr@498 370 bool using_times = spec->refine.flags & ADB_REFINE_DURATION_RATIO;
mas01cr@498 371 double *times_table = NULL;
mas01cr@498 372
mas01cr@498 373
mas01cr@498 374 dbpointers->nvectors = nvectors;
mas01cr@498 375 dbpointers->l2norm_data = new double[nvectors];
mas01cr@498 376
mas01cr@498 377 double *snpp = dbpointers->l2norm_data, *sppp = 0;
mas01cr@596 378 lseek_set_or_goto_error(adb->fd, adb->header->l2normTableOffset);
mas01cr@498 379 read_or_goto_error(adb->fd, dbpointers->l2norm_data, nvectors * sizeof(double));
mas01cr@498 380
mas01cr@498 381 if (using_power) {
mas01cr@509 382 if (!(adb->header->flags & ADB_HEADER_FLAG_POWER)) {
mas01cr@498 383 goto error;
mas01cr@498 384 }
mas01cr@498 385 dbpointers->power_data = new double[nvectors];
mas01cr@498 386 sppp = dbpointers->power_data;
mas01cr@596 387 lseek_set_or_goto_error(adb->fd, adb->header->powerTableOffset);
mas01cr@498 388 read_or_goto_error(adb->fd, dbpointers->power_data, nvectors * sizeof(double));
mas01mc@292 389 }
mas01mc@292 390
mas01cr@498 391 for(unsigned int i = 0; i < adb->header->numFiles; i++){
mas01cr@498 392 size_t track_length = (*adb->track_lengths)[i];
mas01cr@498 393 if(track_length >= sequence_length) {
mas01cr@498 394 audiodb_sequence_sum(snpp, track_length, sequence_length);
mas01cr@498 395 audiodb_sequence_sqrt(snpp, track_length, sequence_length);
mas01cr@498 396 if (using_power) {
mas01cr@498 397 audiodb_sequence_sum(sppp, track_length, sequence_length);
mas01cr@498 398 audiodb_sequence_average(sppp, track_length, sequence_length);
mas01cr@498 399 }
mas01mc@324 400 }
mas01cr@498 401 snpp += track_length;
mas01cr@498 402 if (using_power) {
mas01cr@498 403 sppp += track_length;
mas01mc@324 404 }
mas01mc@292 405 }
mas01mc@292 406
mas01cr@498 407 if (using_times) {
mas01cr@509 408 if(!(adb->header->flags & ADB_HEADER_FLAG_TIMES)) {
mas01cr@498 409 goto error;
mas01cr@498 410 }
mas01mc@292 411
mas01cr@498 412 dbpointers->mean_duration = new double[adb->header->numFiles];
mas01cr@498 413
mas01cr@596 414 malloc_and_fill_or_goto_error(double *, times_table, adb->header->timesTableOffset, 2 * nvectors * sizeof(double));
mas01cr@498 415 for(unsigned int k = 0; k < adb->header->numFiles; k++) {
mas01cr@498 416 size_t track_length = (*adb->track_lengths)[k];
mas01cr@498 417 unsigned int j;
mas01cr@498 418 dbpointers->mean_duration[k] = 0.0;
mas01cr@498 419 for(j = 0; j < track_length; j++) {
mas01cr@498 420 dbpointers->mean_duration[k] += times_table[2*j+1] - times_table[2*j];
mas01mc@292 421 }
mas01cr@498 422 dbpointers->mean_duration[k] /= j;
mas01mc@292 423 }
mas01cr@498 424
mas01cr@498 425 free(times_table);
mas01cr@498 426 times_table = NULL;
mas01mc@292 427 }
mas01cr@498 428
mas01cr@498 429 dbpointers->l2norm = dbpointers->l2norm_data;
mas01cr@498 430 dbpointers->power = dbpointers->power_data;
mas01cr@498 431 return 0;
mas01cr@498 432
mas01cr@498 433 error:
mas01cr@596 434 maybe_delete_array(dbpointers->l2norm_data);
mas01cr@596 435 maybe_delete_array(dbpointers->power_data);
mas01cr@596 436 maybe_delete_array(dbpointers->mean_duration);
mas01cr@596 437 maybe_free(times_table);
mas01cr@498 438 return 1;
mas01cr@498 439
mas01mc@292 440 }
mas01mc@292 441
mas01cr@498 442 int audiodb_query_queue_loop(adb_t *adb, const adb_query_spec_t *spec, adb_qstate_internal_t *qstate, double *query, adb_qpointers_internal_t *qpointers) {
mas01cr@498 443 adb_qpointers_internal_t dbpointers = {0};
mas01mc@292 444
mas01cr@498 445 uint32_t sequence_length = spec->qid.sequence_length;
mas01cr@498 446 bool power_refine = spec->refine.flags & (ADB_REFINE_ABSOLUTE_THRESHOLD|ADB_REFINE_RELATIVE_THRESHOLD);
mas01cr@239 447
mas01cr@498 448 if(qstate->exact_evaluation_queue->size() == 0) {
mas01cr@498 449 return 0;
mas01cr@239 450 }
mas01cr@239 451
mas01cr@498 452 /* We are guaranteed that the order of points is sorted by:
mas01cr@498 453 * {trackID, spos, qpos} so we can be relatively efficient in
mas01cr@498 454 * initialization of track data. We assume that points usually
mas01cr@498 455 * don't overlap, so we will use exhaustive dot product evaluation
mas01cr@498 456 * (instead of memoization of partial sums, as in query_loop()).
mas01cr@498 457 */
mas01cr@498 458 double dist;
mas01cr@498 459 double *dbdata = 0, *dbdata_pointer;
mas01cr@589 460 uint32_t currentTrack = 0x80000000; // KLUDGE: Initialize with a value outside of track index range
mas01cr@589 461 uint32_t npairs = qstate->exact_evaluation_queue->size();
mas01cr@498 462 while(npairs--) {
mas01cr@498 463 PointPair pp = qstate->exact_evaluation_queue->top();
mas01cr@498 464 if(currentTrack != pp.trackID) {
mas01cr@509 465 maybe_delete_array(dbdata);
mas01cr@509 466 maybe_delete_array(dbpointers.l2norm_data);
mas01cr@509 467 maybe_delete_array(dbpointers.power_data);
mas01cr@509 468 maybe_delete_array(dbpointers.mean_duration);
mas01cr@498 469 currentTrack = pp.trackID;
mas01cr@498 470 adb_datum_t d = {0};
mas01cr@498 471 if(audiodb_track_id_datum(adb, pp.trackID, &d)) {
mas01cr@498 472 delete qstate->exact_evaluation_queue;
mas01cr@498 473 return 1;
mas01cr@498 474 }
mas01cr@498 475 if(audiodb_datum_qpointers(&d, sequence_length, &dbdata, &dbdata_pointer, &dbpointers)) {
mas01cr@498 476 delete qstate->exact_evaluation_queue;
mas01cr@580 477 audiodb_really_free_datum(&d);
mas01cr@498 478 return 1;
mas01cr@498 479 }
mas01cr@580 480 audiodb_really_free_datum(&d);
mas01cr@498 481 }
mas01cr@589 482 uint32_t qPos = (spec->qid.flags & ADB_QID_FLAG_EXHAUSTIVE) ? pp.qpos : 0;
mas01cr@589 483 uint32_t sPos = pp.spos; // index into l2norm table
mas01cr@498 484 // Test power thresholds before computing distance
mas01cr@498 485 if( ( (!power_refine) || audiodb_powers_acceptable(&spec->refine, qpointers->power[qPos], dbpointers.power[sPos])) &&
mas01cr@498 486 ( qPos<qpointers->nvectors-sequence_length+1 && sPos<(*adb->track_lengths)[pp.trackID]-sequence_length+1 ) ){
mas01cr@498 487 // Compute distance
mas01mc@768 488 dist = 1.0e9;
mas01mc@768 489 if (spec->params.distance==ADB_DISTANCE_EUCLIDEAN_NORMED || spec->params.distance==ADB_DISTANCE_EUCLIDEAN)
mas01mc@768 490 dist = audiodb_dot_product(query + qPos*adb->header->dim, dbdata + sPos*adb->header->dim, adb->header->dim*sequence_length);
mas01mc@768 491 else if(spec->params.distance==ADB_DISTANCE_KULLBACK_LEIBLER_DIVERGENCE)
mas01mc@768 492 dist = audiodb_kullback_leibler(query + qPos*adb->header->dim, dbdata + sPos*adb->header->dim, adb->header->dim*sequence_length);
mas01cr@498 493 double qn = qpointers->l2norm[qPos];
mas01cr@498 494 double sn = dbpointers.l2norm[sPos];
mas01cr@498 495 switch(spec->params.distance) {
mas01cr@498 496 case ADB_DISTANCE_EUCLIDEAN_NORMED:
mas01cr@498 497 dist = 2 - (2/(qn*sn))*dist;
mas01cr@498 498 break;
mas01cr@498 499 case ADB_DISTANCE_EUCLIDEAN:
mas01cr@498 500 dist = qn*qn + sn*sn - 2*dist;
mas01cr@498 501 break;
mas01cr@498 502 }
mas01cr@498 503 if((!(spec->refine.flags & ADB_REFINE_RADIUS)) ||
mas01cr@509 504 dist <= (spec->refine.radius + ADB_DISTANCE_TOLERANCE)) {
mas01cr@498 505 adb_result_t r;
mas01cr@672 506 r.ikey = (*adb->keys)[pp.trackID].c_str();
mas01cr@672 507 r.qkey = spec->qid.datum->key;
mas01cr@498 508 r.dist = dist;
mas01cr@498 509 r.qpos = pp.qpos;
mas01cr@498 510 r.ipos = pp.spos;
mas01cr@610 511 if(qstate->set->find(r) == qstate->set->end()) {
mas01cr@610 512 qstate->set->insert(r);
mas01cr@610 513 qstate->accumulator->add_point(&r);
mas01cr@610 514 }
mas01cr@239 515 }
mas01cr@239 516 }
mas01cr@498 517 qstate->exact_evaluation_queue->pop();
mas01mc@292 518 }
mas01mc@474 519
mas01mc@315 520 // Cleanup
mas01cr@509 521 maybe_delete_array(dbdata);
mas01cr@509 522 maybe_delete_array(dbpointers.l2norm_data);
mas01cr@509 523 maybe_delete_array(dbpointers.power_data);
mas01cr@509 524 maybe_delete_array(dbpointers.mean_duration);
mas01cr@498 525 delete qstate->exact_evaluation_queue;
mas01cr@498 526 return 0;
mas01mc@292 527 }
mas01mc@292 528
mas01cr@498 529 int audiodb_query_loop(adb_t *adb, const adb_query_spec_t *spec, adb_qstate_internal_t *qstate) {
mas01cr@498 530
mas01cr@498 531 double *query, *query_data;
mas01cr@498 532 adb_qpointers_internal_t qpointers = {0}, dbpointers = {0};
mas01mc@292 533
mas01cr@498 534 bool power_refine = spec->refine.flags & (ADB_REFINE_ABSOLUTE_THRESHOLD|ADB_REFINE_RELATIVE_THRESHOLD);
mas01cr@239 535
mas01cr@509 536 if(adb->header->flags & ADB_HEADER_FLAG_REFERENCES) {
mas01cr@498 537 /* FIXME: actually it would be nice to support this mode of
mas01cr@498 538 * operation, but for now... */
mas01cr@498 539 return 1;
mas01cr@498 540 }
mas01mc@324 541
mas01cr@498 542 if(audiodb_query_spec_qpointers(adb, spec, &query_data, &query, &qpointers)) {
mas01cr@498 543 return 1;
mas01cr@498 544 }
mas01cr@239 545
mas01cr@498 546 if(audiodb_set_up_dbpointers(adb, spec, &dbpointers)) {
mas01cr@498 547 return 1;
mas01cr@498 548 }
mas01cr@239 549
mas01cr@663 550 unsigned j,k,track,trackOffset=0;
mas01cr@498 551 unsigned wL = spec->qid.sequence_length;
mas01cr@239 552 double **D = 0; // Differences query and target
mas01cr@239 553 double **DD = 0; // Matched filter distance
mas01cr@239 554
mas01cr@498 555 D = new double*[qpointers.nvectors]; // pre-allocate
mas01cr@498 556 DD = new double*[qpointers.nvectors];
mas01cr@239 557
mas01cr@675 558 unsigned qhop, ihop;
mas01cr@663 559
mas01cr@663 560 if(spec->refine.flags & ADB_REFINE_HOP_SIZE) {
mas01cr@675 561 qhop = spec->refine.qhopsize;
mas01cr@675 562 qhop = qhop ? qhop : 1;
mas01cr@675 563 ihop = spec->refine.ihopsize;
mas01cr@675 564 ihop = ihop ? ihop : 1;
mas01cr@663 565 } else {
mas01cr@675 566 qhop = 1;
mas01cr@675 567 ihop = 1;
mas01cr@663 568 }
mas01cr@239 569 off_t trackIndexOffset;
mas01cr@239 570
mas01cr@239 571 // Track loop
mas01cr@239 572 size_t data_buffer_size = 0;
mas01cr@239 573 double *data_buffer = 0;
mas01cr@498 574 lseek(adb->fd, adb->header->dataOffset, SEEK_SET);
mas01cr@239 575
mas01cr@498 576 std::set<std::string>::iterator keys_end = qstate->allowed_keys->end();
mas01cr@498 577 for(track = 0; track < adb->header->numFiles; track++) {
mas01cr@498 578 unsigned t = track;
mas01cr@498 579
mas01cr@498 580 while (qstate->allowed_keys->find((*adb->keys)[track]) == keys_end) {
mas01cr@498 581 track++;
mas01cr@498 582 if(track == adb->header->numFiles) {
mas01cr@498 583 goto loop_finish;
mas01cr@239 584 }
mas01cr@239 585 }
mas01cr@498 586 trackOffset = (*adb->track_offsets)[track];
mas01cr@498 587 if(track != t) {
mas01cr@498 588 lseek(adb->fd, adb->header->dataOffset + trackOffset, SEEK_SET);
mas01cr@498 589 }
mas01cr@498 590 trackIndexOffset = trackOffset / (adb->header->dim * sizeof(double)); // dbpointers.nvectors offset
mas01cr@239 591
mas01cr@498 592 if(audiodb_read_data(adb, adb->fd, track, &data_buffer, &data_buffer_size)) {
mas01cr@498 593 return 1;
mas01mc@292 594 }
mas01cr@498 595 if(wL <= (*adb->track_lengths)[track]) { // test for short sequences
mas01cr@498 596
mas01cr@498 597 audiodb_initialize_arrays(adb, spec, track, qpointers.nvectors, query, data_buffer, D, DD);
mas01mc@292 598
mas01cr@498 599 if((!(spec->refine.flags & ADB_REFINE_DURATION_RATIO)) ||
mas01cr@498 600 fabs(dbpointers.mean_duration[track]-qpointers.mean_duration[0]) < qpointers.mean_duration[0]*spec->refine.duration_ratio) {
mas01cr@239 601
mas01cr@239 602 // Search for minimum distance by shingles (concatenated vectors)
mas01cr@675 603 for(j = 0; j <= qpointers.nvectors - wL; j += qhop) {
mas01cr@675 604 for(k = 0; k <= (*adb->track_lengths)[track] - wL; k += ihop) {
mas01cr@498 605 double thisDist = 0;
mas01cr@498 606 double qn = qpointers.l2norm[j];
mas01cr@498 607 double sn = dbpointers.l2norm[trackIndexOffset + k];
mas01cr@498 608 switch(spec->params.distance) {
mas01cr@498 609 case ADB_DISTANCE_EUCLIDEAN_NORMED:
mas01cr@498 610 thisDist = 2-(2/(qn*sn))*DD[j][k];
mas01cr@498 611 break;
mas01cr@498 612 case ADB_DISTANCE_EUCLIDEAN:
mas01cr@498 613 thisDist = qn*qn + sn*sn - 2*DD[j][k];
mas01cr@498 614 break;
mas01cr@498 615 case ADB_DISTANCE_DOT_PRODUCT:
mas01cr@498 616 thisDist = DD[j][k];
mas01cr@498 617 break;
mas01cr@498 618 }
mas01cr@239 619 // Power test
mas01cr@498 620 if ((!power_refine) || audiodb_powers_acceptable(&spec->refine, qpointers.power[j], dbpointers.power[trackIndexOffset + k])) {
mas01cr@239 621 // radius test
mas01cr@498 622 if((!(spec->refine.flags & ADB_REFINE_RADIUS)) ||
mas01cr@509 623 thisDist <= (spec->refine.radius + ADB_DISTANCE_TOLERANCE)) {
mas01cr@498 624 adb_result_t r;
mas01cr@672 625 r.ikey = (*adb->keys)[track].c_str();
mas01cr@672 626 r.qkey = spec->qid.datum->key;
mas01cr@498 627 r.dist = thisDist;
mas01cr@498 628 if(spec->qid.flags & ADB_QID_FLAG_EXHAUSTIVE) {
mas01cr@498 629 r.qpos = j;
mas01cr@498 630 } else {
mas01cr@498 631 r.qpos = spec->qid.sequence_start;
mas01cr@498 632 }
mas01cr@498 633 r.ipos = k;
mas01cr@498 634 qstate->accumulator->add_point(&r);
mas01cr@239 635 }
mas01cr@239 636 }
mas01cr@239 637 }
mas01cr@239 638 }
mas01cr@239 639 } // Duration match
mas01cr@498 640 audiodb_delete_arrays(track, qpointers.nvectors, D, DD);
mas01cr@239 641 }
mas01cr@239 642 }
mas01cr@239 643
mas01cr@498 644 loop_finish:
mas01cr@498 645
mas01cr@239 646 free(data_buffer);
mas01cr@596 647 maybe_delete_array(query_data);
mas01cr@596 648 maybe_delete_array(qpointers.power_data);
mas01cr@596 649 maybe_delete_array(qpointers.l2norm_data);
mas01cr@596 650 maybe_delete_array(qpointers.mean_duration);
mas01cr@596 651 maybe_delete_array(dbpointers.power_data);
mas01cr@596 652 maybe_delete_array(dbpointers.l2norm_data);
mas01cr@596 653 maybe_delete_array(dbpointers.mean_duration);
mas01cr@596 654 maybe_delete_array(D);
mas01cr@596 655 maybe_delete_array(DD);
mas01cr@498 656
mas01cr@498 657 return 0;
mas01cr@239 658 }