annotate bayesianArraySrc/BayesianArrayStructure.cpp @ 37:9806a4f22fd0

Fixed bugs in the likelohoods that caused some to zero when no events found.
author Andrew N Robertson <andrew.robertson@eecs.qmul.ac.uk>
date Tue, 24 Apr 2012 14:16:01 +0100
parents eb43b2a007ea
children f5de07b4d733
rev   line source
andrew@0 1 /*
andrew@0 2 * BayesianArrayStructure.cpp
andrew@0 3 * midiCannamReader
andrew@0 4 *
andrew@0 5 * Created by Andrew on 17/07/2011.
andrew@0 6 * Copyright 2011 QMUL. All rights reserved.
andrew@0 7 *
andrew@0 8 */
andrew@0 9
andrew@0 10 //look at reset speed to one - what does this do? - get rid of?
andrew@0 11
andrew@0 12
andrew@0 13 #include "BayesianArrayStructure.h"
andrew@0 14
andrew@0 15 BayesianArrayStructure::BayesianArrayStructure(){
andrew@0 16 printf("Bayesian structure: DeFault constructor called");
andrew@0 17
andrew@4 18 usingIntegratedTempoEstimate = false;// use max index
andrew@7 19 updatingSpeedDistribution = false;
andrew@0 20
andrew@0 21 relativeSpeedLikelihoodStdDev = 5.0;
andrew@0 22
andrew@0 23 prior.createVector(1);
andrew@0 24 likelihood.createVector(1);
andrew@0 25 posterior.createVector(1);
andrew@0 26
andrew@0 27 speedPriorValue = 1.0;//default value for the speed prior
andrew@0 28 speedEstimate = speedPriorValue;
andrew@0 29
andrew@0 30 lastEventTime = 0;//ofGetElapsedTimeMillis();
andrew@0 31
andrew@0 32 tmpBestEstimate = 0;
andrew@5 33 crossUpdateTimeThreshold = 60;
andrew@37 34 priorWidth = 20;
andrew@36 35
andrew@36 36
andrew@0 37
andrew@0 38 }
andrew@0 39
andrew@0 40 BayesianArrayStructure::BayesianArrayStructure(int length){
andrew@0 41 printf("BAYESIAN STURCTURE CREATED LENGTH: %i\n", length);
andrew@0 42 //this constructor isnt called it seems
andrew@0 43 prior.createVector(length);
andrew@0 44 likelihood.createVector(length);
andrew@0 45 posterior.createVector(length);
andrew@0 46
andrew@0 47 lastEventTime = 0;
andrew@0 48
andrew@0 49
andrew@0 50 }
andrew@0 51
andrew@0 52
andrew@0 53
andrew@0 54 void BayesianArrayStructure::resetSize(int length){
andrew@0 55 printf("BAYESIAN STRUCTURE size is : %i\n", length);
andrew@0 56
andrew@0 57 prior.createVector(length);
andrew@0 58 likelihood.createVector(length);
andrew@0 59 posterior.createVector(length);
andrew@0 60
andrew@0 61 acceleration.createVector(length);
andrew@36 62 printf("scalar %f\n", posterior.scalar);
andrew@0 63
andrew@0 64 }
andrew@0 65
andrew@0 66 void BayesianArrayStructure::resetSpeedSize(int length){
andrew@0 67 printf("BAYESIAN reset SPEED size is : %i\n", length);
andrew@0 68
andrew@0 69 relativeSpeedPrior.createVector(length);
andrew@0 70 relativeSpeedLikelihood.createVector(length);
andrew@0 71 relativeSpeedPosterior.createVector(length);
andrew@0 72 tmpPosteriorForStorage.createVector(length);
andrew@0 73
andrew@0 74
andrew@0 75
andrew@0 76 }
andrew@0 77
andrew@0 78 void BayesianArrayStructure::setRelativeSpeedScalar(double f){
andrew@0 79 relativeSpeedPrior.scalar = f;
andrew@0 80 relativeSpeedPosterior.scalar = f;
andrew@0 81 relativeSpeedLikelihood.scalar = f;
andrew@0 82 }
andrew@0 83
andrew@0 84 void BayesianArrayStructure::resetSpeedToOne(){
andrew@0 85 relativeSpeedPrior.zero();
andrew@0 86 relativeSpeedPosterior.zero();
andrew@0 87 relativeSpeedLikelihood.zero();
andrew@0 88
andrew@0 89
andrew@20 90 // relativeSpeedPosterior.addGaussianShape(100, 8, 0.1);
andrew@20 91 // relativeSpeedPosterior.renormalise();
andrew@20 92 // relativeSpeedPosterior.getMaximum();
andrew@0 93
andrew@0 94 setSpeedPrior(speedPriorValue);
andrew@0 95 speedEstimate = speedPriorValue;
andrew@0 96
andrew@0 97 prior.zero();
andrew@0 98 posterior.zero();
andrew@0 99
andrew@0 100 posterior.addToIndex(0, 1);
andrew@0 101 posterior.renormalise();
andrew@0 102
andrew@0 103 }
andrew@0 104
andrew@0 105 void BayesianArrayStructure::setSpeedPrior(double f){
andrew@0 106 speedPriorValue = f;
andrew@0 107 int index = relativeSpeedPosterior.getRealTermsAsIndex(speedPriorValue);
andrew@0 108 relativeSpeedPosterior.zero();
andrew@37 109
andrew@37 110 //changed these a bit - still need to figure tempo process out properly
andrew@37 111 relativeSpeedPosterior.addGaussianShape(index, priorWidth, 1.9);
andrew@36 112 relativeSpeedPosterior.addGaussianShape(relativeSpeedPosterior.getRealTermsAsIndex(1), 4, 1);
andrew@36 113
andrew@4 114 printf("speed adding to index for 1 = %f\n", relativeSpeedPosterior.getRealTermsAsIndex(1));
andrew@36 115 // relativeSpeedPosterior.addToIndex(relativeSpeedPosterior.getRealTermsAsIndex(1), 0.5);
andrew@20 116 //relativeSpeedPosterior.addGaussianShapeFromRealTime(1, 3, 0.5);
andrew@11 117
andrew@11 118
andrew@11 119 //tmp debug test
andrew@11 120 //relativeSpeedPosterior.addToIndex(relativeSpeedPosterior.getRealTermsAsIndex(1.8), 0.15);
andrew@4 121
andrew@0 122 relativeSpeedPosterior.renormalise();
andrew@0 123 relativeSpeedPosterior.getMaximum();
andrew@11 124
andrew@11 125
andrew@0 126 relativeSpeedPrior.copyFromDynamicVector(relativeSpeedPosterior);
andrew@11 127
andrew@11 128
andrew@11 129
andrew@0 130 printf("BAYES STRUCTU ' SPEED PRIOR %f . index %i\n", speedPriorValue, index);
andrew@0 131
andrew@0 132 }
andrew@0 133
andrew@0 134
andrew@0 135 void BayesianArrayStructure::setPositionDistributionScalar(double f){
andrew@0 136 if (f > 0){
andrew@0 137 prior.scalar = f;
andrew@0 138 posterior.scalar = f;
andrew@0 139 likelihood.scalar = f;
andrew@36 140 printf("SET POS DISTBN SCALAR %f\n", f);
andrew@0 141 }
andrew@0 142 }
andrew@0 143
andrew@0 144 void BayesianArrayStructure::simpleExample(){
andrew@0 145 relativeSpeedPosterior.getMaximum();
andrew@0 146 relativeSpeedPrior.copyFromDynamicVector(relativeSpeedPosterior);
andrew@0 147 }
andrew@0 148
andrew@0 149 void BayesianArrayStructure::copyPriorToPosterior(){
andrew@0 150
andrew@0 151 for (int i = 0;i < prior.arraySize;i++){
andrew@0 152 posterior.array[i] = prior.array[i];
andrew@0 153 }
andrew@0 154 }
andrew@0 155
andrew@0 156 void BayesianArrayStructure::setStartPlaying(){
andrew@0 157
andrew@0 158 lastEventTime = 0;
andrew@0 159 bestEstimate = 0;
andrew@0 160 lastBestEstimateUpdateTime = 0;
andrew@3 161
andrew@0 162 if (*realTimeMode)
andrew@0 163 lastBestEstimateUpdateTime = ofGetElapsedTimeMillis();
andrew@0 164 //cannot just be zero - offline bug
andrew@0 165 //printf("start playing - best estimate %f\n", lastBestEstimateUpdateTime);
andrew@0 166
andrew@0 167 resetArrays();
andrew@0 168 }
andrew@0 169
andrew@0 170 void BayesianArrayStructure::resetArrays(){
andrew@0 171 //called when we start playing
andrew@0 172
andrew@0 173 prior.zero();
andrew@0 174 likelihood.zero();
andrew@0 175 posterior.zero();
andrew@0 176
andrew@0 177 updateCounter = 0;
andrew@0 178
andrew@0 179 posterior.offset = -1000;
andrew@0 180 setNewDistributionOffsets(0);
andrew@0 181
andrew@0 182 int zeroIndex = posterior.getRealTermsAsIndex(0);
andrew@3 183 printf("ZERO INDEX %i\n", zeroIndex);
andrew@0 184
andrew@31 185 posterior.addGaussianShapeFromRealTime(0, 50, 1);//one way to add at x msec
andrew@7 186 // posterior.addGaussianShape(posterior.getRealTermsAsIndex(10), 50, 1);//alternative way
andrew@3 187
andrew@0 188 likelihood.addConstant(1);
andrew@0 189
andrew@0 190 updateCounter = 0;
andrew@0 191
andrew@0 192
andrew@0 193 printf("bayes reset arrays - best estimate %f\n", lastBestEstimateUpdateTime);
andrew@0 194
andrew@0 195 setSpeedPrior(speedPriorValue);
andrew@3 196 relativeSpeedPosterior.copyFromDynamicVector(relativeSpeedPrior);
andrew@3 197
andrew@0 198 }
andrew@0 199
andrew@0 200
andrew@0 201 void BayesianArrayStructure::zeroArrays(){
andrew@0 202 prior.zero();
andrew@0 203 likelihood.zero();
andrew@0 204 posterior.zero();
andrew@0 205
andrew@0 206 relativeSpeedPrior.zero();
andrew@0 207 relativeSpeedPosterior.zero();
andrew@0 208 relativeSpeedLikelihood.zero();
andrew@0 209
andrew@0 210 }
andrew@0 211
andrew@4 212 /*
andrew@0 213 void BayesianArrayStructure::updateTmpBestEstimate(const double& timeDifference){
andrew@0 214 //input is the time since the start of playing
andrew@0 215 // double timeDiff = ofGetElapsedTimeMillis() - lastEventTime;//lastBestEstimateUpdateTime;
andrew@0 216 double timeDiff = timeDifference;
andrew@0 217 if (*realTimeMode)
andrew@0 218 timeDiff = ofGetElapsedTimeMillis() - lastBestEstimateUpdateTime;
andrew@0 219
andrew@0 220 double tmp = relativeSpeedPosterior.getIntegratedEstimate();
andrew@0 221 tmpBestEstimate = posterior.getIndexInRealTerms(posterior.MAPestimate) + timeDiff*relativeSpeedPosterior.getIndexInRealTerms(relativeSpeedPosterior.integratedEstimate);
andrew@0 222 //
andrew@0 223 //printf("tmp best %f and best %f time diff %f posterior MAP %f at speed %f\n", 0Estimate, bestEstimate, timeDifference, posterior.getIndexInRealTerms(posterior.MAPestimate), relativeSpeedPosterior.getIndexInRealTerms(relativeSpeedPosterior.integratedEstimate));
andrew@0 224 //lastBestEstimateUpdateTime = ofGetElapsedTimeMillis();
andrew@0 225 }
andrew@4 226 */
andrew@4 227
andrew@0 228 void BayesianArrayStructure::updateBestEstimate(const double& timeDifference){
andrew@0 229 // double timeDiff = ofGetElapsedTimeMillis() - lastEventTime;//
andrew@3 230 double tmp = bestEstimate;
andrew@10 231 // printf("best est routine: posterior offset %f\n", posterior.offset);
andrew@0 232
andrew@0 233 double timeDiff = timeDifference;
andrew@0 234
andrew@0 235 //Using timedifferencfe here will make it go wrong. Is time since beginning of playing
andrew@0 236
andrew@0 237 if (*realTimeMode)
andrew@0 238 timeDiff = ofGetElapsedTimeMillis() - lastBestEstimateUpdateTime;
andrew@0 239
andrew@0 240 //lastbest is time we started playing
andrew@4 241 /*
andrew@0 242 if (usingIntegratedTempoEstimate)
andrew@0 243 speedEstimateIndex = relativeSpeedPosterior.getIntegratedEstimate();
andrew@0 244 else
andrew@4 245 speedEstimateIndex = relativeSpeedPosterior.getMAPestimate();
andrew@4 246 */
andrew@4 247 speedEstimateIndex = getSpeedEstimateIndex();
andrew@0 248
andrew@0 249 speedEstimate = relativeSpeedPosterior.getIndexInRealTerms(speedEstimateIndex);
andrew@0 250 bestEstimate = posterior.getIndexInRealTerms(posterior.MAPestimate) + timeDiff*speedEstimate;
andrew@0 251
andrew@10 252 // printf("best estimate update from %f to %f; time diff %f MAP %i = %f ms speed index %f est %f SpeedxTime %f\n", tmp, bestEstimate, timeDiff,
andrew@10 253 // posterior.MAPestimate, posterior.getIndexInRealTerms(posterior.MAPestimate), speedEstimateIndex, speedEstimate, timeDiff*speedEstimate);
andrew@17 254
andrew@17 255 printf("BEST ESTIMATE %f\n", bestEstimate);
andrew@0 256 }
andrew@0 257
andrew@10 258
andrew@0 259 void BayesianArrayStructure::calculatePosterior(){
andrew@0 260 //posterior.doProduct(prior, likelihood);
andrew@4 261 assert(posterior.offset == prior.offset);
andrew@0 262
andrew@0 263 int i;
andrew@0 264 for (i = 0;i < posterior.length;i++){
andrew@0 265 posterior.array[i] = likelihood.array[i] * prior.array[i];
andrew@0 266 }
andrew@0 267
andrew@0 268 posterior.renormalise();
andrew@0 269
andrew@17 270 // printf("After CALC");
andrew@17 271 // printPostOffset();
andrew@6 272
andrew@0 273 }
andrew@0 274
andrew@0 275
andrew@4 276 double BayesianArrayStructure::getSpeedEstimateIndex(){
andrew@4 277 if (usingIntegratedTempoEstimate)
andrew@4 278 return relativeSpeedPosterior.getIntegratedEstimate();
andrew@4 279 else
andrew@4 280 return relativeSpeedPosterior.getMAPestimate();
andrew@4 281 }
andrew@0 282
andrew@0 283
andrew@0 284 void BayesianArrayStructure::setNewDistributionOffsets(const double& newOffset){
andrew@4 285
andrew@4 286 printf("prior offset was %f now %f\n", prior.offset, newOffset);
andrew@4 287
andrew@0 288 prior.offset = newOffset;
andrew@0 289 likelihood.offset = newOffset;
andrew@4 290
andrew@0 291 //posterior.offset = newOffset;
andrew@0 292 }
andrew@0 293
andrew@0 294
andrew@6 295 void BayesianArrayStructure::updateBayesianDistributions(const double& newEventTime){
andrew@6 296
andrew@6 297 //NEED TO CHECK HERE THAT THEY HAVE THE SAME OFFSETS
andrew@6 298 prior.copyFromDynamicVector(posterior);//try the otehr way
andrew@6 299
andrew@6 300 //bayesianStruct.copyPriorToPosterior();
andrew@6 301 //need to get new MAP position and set the offset of the arrays
andrew@6 302 //currently bestEstimate is the approx for the new MAP position
andrew@6 303 // int tmpMap = bayesianStruct.posterior.getMAPestimate();
andrew@6 304
andrew@6 305 double timeDifference = newEventTime - lastEventTime;
andrew@17 306 // printf("updating distributions at time %f diff %f offset %f\n", newEventTime, timeDifference, posterior.offset);
andrew@6 307
andrew@6 308 //addnoise to the tempo distribution
andrew@6 309 //bayesianStruct.decaySpeedDistribution(timeDifference);
andrew@6 310
andrew@7 311 if (timeDifference > 50 && updatingSpeedDistribution){
andrew@6 312 addGaussianNoiseToSpeedPosterior(timeDifference * 10.0 / 100.);
andrew@6 313 }
andrew@6 314
andrew@6 315 printPostOffset();
andrew@6 316
andrew@6 317 updateBestEstimate(timeDifference);
andrew@6 318 lastBestEstimateUpdateTime = newEventTime;//getTimeNow(timePlayed);
andrew@6 319
andrew@8 320 //set TARGETS - commented tenmporarily
andrew@8 321 setNewDistributionOffsets(max(0., bestEstimate - (prior.scalar*prior.arraySize/2)));
andrew@6 322
andrew@6 323 crossUpdateArrays(posterior, relativeSpeedPosterior, timeDifference);
andrew@6 324
andrew@6 325 //i.e. using the same offset as prior
andrew@6 326 posterior.offset = prior.offset;//
andrew@6 327
andrew@6 328 // float tmpPrior = max(0., bestEstimate - (prior.scalar*prior.arraySize/2));// prior.offset = max(0., bestEstimate - (prior.scalar*prior.arraySize/2));
andrew@6 329 // printf("Using prior offset of %f not %f\n", tmpPrior, prior.offset);
andrew@6 330
andrew@6 331 lastEventTime = newEventTime;//lastEventTime = ofGetElapsedTimeMillis();
andrew@6 332
andrew@6 333 }
andrew@6 334
andrew@6 335
andrew@10 336
andrew@10 337
andrew@0 338 void BayesianArrayStructure::crossUpdateArrays(DynamicVector& position, DynamicVector& speed, double timeDifference){
andrew@8 339
andrew@0 340 //set the cutoff for offset of position first! XXX
andrew@0 341
andrew@0 342 double timeDifferenceInPositionVectorUnits = timeDifference / prior.scalar;
andrew@0 343
andrew@17 344 //printf("CROSS UPDATE time diff %f ms is %f units; ", timeDifference, timeDifferenceInPositionVectorUnits);
andrew@0 345 prior.zero();//kill prior
andrew@4 346
andrew@4 347 // calculateNewPriorOffset(timeDifference);//dioesnt do anything
andrew@4 348
andrew@5 349 // printf("new prior offset %f and post offset %f\n", prior.offset, posterior.offset);
andrew@0 350
andrew@0 351 if (timeDifferenceInPositionVectorUnits > crossUpdateTimeThreshold)
andrew@0 352 complexCrossUpdate(timeDifferenceInPositionVectorUnits);
andrew@0 353 else
andrew@0 354 translateByMaximumSpeed(timeDifferenceInPositionVectorUnits);
andrew@0 355
andrew@0 356
andrew@0 357 updateCounter++;
andrew@8 358 prior.renormalise();//not strictly necessary??
andrew@0 359
andrew@0 360 }
andrew@0 361
andrew@0 362 void BayesianArrayStructure::complexCrossUpdate(const double& timeDifferenceInPositionVectorUnits){
andrew@4 363
andrew@5 364
andrew@17 365 // printf("before cross c : posterior map is %i = %f ms time diff pos vec %f\n", posterior.getMAPestimate(), posterior.getIndexInRealTerms(prior.getMAPestimate()), timeDifferenceInPositionVectorUnits);
andrew@5 366
andrew@0 367 int distanceMoved, newPriorIndex;
andrew@0 368
andrew@0 369 double speedValue = relativeSpeedPosterior.offset;
andrew@0 370
andrew@0 371 for (int i = 0;i < relativeSpeedPosterior.arraySize;i++){
andrew@0 372
andrew@0 373 // double speedValue = relativeSpeedPosterior.getIndexInRealTerms(i);//so for scalar 0.01, 50 -> speed value of 0.5
andrew@0 374
andrew@0 375 //so we have moved
andrew@0 376 distanceMoved = round(timeDifferenceInPositionVectorUnits * speedValue);//round the value
andrew@5 377 // printf("Speed value %f time %f gives distance %i\n", speedValue, timeDifferenceInPositionVectorUnits, distanceMoved);
andrew@5 378
andrew@0 379 if (relativeSpeedPosterior.array[i] != 0){
andrew@4 380
andrew@0 381 double speedContribution = relativeSpeedPosterior.array[i];
andrew@0 382
andrew@4 383 // printf("speed [%i](val[%f]) gives %f moved %i in %f units \n", i, relativeSpeedPosterior.array[i], speedValue, distanceMoved, timeDifferenceInPositionVectorUnits);
andrew@0 384
andrew@4 385 //1/2/12 deleted line
andrew@15 386 newPriorIndex = (int)posterior.millisToVectorUnits(posterior.offset - prior.offset) + distanceMoved;//i.e. where post[0] goes to in terms of prior at this speed
andrew@5 387 int postIndex = 0;//index of posterior that will contribute
andrew@5 388
andrew@4 389 while (postIndex < posterior.arraySize && newPriorIndex < prior.arraySize){
andrew@4 390
andrew@4 391 //did use a for loop
andrew@4 392 // for (postIndex = 0;postIndex < posterior.arraySize;postIndex++){
andrew@0 393 //old posterior contributing to new prior
andrew@4 394
andrew@4 395 //would use this method
andrew@4 396 //newPriorIndex = postIndex + posterior.offset - prior.offset + distanceMoved;
andrew@0 397
andrew@4 398 if (newPriorIndex >= 0){
andrew@0 399 prior.addToIndex(newPriorIndex, posterior.array[postIndex]*speedContribution);
andrew@4 400 // printf("speed index %i new prior index %i post val %f speed contrib %f dist %i\n", i, newPriorIndex, posterior.array[postIndex], speedContribution, distanceMoved);
andrew@0 401 }
andrew@4 402 //but we actually do this for simplicity
andrew@4 403 newPriorIndex++;
andrew@4 404 postIndex++;
andrew@4 405 }//end for. now while
andrew@0 406
andrew@0 407
andrew@0 408 }//if not zero
andrew@5 409
andrew@5 410 speedValue += relativeSpeedPosterior.scalar;//optimised line
andrew@0 411 //as we wanted:
andrew@0 412 // double speedValue = relativeSpeedPosterior.getIndexInRealTerms(i);//so for scalar 0.01, 50 -> speed value of 0.5
andrew@5 413
andrew@0 414 }//end speed
andrew@5 415
andrew@5 416
andrew@17 417 // printf("after cross c : prior map is %i = %f ms\n", prior.getMAPestimate(), prior.getIndexInRealTerms(prior.getMAPestimate()));
andrew@0 418 }
andrew@0 419
andrew@0 420
andrew@0 421
andrew@0 422 void BayesianArrayStructure::translateByMaximumSpeed(const double& timeDifferenceInPositionVectorUnits){
andrew@0 423
andrew@8 424
andrew@0 425 int distanceMoved, newPriorIndex;
andrew@8 426 double speedIndex = getSpeedEstimateIndex();
andrew@8 427 double speedValue = relativeSpeedPosterior.getIndexInRealTerms(speedIndex);
andrew@0 428
andrew@8 429 // double speedValue = relativeSpeedPosterior.getIndexInRealTerms(relativeSpeedPosterior.integratedEstimate);
andrew@8 430
andrew@0 431 //so for scalar 0.01, 50 -> speed value of 0.5
andrew@8 432 double speedContribution = relativeSpeedPosterior.array[(int)round(speedIndex)];
andrew@0 433 //so we have moved
andrew@0 434 distanceMoved = round(timeDifferenceInPositionVectorUnits * speedValue);//round the value
andrew@8 435
andrew@8 436 // printf("speed [%i] gives %f moved %i in %f units \n", i, speedValue, distanceMoved, timeDifferenceInPositionVectorUnits);
andrew@0 437
andrew@0 438 for (int postIndex = 0;postIndex < posterior.arraySize;postIndex++){
andrew@0 439 //old posterior contributing to new prior
andrew@15 440 newPriorIndex = (int)round(posterior.millisToVectorUnits(posterior.offset - prior.offset)) + postIndex + distanceMoved;
andrew@0 441 if (newPriorIndex >= 0 && newPriorIndex < prior.arraySize){
andrew@0 442 prior.addToIndex(newPriorIndex, posterior.array[postIndex]*speedContribution);
andrew@0 443 }
andrew@0 444
andrew@0 445 }
andrew@0 446
andrew@0 447 }
andrew@0 448
andrew@0 449 void BayesianArrayStructure::addGaussianNoiseToSpeedPosterior(const double& std_dev){
andrew@0 450 tmpPosteriorForStorage.copyFromDynamicVector(relativeSpeedPosterior);
andrew@0 451
andrew@0 452 for (int i = 0;i < relativeSpeedPosterior.length;i++){
andrew@0 453 tmpPosteriorForStorage.addGaussianShape(i, std_dev, relativeSpeedPosterior.array[i]);
andrew@0 454 }
andrew@0 455
andrew@0 456 tmpPosteriorForStorage.renormalise();
andrew@0 457
andrew@0 458 relativeSpeedPosterior.copyFromDynamicVector(tmpPosteriorForStorage);
andrew@0 459 }
andrew@0 460
andrew@0 461
andrew@0 462 void BayesianArrayStructure::addTriangularNoiseToSpeedPosterior(const double& std_dev){
andrew@0 463 tmpPosteriorForStorage.copyFromDynamicVector(relativeSpeedPosterior);
andrew@0 464
andrew@0 465 for (int i = 0;i < relativeSpeedPosterior.length;i++){
andrew@0 466 //adding a linear amount depending on distance
andrew@0 467 tmpPosteriorForStorage.addTriangularShape(i, std_dev*2.0, relativeSpeedPosterior.array[i]);
andrew@0 468 }
andrew@0 469
andrew@0 470 tmpPosteriorForStorage.renormalise();
andrew@0 471
andrew@0 472 relativeSpeedPosterior.copyFromDynamicVector(tmpPosteriorForStorage);
andrew@0 473 }
andrew@0 474
andrew@0 475 void BayesianArrayStructure::calculateNewPriorOffset(const double& timeDifference){
andrew@0 476
andrew@4 477 // double maxSpeed = relativeSpeedPosterior.getIndexInRealTerms(relativeSpeedPosterior.integratedEstimate);
andrew@4 478
andrew@4 479 double maxSpeed = relativeSpeedPosterior.getIndexInRealTerms(getSpeedEstimateIndex());//either integrated or MAP
andrew@0 480 // printf("Maxspeed is %f\n", maxSpeed);
andrew@0 481
andrew@0 482 double priorMax = posterior.getMaximum();
andrew@0 483 double distanceTravelled = maxSpeed * (timeDifference / prior.scalar);
andrew@0 484 double newMaxLocation = posterior.MAPestimate + distanceTravelled;
andrew@0 485 // printf("MAP: %i, tim df %f, distance %f, new location %f\n", posterior.MAPestimate, timeDifference, distanceTravelled, newMaxLocation);
andrew@0 486
andrew@0 487 }
andrew@0 488
andrew@0 489
andrew@0 490 void BayesianArrayStructure::decaySpeedDistribution(double timeDifference){
andrew@0 491
andrew@0 492 // commented for the moment
andrew@0 493 double relativeAmount = max(1.0, timeDifference/1000.);
andrew@0 494 // printf("decay %f around %i \n", timeDifference, relativeSpeedPosterior.MAPestimate);
andrew@0 495 relativeAmount *= speedDecayAmount;
andrew@0 496 relativeSpeedPosterior.renormalise();
andrew@0 497 relativeSpeedPosterior.addGaussianShape(relativeSpeedPosterior.MAPestimate, speedDecayWidth, relativeAmount);
andrew@0 498
andrew@0 499 relativeSpeedPosterior.renormalise();
andrew@0 500 double newMax = relativeSpeedPosterior.getMaximum();
andrew@0 501
andrew@0 502 //old code
andrew@0 503 // relativeSpeedPosterior.addGaussianShape(relativeSpeedPosterior.MAPestimate, speedDecayWidth, 10);
andrew@0 504 //relativeSpeedPosterior.addConstant(1);
andrew@0 505
andrew@0 506 /*
andrew@0 507 relativeSpeedPrior.copyFromDynamicVector(relativeSpeedPosterior);
andrew@0 508 relativeSpeedLikelihood.zero();
andrew@0 509 relativeSpeedLikelihood.addConstant(0.2);
andrew@0 510 relativeSpeedLikelihood.addGaussianShape(relativeSpeedPosterior.maximumValue, speedDecayWidth, relativeAmount);
andrew@0 511 relativeSpeedPosterior.doProduct(relativeSpeedPrior, relativeSpeedLikelihood);
andrew@0 512 relativeSpeedPosterior.renormalise();
andrew@0 513 */
andrew@0 514
andrew@0 515
andrew@0 516
andrew@0 517 }
andrew@0 518
andrew@0 519 void BayesianArrayStructure::setLikelihoodToConstant(){
andrew@0 520 //set new likelihood
andrew@0 521 relativeSpeedLikelihood.zero();
andrew@0 522 relativeSpeedLikelihood.addConstant(speedLikelihoodNoise);
andrew@0 523 }
andrew@0 524
andrew@0 525
andrew@0 526 void BayesianArrayStructure::updateTempoLikelihood(const double& speedRatio, const double& matchFactor){
andrew@0 527
andrew@0 528 //speedratio is speed of played relative to the recording
andrew@0 529
andrew@0 530 double index = relativeSpeedLikelihood.getRealTermsAsIndex(speedRatio);
andrew@0 531 // printf("index of likelihood would be %f for ratio %f\n", index, speedRatio);
andrew@0 532 if (index >= 0 && index < relativeSpeedPrior.length){
andrew@0 533 relativeSpeedLikelihood.addGaussianShape(index , relativeSpeedLikelihoodStdDev, matchFactor);
andrew@0 534 }
andrew@0 535 }
andrew@0 536
andrew@0 537
andrew@0 538
andrew@0 539 void BayesianArrayStructure::updateTempoDistribution(){
andrew@0 540
andrew@0 541 //copy posterior to prior
andrew@0 542 relativeSpeedPrior.copyFromDynamicVector(relativeSpeedPosterior);
andrew@0 543
andrew@0 544 //update
andrew@0 545 relativeSpeedPosterior.doProduct(relativeSpeedPrior, relativeSpeedLikelihood);
andrew@0 546
andrew@0 547 //normalise
andrew@0 548 relativeSpeedPosterior.renormalise();
andrew@0 549
andrew@0 550 relativeSpeedPosterior.getMaximum();
andrew@0 551
andrew@0 552 //relativeSpeedPosterior.updateIntegratedEstimate(); - could be swayed when off to side
andrew@0 553 //so now limited to where there is room sround the MAP estimate
andrew@0 554 relativeSpeedPosterior.updateLimitedIntegratedEstimate();
andrew@0 555
andrew@0 556 speedEstimate = relativeSpeedPosterior.getIndexInRealTerms(relativeSpeedPosterior.integratedEstimate);
andrew@0 557 }
andrew@0 558
andrew@0 559
andrew@0 560 void BayesianArrayStructure::calculateTempoUpdate(){
andrew@0 561 //copy posterior to prior
andrew@0 562 relativeSpeedPrior.copyFromDynamicVector(relativeSpeedPosterior);
andrew@0 563
andrew@0 564 //update
andrew@0 565 relativeSpeedPosterior.doProduct(relativeSpeedPrior, relativeSpeedLikelihood);
andrew@0 566
andrew@0 567 //normalise
andrew@0 568 relativeSpeedPosterior.renormalise();
andrew@0 569
andrew@0 570 relativeSpeedPosterior.getMaximum();
andrew@0 571
andrew@0 572 }
andrew@0 573
andrew@0 574
andrew@0 575 void BayesianArrayStructure::drawArrays(){
andrew@0 576
andrew@0 577 //bayesArray.drawFloatArray(&bayesArray.prior[0], 0, 200);
andrew@0 578 //bayesArray.drawFloatArray(&bayesArray.prior[0], 0, 200);
andrew@0 579
andrew@0 580 int displaySize = prior.arraySize;
andrew@0 581 ofSetColor(0,0,255);
andrew@0 582 prior.drawVector(0, displaySize);
andrew@0 583 ofSetColor(0,255,0);
andrew@0 584 likelihood.drawVector(0, displaySize);
andrew@0 585 ofSetColor(255,0,255);
andrew@0 586 posterior.drawVector(0, displaySize);
andrew@0 587
andrew@0 588
andrew@0 589
andrew@0 590 }
andrew@0 591
andrew@0 592
andrew@0 593 void BayesianArrayStructure::drawTempoArrays(){
andrew@0 594 ofSetColor(0,0,255);
andrew@0 595 // relativeSpeedPrior.drawVector(0, relativeSpeedPrior.arraySize);
andrew@0 596
andrew@0 597 ofSetColor(0,150,255);
andrew@0 598 relativeSpeedLikelihood.drawVector(0, relativeSpeedLikelihood.arraySize);
andrew@0 599
andrew@0 600 // relativeSpeedLikelihood.drawConstrainedVector(0, 199, 0, 1000);// relativeSpeedLikelihood.arraySize);
andrew@0 601 ofSetColor(255,0,0);
andrew@0 602 relativeSpeedPosterior.drawVector(0, relativeSpeedPosterior.arraySize);
andrew@0 603
andrew@0 604 // ofSetColor(0,0,255);
andrew@0 605 // tmpPosteriorForStorage.drawVector(0, tmpPosteriorForStorage.arraySize);
andrew@0 606
andrew@0 607 ofSetColor(255,255, 255);
andrew@0 608 ofLine(screenWidth/2, 0, screenWidth/2, ofGetHeight());//middle of screen
andrew@0 609
andrew@0 610 ofSetColor(25, 0, 250);
andrew@0 611 double fractionOfScreen = ((double)relativeSpeedPosterior.integratedEstimate / relativeSpeedPosterior.length);
andrew@0 612 ofLine(screenWidth * fractionOfScreen, 0, screenWidth * fractionOfScreen, ofGetHeight());
andrew@0 613
andrew@0 614 ofSetColor(255, 0, 20);
andrew@0 615 fractionOfScreen = ((double)relativeSpeedPosterior.MAPestimate / relativeSpeedPosterior.length);
andrew@0 616 ofLine(screenWidth * fractionOfScreen, 0, screenWidth * fractionOfScreen, ofGetHeight());
andrew@0 617
andrew@0 618
andrew@0 619 }
andrew@0 620
andrew@0 621
andrew@0 622 void BayesianArrayStructure::drawArraysRelativeToTimeframe(const double& startTimeMillis, const double& endTimeMillis){
andrew@0 623
andrew@0 624 screenWidth = ofGetWidth();
andrew@0 625
andrew@0 626 int startArrayIndex = 0;
andrew@0 627
andrew@0 628 if (prior.getIndexInRealTerms(prior.arraySize-1) > startTimeMillis){
andrew@0 629 //i.e. the array is on the page
andrew@0 630
andrew@0 631 while (prior.getIndexInRealTerms(startArrayIndex) < startTimeMillis){
andrew@0 632 startArrayIndex++;
andrew@0 633 }
andrew@0 634 int endArrayIndex = prior.arraySize-1;
andrew@0 635 //could find constraints here
andrew@0 636 if (prior.getIndexInRealTerms(prior.arraySize-1) > endTimeMillis)
andrew@0 637 endArrayIndex = (floor)((endTimeMillis - prior.offset)/prior.scalar);
andrew@0 638
andrew@0 639 //so we need to figure where start and end array are on screen
andrew@0 640 int startScreenPosition, endScreenPosition;
andrew@0 641 double screenWidthMillis = endTimeMillis - startTimeMillis;
andrew@0 642
andrew@0 643 startScreenPosition = (prior.getIndexInRealTerms(startArrayIndex) - startTimeMillis)*screenWidth/screenWidthMillis;
andrew@0 644 endScreenPosition = (double)(prior.getIndexInRealTerms(endArrayIndex) - startTimeMillis)*screenWidth/screenWidthMillis;
andrew@0 645
andrew@0 646 ofSetColor(0,0,100);
andrew@0 647 string relativeString = " offset "+ofToString(prior.offset, 1);//starttimes("+ofToString(startTimeMillis)+", "+ofToString(endTimeMillis);
andrew@0 648 relativeString += ": index "+ofToString(startArrayIndex)+" , "+ofToString(endArrayIndex)+" [";
andrew@0 649 // relativeString += ofToString(prior.getIndexInRealTerms(endArrayIndex), 3)+"] (sc-width:"+ofToString(screenWidthMillis, 1)+") ";
andrew@0 650 relativeString += " mapped to screen "+ofToString(startScreenPosition)+" , "+ofToString(endScreenPosition);
andrew@0 651 // ofDrawBitmapString(relativeString, 100, 180);
andrew@0 652
andrew@0 653
andrew@0 654
andrew@0 655 ofSetColor(100,100,100);//255, 255, 0);
andrew@0 656 likelihood.drawConstrainedVector(startArrayIndex, endArrayIndex, startScreenPosition, endScreenPosition);
andrew@0 657
andrew@0 658 // ofSetColor(0,0,200);
andrew@0 659 ofSetColor(170,170,170);//00,200);
andrew@0 660 prior.drawConstrainedVector(startArrayIndex, endArrayIndex, startScreenPosition, endScreenPosition);
andrew@0 661
andrew@0 662 ofSetColor(0,0,150);
andrew@0 663 // ofSetColor(200, 0, 0);
andrew@0 664 posterior.drawConstrainedVector(startArrayIndex, endArrayIndex, startScreenPosition, endScreenPosition);
andrew@0 665
andrew@0 666
andrew@0 667 // ofSetColor(0, 200, 255);
andrew@0 668 // acceleration.drawConstrainedVector(startArrayIndex, endArrayIndex, startScreenPosition, endScreenPosition);
andrew@0 669
andrew@0 670
andrew@0 671 }
andrew@0 672
andrew@0 673 }
andrew@0 674
andrew@0 675
andrew@6 676 void BayesianArrayStructure::printPostOffset(){
andrew@6 677 double tmp = posterior.getMAPestimate();
andrew@6 678 printf(" MAP index %i post offset %f == %f ms\n", posterior.MAPestimate, posterior.offset, posterior.getIndexInRealTerms(posterior.MAPestimate));
andrew@6 679 }
andrew@6 680
andrew@10 681 //PROJECT PRIOR CODE
andrew@11 682
andrew@11 683 void BayesianArrayStructure::projectDistribution(const double& newEventTime, const double& newAlignmentPosition, DynamicVector& projectedPrior){
andrew@10 684
andrew@10 685 projectedPrior.copyFromDynamicVector(posterior);
andrew@10 686
andrew@10 687 double timeDifference = newEventTime - lastEventTime;
andrew@10 688
andrew@10 689 // if (timeDifference > 50 && updatingSpeedDistribution){
andrew@10 690 // addGaussianNoiseToSpeedPosterior(timeDifference * 10.0 / 100.);
andrew@10 691 // }
andrew@10 692
andrew@10 693
andrew@11 694 // updateBestEstimate(timeDifference);
andrew@11 695 // lastBestEstimateUpdateTime = newEventTime;//getTimeNow(timePlayed);
andrew@10 696
andrew@10 697 //set TARGETS - commented tenmporarily
andrew@11 698
andrew@11 699 //setNewDistributionOffsets(max(0., bestEstimate - (prior.scalar*prior.arraySize/2)));
andrew@11 700 projectedPrior.offset = posterior.offset;//max(0., newAlignmentPosition - (projectedPrior.scalar*projectedPrior.arraySize/2));
andrew@11 701
andrew@11 702 // int timeDifference = newEventTime - lastEventTime;
andrew@11 703
andrew@11 704 double timeDifferenceInPositionVectorUnits = timeDifference / projectedPrior.scalar;
andrew@11 705
andrew@11 706 // printf("CROSS UPDATE time diff %f ms is %f units; ", timeDifference, timeDifferenceInPositionVectorUnits);
andrew@11 707 projectedPrior.zero();//kill prior
andrew@11 708
andrew@11 709 // calculateNewPriorOffset(timeDifference);//dioesnt do anything
andrew@11 710
andrew@11 711 // printf("new prior offset %f and post offset %f\n", prior.offset, posterior.offset);
andrew@11 712
andrew@11 713 if (timeDifferenceInPositionVectorUnits > crossUpdateTimeThreshold)
andrew@11 714 complexCrossUpdateProjection(projectedPrior, timeDifferenceInPositionVectorUnits);
andrew@11 715 else
andrew@11 716 translatePosteriorByMaximumSpeed(projectedPrior, timeDifferenceInPositionVectorUnits);
andrew@11 717
andrew@11 718
andrew@11 719 // updateCounter++;
andrew@11 720 projectedPrior.renormalise();//not strictly necessary??
andrew@11 721
andrew@11 722
andrew@10 723 //i.e. using the same offset as prior
andrew@11 724 // posterior.offset = prior.offset;//
andrew@10 725
andrew@10 726 // float tmpPrior = max(0., bestEstimate - (prior.scalar*prior.arraySize/2));// prior.offset = max(0., bestEstimate - (prior.scalar*prior.arraySize/2));
andrew@10 727 // printf("Using prior offset of %f not %f\n", tmpPrior, prior.offset);
andrew@10 728
andrew@11 729 // lastEventTime = newEventTime;//lastEventTime = ofGetElapsedTimeMillis();
andrew@10 730
andrew@10 731 }
andrew@10 732
andrew@10 733
andrew@10 734
andrew@11 735 void BayesianArrayStructure::complexCrossUpdateProjection(DynamicVector& projectedPrior, const double& timeDifferenceInPositionVectorUnits){
andrew@10 736
andrew@10 737 int distanceMoved, newPriorIndex;
andrew@11 738
andrew@10 739 double speedValue = relativeSpeedPosterior.offset;
andrew@10 740
andrew@10 741 for (int i = 0;i < relativeSpeedPosterior.arraySize;i++){
andrew@10 742
andrew@10 743 distanceMoved = round(timeDifferenceInPositionVectorUnits * speedValue);//round the value
andrew@10 744
andrew@10 745 if (relativeSpeedPosterior.array[i] != 0){
andrew@10 746
andrew@10 747 double speedContribution = relativeSpeedPosterior.array[i];
andrew@10 748
andrew@11 749 newPriorIndex = posterior.offset - projectedPrior.offset + distanceMoved;//i.e. where post[0] goes to in terms of prior at this speed
andrew@10 750 int postIndex = 0;//index of posterior that will contribute
andrew@10 751
andrew@11 752 while (postIndex < posterior.arraySize && newPriorIndex < projectedPrior.arraySize){
andrew@11 753
andrew@10 754 //would use this method
andrew@10 755 //newPriorIndex = postIndex + posterior.offset - prior.offset + distanceMoved;
andrew@10 756
andrew@10 757 if (newPriorIndex >= 0){
andrew@11 758 projectedPrior.addToIndex(newPriorIndex, posterior.array[postIndex]*speedContribution);
andrew@10 759 // printf("speed index %i new prior index %i post val %f speed contrib %f dist %i\n", i, newPriorIndex, posterior.array[postIndex], speedContribution, distanceMoved);
andrew@10 760 }
andrew@10 761 //but we actually do this for simplicity
andrew@10 762 newPriorIndex++;
andrew@10 763 postIndex++;
andrew@10 764 }//end for. now while
andrew@10 765
andrew@10 766
andrew@10 767 }//if not zero
andrew@10 768
andrew@10 769 speedValue += relativeSpeedPosterior.scalar;//optimised line
andrew@10 770 //as we wanted:
andrew@10 771 // double speedValue = relativeSpeedPosterior.getIndexInRealTerms(i);//so for scalar 0.01, 50 -> speed value of 0.5
andrew@10 772
andrew@10 773 }//end speed
andrew@11 774
andrew@10 775 }
andrew@10 776
andrew@10 777
andrew@10 778
andrew@11 779 void BayesianArrayStructure::translatePosteriorByMaximumSpeed(DynamicVector& translatedPosterior, const double& timeDifferenceInPositionVectorUnits){
andrew@10 780
andrew@10 781
andrew@10 782 int distanceMoved, newPriorIndex;
andrew@10 783 double speedIndex = getSpeedEstimateIndex();
andrew@10 784 double speedValue = relativeSpeedPosterior.getIndexInRealTerms(speedIndex);
andrew@10 785
andrew@10 786 // double speedValue = relativeSpeedPosterior.getIndexInRealTerms(relativeSpeedPosterior.integratedEstimate);
andrew@10 787
andrew@10 788 //so for scalar 0.01, 50 -> speed value of 0.5
andrew@10 789 double speedContribution = relativeSpeedPosterior.array[(int)round(speedIndex)];
andrew@10 790 //so we have moved
andrew@10 791 distanceMoved = round(timeDifferenceInPositionVectorUnits * speedValue);//round the value
andrew@10 792
andrew@10 793 // printf("speed [%i] gives %f moved %i in %f units \n", i, speedValue, distanceMoved, timeDifferenceInPositionVectorUnits);
andrew@10 794
andrew@10 795 for (int postIndex = 0;postIndex < posterior.arraySize;postIndex++){
andrew@10 796 //old posterior contributing to new prior
andrew@10 797 newPriorIndex = postIndex + posterior.offset - prior.offset + distanceMoved;
andrew@10 798 if (newPriorIndex >= 0 && newPriorIndex < prior.arraySize){
andrew@11 799 translatedPosterior.addToIndex(newPriorIndex, posterior.array[postIndex]*speedContribution);
andrew@10 800 }
andrew@10 801
andrew@10 802 }
andrew@10 803
andrew@10 804 }
andrew@11 805
andrew@10 806
andrew@10 807 //END PROJECT PRIOR CODE
andrew@10 808
andrew@10 809
andrew@11 810
andrew@0 811 /*
andrew@0 812
andrew@0 813 void BayesianArrayStructure::updateTempoDistribution(const double& speedRatio, const double& matchFactor){
andrew@0 814 //speedratio is speed of played relative to the recording
andrew@0 815
andrew@0 816 double index = relativeSpeedLikelihood.getRealTermsAsIndex(speedRatio);
andrew@0 817 // printf("\nindex of likelihood would be %f\n", index);
andrew@0 818 if (index >= 0 && index < relativeSpeedPrior.length){
andrew@0 819 //then we can do update
andrew@0 820
andrew@0 821 //set new likelihood
andrew@0 822 relativeSpeedLikelihood.zero();
andrew@0 823 relativeSpeedLikelihood.addConstant(speedLikelihoodNoise);
andrew@0 824
andrew@0 825 relativeSpeedLikelihood.addGaussianShape(index , 5, 0.5*matchFactor);
andrew@0 826
andrew@0 827
andrew@0 828 //copy posterior to prior
andrew@0 829 relativeSpeedPrior.copyFromDynamicVector(relativeSpeedPosterior);
andrew@0 830
andrew@0 831 //update
andrew@0 832 relativeSpeedPosterior.doProduct(relativeSpeedPrior, relativeSpeedLikelihood);
andrew@0 833
andrew@0 834 //normalise
andrew@0 835 relativeSpeedPosterior.renormalise();
andrew@0 836
andrew@0 837 relativeSpeedPosterior.getMaximum();
andrew@0 838 }//end if within range
andrew@0 839
andrew@0 840
andrew@0 841 }
andrew@0 842
andrew@0 843 */