libcm is a C development framework with an emphasis on audio signal processing applications.
You can not select more than 25 topics Topics must start with a letter or number, can include dashes ('-') and can be up to 35 characters long.

cmProcTest.c 39KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520521522523524525526527528529530531532533534535536537538539540541542543544545546547548549550551552553554555556557558559560561562563564565566567568569570571572573574575576577578579580581582583584585586587588589590591592593594595596597598599600601602603604605606607608609610611612613614615616617618619620621622623624625626627628629630631632633634635636637638639640641642643644645646647648649650651652653654655656657658659660661662663664665666667668669670671672673674675676677678679680681682683684685686687688689690691692693694695696697698699700701702703704705706707708709710711712713714715716717718719720721722723724725726727728729730731732733734735736737738739740741742743744745746747748749750751752753754755756757758759760761762763764765766767768769770771772773774775776777778779780781782783784785786787788789790791792793794795796797798799800801802803804805806807808809810811812813814815816817818819820821822823824825826827828829830831832833834835836837838839840841842843844845846847848849850851852853854855856857858859860861862863864865866867868869870871872873874875876877878879880881882883884885886887888889890891892893894895896897898899900901902903904905906907908909910911912913914915916917918919920921922923924925926927928929930931932933934935936937938939940941942943944945946947948949950951952953954955956957958959960961962963964965966967968969970971972973974975976977978979980981982983984985986987988989990991992993994995996997998999100010011002100310041005100610071008100910101011101210131014101510161017101810191020102110221023102410251026102710281029103010311032103310341035103610371038103910401041104210431044104510461047104810491050105110521053105410551056105710581059106010611062106310641065106610671068106910701071107210731074107510761077107810791080108110821083108410851086108710881089109010911092109310941095109610971098109911001101110211031104110511061107110811091110111111121113111411151116111711181119112011211122112311241125112611271128112911301131113211331134113511361137113811391140114111421143114411451146114711481149115011511152115311541155115611571158115911601161116211631164116511661167116811691170117111721173117411751176117711781179118011811182118311841185118611871188118911901191119211931194119511961197119811991200120112021203120412051206120712081209121012111212121312141215121612171218121912201221122212231224122512261227122812291230123112321233123412351236123712381239124012411242124312441245124612471248124912501251125212531254125512561257125812591260126112621263126412651266126712681269127012711272127312741275127612771278127912801281128212831284128512861287128812891290129112921293129412951296129712981299130013011302130313041305130613071308130913101311131213131314131513161317131813191320132113221323132413251326132713281329133013311332133313341335133613371338133913401341134213431344134513461347134813491350135113521353135413551356135713581359136013611362136313641365136613671368136913701371137213731374137513761377137813791380138113821383138413851386138713881389139013911392139313941395139613971398139914001401140214031404140514061407140814091410141114121413141414151416141714181419142014211422142314241425142614271428142914301431143214331434143514361437143814391440144114421443144414451446144714481449145014511452145314541455145614571458145914601461146214631464146514661467146814691470147114721473147414751476147714781479148014811482148314841485148614871488148914901491149214931494149514961497149814991500150115021503150415051506150715081509151015111512151315141515151615171518151915201521152215231524152515261527152815291530153115321533153415351536153715381539154015411542154315441545154615471548154915501551155215531554155515561557
  1. #include "cmPrefix.h"
  2. #include "cmGlobal.h"
  3. #include "cmFloatTypes.h"
  4. #include "cmComplexTypes.h"
  5. #include "cmRpt.h"
  6. #include "cmErr.h"
  7. #include "cmCtx.h"
  8. #include "cmMem.h"
  9. #include "cmMallocDebug.h"
  10. #include "cmLinkedHeap.h"
  11. #include "cmSymTbl.h"
  12. #include "cmAudioFile.h"
  13. #include "cmTime.h"
  14. #include "cmMidi.h"
  15. #include "cmFile.h"
  16. #include "cmFileSys.h"
  17. #include "cmProcObj.h"
  18. #include "cmProcTemplate.h"
  19. #include "cmProc.h"
  20. #include "cmProc2.h"
  21. #include "cmMath.h"
  22. #include "cmVectOps.h"
  23. #include "cmKeyboard.h"
  24. #include "cmGnuPlot.h"
  25. #include "cmStack.h"
  26. #include "cmRbm.h"
  27. #include <time.h> // time()
  28. void cmTestPrint( cmRpt_t* rpt, const char* fmt, ... )
  29. {
  30. va_list vl;
  31. va_start(vl,fmt);
  32. if(rpt != NULL )
  33. cmRptVPrintf(rpt,fmt,vl);
  34. else
  35. vprintf(fmt,vl);
  36. va_end(vl);
  37. }
  38. void cmStandardizeTest(cmCtx_t* ctx )
  39. {
  40. cmRpt_t* rpt = ctx->err.rpt;
  41. unsigned rn = 3;
  42. unsigned cn = 2;
  43. double m[] = { 4, 5, 6, 7, 8, 9 };
  44. double uV[rn];
  45. double sdV[rn];
  46. cmVOD_PrintL("m", rpt, rn, cn, m );
  47. cmVOD_StandardizeRows(m, rn, cn, uV, sdV );
  48. cmVOD_PrintL("uV", rpt, 1, rn, uV);
  49. cmVOD_PrintL("sdV", rpt, 1, rn, sdV );
  50. cmVOD_PrintL("m", rpt, rn, cn, m );
  51. }
  52. void cmBinMtxFileTest(cmCtx_t* ctx )
  53. {
  54. const char* dataFn = "/home/kevin/temp/cmRbmData0.mtx";
  55. unsigned pointsN = 6; // count of data points in mtx (columns)
  56. unsigned dimN = 4; // dim of binary data matrix (rows)
  57. double probV[] = {.2,.5,.8,.6}; // probabilities of generating a 1.0
  58. assert(sizeof(probV)/sizeof(probV[0]) == dimN);
  59. // alloc the data matrix
  60. double* data0M = cmMemAllocZ( double, dimN*pointsN );
  61. unsigned i,j;
  62. // generate a stochastic binary data matrix according to prob's in probV[]
  63. for(i=0; i<pointsN; ++i)
  64. for(j=0; j<dimN; ++j)
  65. data0M[ i*dimN + j ] = rand() < (probV[j] * RAND_MAX);
  66. // write the binary matrix file
  67. cmBinMtxFileWrite(dataFn, dimN, pointsN, NULL, data0M, NULL, ctx->err.rpt );
  68. // print the original data
  69. cmVOD_PrintL("data0M",&ctx->rpt,dimN,pointsN,data0M);
  70. // read the binary matrix file
  71. unsigned rn,cn,en;
  72. cmRC_t rc0 = cmBinMtxFileSize(ctx,dataFn,&rn,&cn,&en);
  73. double* data1M = cmMemAllocZ(double,rn*cn);
  74. unsigned* colCntV = cmMemAllocZ(unsigned,rn);
  75. cmRC_t rc1 = cmBinMtxFileRead(ctx,dataFn,rn,cn,en,data1M,colCntV);
  76. if( rc0 != cmOkRC || rc1 != cmOkRC )
  77. cmRptPrintf(&ctx->rpt,"Error reading binary matrix file:%s",cmStringNullGuard(dataFn));
  78. // print the binary matrix file
  79. cmVOU_PrintL("colCntV",&ctx->rpt,1,rn,colCntV);
  80. cmVOD_PrintL("data1M", &ctx->rpt,rn,cn,data1M);
  81. cmMemFree(colCntV);
  82. cmMemFree(data0M);
  83. cmMemFree(data1M);
  84. }
  85. void cmFileGetLineTest(cmCtx_t* ctx )
  86. {
  87. const cmChar_t* fn = "/home/kevin/temp/labels.txt";
  88. cmFileH_t f = cmFileNullHandle;
  89. if( cmFileOpen(&f,fn,kReadFileFl,&ctx->rpt) == kOkFileRC )
  90. {
  91. unsigned bufByteCnt = 0;
  92. cmChar_t* buf = NULL;
  93. bool fl = true;
  94. while(fl)
  95. {
  96. switch( cmFileGetLine(f,buf,&bufByteCnt) )
  97. {
  98. case kOkFileRC:
  99. cmRptPrintf(&ctx->rpt,"%i %i %s\n",bufByteCnt,strlen(buf),buf);
  100. break;
  101. case kBufTooSmallFileRC:
  102. buf = cmMemResizeZ(cmChar_t,buf,bufByteCnt);
  103. break;
  104. default:
  105. fl = false;
  106. break;
  107. }
  108. }
  109. cmMemFree(buf);
  110. cmFileClose(&f);
  111. }
  112. }
  113. // Test the cmPvAnl API.
  114. // See cmPvAnlProcTest.m for the equivalent octave code.
  115. void cmPvAnlTest(cmCtx* c )
  116. {
  117. const char* ifn = "/home/kevin/temp/onset0.wav";
  118. const char* ofn = "/home/kevin/temp/test0.aif";
  119. unsigned afChIdx = 0;
  120. unsigned afChCnt = 1;
  121. unsigned afBegSmpIdx = 1000000;
  122. unsigned afEndSmpIdx = 1010000;
  123. unsigned wndSmpCnt = 4096;
  124. unsigned hopSmpCnt = 1024;
  125. unsigned procSmpCnt = hopSmpCnt;
  126. cmAudioFileRd* afRd = cmAudioFileRdAlloc(c,NULL,procSmpCnt,ifn, afChIdx, afBegSmpIdx, afEndSmpIdx );
  127. assert(afRd != NULL );
  128. cmAudioFileWr* afWr = cmAudioFileWrAlloc(c,NULL,procSmpCnt,ofn, afRd->info.srate,afChCnt,afRd->info.bits );
  129. assert(afWr != NULL );
  130. cmPvAnl* pvAnl = cmPvAnlAlloc(c,NULL,procSmpCnt,afRd->info.srate,wndSmpCnt,hopSmpCnt, kNoCalcHzPvaFl );
  131. assert(pvAnl != NULL);
  132. while( cmAudioFileRdRead(afRd) != cmEofRC )
  133. {
  134. while( cmPvAnlExec(pvAnl,afRd->outV, afRd->outN ) )
  135. {
  136. printf(" : %f %i\n",cmVOR_Sum(pvAnl->magV,pvAnl->binCnt),pvAnl->binCnt);
  137. }
  138. cmAudioFileWrExec(afWr,0,afRd->outV,afRd->outN);
  139. }
  140. cmPvAnlFree(&pvAnl);
  141. cmAudioFileWrFree(&afWr);
  142. cmAudioFileRdFree(&afRd);
  143. }
  144. void cmAudioFileProcTest(cmCtx_t* ctx)
  145. {
  146. //const cmChar_t* aifFn = "/home/kevin/media/audio/sourcetone/00-11-060-I-Shapeshifter-TranquilVapor.aiff";
  147. //const cmChar_t* wavFn = "/home/kevin/temp/mas/onset_conv/Piano 3_01.aif";
  148. const cmChar_t* wavFn = "/home/kevin/temp/mas/onsets/Piano 3_01.wav";
  149. //const cmChar_t* wavFn = "/home/kevin/temp/onsetsConv0.aif";
  150. //const cmChar_t* wavFn = "/home/kevin/media/audio/20110723-Kriesberg/Audio Files/Piano 3_01.aif";
  151. const cmChar_t* fn = wavFn;
  152. cmAudioFileInfo_t afInfo;
  153. cmRC_t cmRC;
  154. cmAudioFileH_t afH = cmAudioFileNewOpen( fn, &afInfo, &cmRC, &ctx->rpt );
  155. if( cmRC != kOkAfRC )
  156. printf("Unable to open the audio file:%s\n",fn);
  157. else
  158. {
  159. //cmAudioFileReport( afH, &ctx->rpt, 9785046, 100);
  160. //cmAudioFileReport( afH, &ctx->rpt, 15134420, 100); // onset_conv/Piano 3_01.aif
  161. //cmAudioFileReport( afH, &ctx->rpt, 12862654, 100); // Audio Files/Piano 3_01.wav
  162. cmAudioFileReport( afH, &ctx->rpt, 96092658, 100); // onsets/Piano 3_01.wav
  163. cmAudioFileDelete(&afH);
  164. }
  165. }
  166. // This code test the af
  167. void cmAudioFileReadWriteTest(cmCtx_t* ctx )
  168. {
  169. const cmChar_t* inDir = "/home/kevin/src/cm/src/data/audio_file_format_test";
  170. const cmChar_t* outDir = "/home/kevin/temp/af1/out";
  171. unsigned dirCnt = 0;
  172. unsigned i = 0;
  173. if( cmPlotInitialize(NULL) != cmOkRC )
  174. return;
  175. if( cmFsIsDir(outDir) == false )
  176. cmFsMkDir(outDir);
  177. // get the files in the input directory
  178. cmFileSysDirEntry_t* dep = cmFsDirEntries( inDir, kFileFsFl | kFullPathFsFl, &dirCnt );
  179. for(i=0; dep != NULL && i<dirCnt; ++i)
  180. {
  181. cmAudioFileH_t afH;
  182. cmAudioFileInfo_t afInfo;
  183. cmRC_t cmRC;
  184. // open the ith file in the input directory
  185. if( cmAudioFileIsValid( afH = cmAudioFileNewOpen( dep[i].name, &afInfo, &cmRC, &ctx->rpt)) == false)
  186. {
  187. cmRptPrintf(&ctx->rpt,"Audio file open error occurred on %s\n",dep[i].name);
  188. continue;
  189. }
  190. cmFileSysPathPart_t* pp;
  191. unsigned smpCnt = afInfo.frameCnt * afInfo.chCnt; // count of samples to read and write
  192. cmSample_t* buf = cmMemAlloc( cmSample_t, smpCnt ); // allocate the sample buffer
  193. cmSample_t* bp[ afInfo.chCnt];
  194. unsigned actualFrmCnt;
  195. unsigned chIdx = 0;
  196. unsigned j;
  197. // initialize the audio channel buffer
  198. for(j=0; j<afInfo.chCnt; ++j)
  199. bp[j] = buf + (j*afInfo.frameCnt);
  200. // parse the input file name
  201. if((pp = cmFsPathParts(dep[i].name)) == NULL )
  202. cmRptPrintf(&ctx->rpt,"Unable to locate the file parts for '%s'.",dep[i].name);
  203. else
  204. {
  205. // use the input file name to form the output file name
  206. const cmChar_t* outFn = cmFsMakeFn(outDir,pp->fnStr,"bin",NULL);
  207. const cmChar_t* audOutFn = cmFsMakeFn(outDir,pp->fnStr,"aif",NULL);
  208. cmAudioFileH_t aofH = cmAudioFileNewCreate(audOutFn, afInfo.srate, afInfo.bits, afInfo.chCnt, &cmRC, &ctx->rpt);
  209. // read the entire audio file into the sample buffer
  210. if( cmAudioFileReadSample(afH,afInfo.frameCnt,chIdx,afInfo.chCnt,bp, &actualFrmCnt ) != kOkAfRC )
  211. cmRptPrintf(&ctx->rpt,"Audio file read error occurred on %s\n",dep[i].name);
  212. // write the audio file out as a binary matrix file
  213. if( cmBinMtxFileWrite(outFn, actualFrmCnt, afInfo.chCnt, buf, NULL, NULL, &ctx->rpt ) != cmOkRC )
  214. cmRptPrintf(&ctx->rpt,"Binary matrix write failed on '%s'\n", outFn );
  215. // write the audio output file
  216. if( cmAudioFileIsValid(aofH) )
  217. if( cmAudioFileWriteSample(aofH,afInfo.frameCnt,afInfo.chCnt, bp ) != cmOkRC)
  218. cmRptPrintf(&ctx->rpt,"Audio output file write failed on '%s'.", audOutFn);
  219. // plot the audio file signals
  220. cmPlotSetup(pp->fnStr,1,1);
  221. for(j=0; j<afInfo.chCnt; ++j)
  222. cmPlotLineS(NULL,NULL,bp[j],NULL,afInfo.frameCnt,NULL,kSolidPlotLineId);
  223. cmPlotDraw();
  224. //cmKeyPress(NULL);
  225. // close the output audio file
  226. if( cmAudioFileIsValid(aofH) )
  227. cmAudioFileDelete(&aofH);
  228. cmFsFreeFn(audOutFn); // release the output audio file name
  229. cmFsFreeFn(outFn); // release the output file name
  230. cmFsFreePathParts(pp); // release the parse recd
  231. }
  232. cmAudioFileDelete(&afH); // release the audio file
  233. cmMemPtrFree(&buf); // release the audio buffer
  234. }
  235. cmFsDirFreeEntries(dep);
  236. cmPlotFinalize();
  237. }
  238. void cmZeroCrossTest( cmRpt_t* rpt )
  239. {
  240. double srate = 32;
  241. unsigned vn = srate * 2;
  242. cmSample_t v6[ vn ];
  243. unsigned impPhs = 0;
  244. double hz = 2;
  245. cmSample_t d = 0;
  246. impPhs = cmVOS_SynthSine( v6, vn, impPhs, srate, hz );
  247. cmVOS_MultVS( v6, vn, -1 );
  248. unsigned zc = cmVOS_ZeroCrossCount( v6, vn, &d);
  249. cmTestPrint(rpt,"zero cross: %i \n",zc);
  250. cmVOS_Print( rpt, 1, vn, v6 );
  251. }
  252. void cmMelTest(cmRpt_t* rpt)
  253. {
  254. double srate = 44100;
  255. unsigned binCnt = 513;
  256. unsigned bandCnt = 36;
  257. cmSample_t t[ binCnt * bandCnt ];
  258. cmCtx* c = cmCtxAlloc(NULL,rpt,cmLHeapNullHandle,cmSymTblNullHandle);
  259. cmMatrixBuf* m = cmMatrixBufAlloc( c, NULL, binCnt, bandCnt );
  260. cmVOS_MelMask( m->bufPtr, bandCnt, binCnt, srate, kShiftMelFl );
  261. cmVOS_Transpose( t, m->bufPtr, bandCnt, binCnt );
  262. cmPlotSetup("Test Proc Impl",1,1);
  263. unsigned i;
  264. double sum = 0;
  265. for(i=0; i<bandCnt; ++i )
  266. {
  267. sum += cmVOS_Sum( t + (i*binCnt), binCnt );
  268. cmPlotLineS( NULL, NULL, t + (i*binCnt), NULL, 35, NULL, kSolidPlotLineId );
  269. }
  270. printf("sum:%f\n",sum);
  271. cmPlotDraw();
  272. cmKeyPress(NULL);
  273. cmCtxFree(&c);
  274. cmMatrixBufFree( &m);
  275. }
  276. void cmDctTest()
  277. {
  278. unsigned coeffCnt = 20;
  279. unsigned filtCnt = 36;
  280. cmSample_t m[ coeffCnt * filtCnt ];
  281. cmSample_t t[ coeffCnt * filtCnt ];
  282. cmVOS_DctMatrix( m, coeffCnt, filtCnt );
  283. cmVOS_Transpose( t, m, coeffCnt, filtCnt );
  284. cmPlotSetup("Test",1,1);
  285. unsigned i;
  286. for(i=0; i<coeffCnt; ++i )
  287. {
  288. cmPlotLineS( NULL, NULL, t+(i*filtCnt), NULL, filtCnt, NULL, kSolidPlotLineId );
  289. }
  290. cmPlotDraw();
  291. cmKeyPress(NULL);
  292. }
  293. void cmMtxMultTest(cmRpt_t* rpt)
  294. {
  295. unsigned mrn = 3;
  296. unsigned mcn = 2;
  297. unsigned vn = mcn;
  298. unsigned on = mrn;
  299. double m[] = { 1, 2, 3, 4, 5, 6 };
  300. double v[] = { 1, 2 };
  301. double o[ on ];
  302. cmVOD_MultVMV( o, mrn, m, mcn, v );
  303. cmVOD_PrintL( "o:\n", rpt, 1, on, o );
  304. cmVOD_MultVVM( o, on, v, vn, m );
  305. cmVOD_PrintL( "o:\n", rpt, 1, on, o );
  306. cmReal_t A[] = { 8, 3, 4, 1, 5, 9, 6, 7, 2 }; // magic(3)
  307. cmReal_t B[] = { 1, 2, 3, 4, 5, 6};
  308. unsigned arn = 3;
  309. unsigned cmn = 3;
  310. unsigned brn = cmn;
  311. unsigned bcn = 2;
  312. cmReal_t D[arn*bcn];
  313. cmVOR_MultMMM( D, arn, bcn, A, B, brn );
  314. cmVOR_PrintL("D:\n",rpt,arn,bcn,D);
  315. // D = B*At
  316. cmVOR_MultMMMt( D, bcn, arn, B, A, cmn );
  317. cmVOR_PrintL("Dt:\n",rpt,bcn,arn,D);
  318. // D += B*At
  319. cmVOR_MultMMM1( D, bcn, arn, 1.0, B, A, cmn, 1.0, kTransposeM1Fl );
  320. cmVOR_PrintL("Dt:\n",rpt,bcn,arn,D);
  321. // D = B*At - with explicit physical row counts for each matrix
  322. cmVOR_MultMMM2( D, bcn, arn, 1.0, B, A, cmn, 0.0, kTransposeM1Fl, bcn, bcn, cmn );
  323. cmVOR_PrintL("Dt:\n",rpt,bcn,arn,D);
  324. // D = 3*B*At - with explicit physical row counts for each matrix
  325. cmVOR_MultMMM2( D, bcn, arn, 3.0, B, A, cmn, 0.0, kTransposeM1Fl, bcn, bcn, cmn );
  326. cmVOR_PrintL("Dt:\n",rpt,bcn,arn,D);
  327. /// dpb[dn] = mp[mrn,dn] * vp[mrn]
  328. double v1[] = {1, 2, 3};
  329. cmVOD_MultVMtV( o, mcn, m, mrn, v1 );
  330. cmVOD_PrintL( "o:\n", rpt, 1, on, o );
  331. cmTestPrint(rpt,"multsum: %f\n", cmVOR_MultSumVV( A, B, 6 ));
  332. }
  333. void cmShiftRotateTest(cmRpt_t* rpt)
  334. {
  335. unsigned vn = 5;
  336. double v[] = { 0, 1, 2, 3, 4 };
  337. cmVOD_Print( rpt, 1, vn, v );
  338. cmVOD_Rotate( v, vn, 2 );
  339. cmVOD_Print( rpt, 1, vn, v );
  340. cmVOD_Rotate( v, vn, -3 );
  341. cmVOD_Print( rpt, 1, vn, v );
  342. cmVOD_Shift( v, vn, 2, 9 );
  343. cmVOD_Print( rpt, 1, vn, v );
  344. cmVOD_Shift( v, vn, -2, 9 );
  345. cmVOD_Print( rpt, 1, vn, v );
  346. unsigned pn = 6;
  347. double p[] = { 0, 1, 0, 1, 0, 1 };
  348. unsigned pi[5];
  349. unsigned pm = cmVOD_PeakIndexes( pi, 5, p, pn, 1.1 );
  350. printf("%i : ",pm);
  351. cmVOU_Print( rpt, 1, pm, pi );
  352. cmReal_t m0[] = { 1, 2, 3, 4, 5, 6, 7, 8, 9 };
  353. cmReal_t m1[9];
  354. cmVOD_PrintL("m0:\n",rpt,3,3,m0);
  355. cmVOR_RotateM( m1, 3,3, m0, 3, -2 );
  356. cmVOD_PrintL("m1:\n",rpt,3,3,m1);
  357. }
  358. void cmRandomTest(cmRpt_t* rpt )
  359. {
  360. cmPlotSetup("Random",1,1);
  361. unsigned i;
  362. unsigned yn = 1000;
  363. unsigned y[yn];
  364. cmReal_t w[] = { .2, .3, .5};
  365. // seed the random number generator with a clock with seconds resolution
  366. srand((unsigned)time(NULL));
  367. unsigned wn = sizeof(w)/sizeof(w[0]);
  368. unsigned h[wn];
  369. cmVOR_WeightedRandInt(y,yn,w,wn);
  370. cmVOU_Hist(h,wn,y,yn);
  371. cmVOU_Print(rpt,1,wn,h);
  372. for(i=0; i<wn; ++i)
  373. printf("%f ",(float)h[i]/yn);
  374. printf("\n");
  375. }
  376. void cmTestGaussWin(cmRpt_t* rpt)
  377. {
  378. unsigned vn = 10;
  379. double v[ vn ];
  380. unsigned n = 4;
  381. cmVOD_GaussWin(v,n,5.0/3.0);
  382. cmVOD_Print(rpt,1,n,v);
  383. }
  384. // meanV[K] and varV[K] identify the center and variance of K clusters
  385. // kV[xcn] identifes the cluster associated with each data point
  386. // Each column of xM[] contains a single 2D data point (xrn is therefore always 2)
  387. void cmPlotGauss2D( cmRpt_t* rpt, const cmReal_t* meanV, const cmReal_t* varV, const cmReal_t* xM, unsigned xD, unsigned xN, const char* colorStr )
  388. {
  389. cmReal_t txM[ xD*xN ];
  390. cmVOR_Transpose(txM,xM,xD,xN);
  391. cmPlotLineR( NULL, txM, txM + xN, NULL, xN, colorStr, kCirclePlotPtId );
  392. unsigned cn = 100;
  393. cmReal_t cV[ cn*2 ];
  394. cmVOR_CircleCoords( cV, cn, meanV[0], meanV[1], varV[0], varV[1] );
  395. cmPlotLineR( NULL, cV, cV+cn, NULL, cn, colorStr, kSolidPlotLineId );
  396. cmPlotLineR( NULL, meanV, meanV+1, NULL, 1, colorStr, kAsteriskPlotPtId );
  397. printf("%s\n",colorStr);
  398. cmVOR_Print(rpt,xN,xD,txM);
  399. }
  400. cmReal_t _cmProcTestKmeansDist( void* userPtr, const cmReal_t* v0, const cmReal_t* v1, unsigned vn )
  401. { return cmVOR_EuclidDistance(v0,v1,vn); }
  402. void cmGaussTest(cmRpt_t* rpt )
  403. {
  404. cmReal_t minV = -1;
  405. cmReal_t incr = .1;
  406. unsigned xn = ((fabs(minV)-minV)/incr)+2;
  407. cmReal_t x0[xn],y0[xn];
  408. cmVOR_Seq(x0,xn,minV,incr);
  409. //cmVOR_Print(rpt,1,xn,x0);
  410. cmVOR_GaussPDF( y0, xn, x0, 0, 1 );
  411. cmPlotSetup("Gauss",1,1);
  412. //cmPlotLineR( NULL, NULL, y0, NULL, xn, NULL, kSolidPlotLineId );
  413. unsigned x1N = 10;
  414. unsigned x1D = 2;
  415. cmReal_t m1V[] = { 0, 0 };
  416. cmReal_t v1V[] = { .5,.5 };
  417. cmReal_t x1[ x1D * x1N ];
  418. cmVOR_RandomGaussM( x1, x1D, x1N, m1V, v1V );
  419. //cmPlotGauss2D(rpt,m1V,v1V,x1,x1D,x1N,"magenta");
  420. unsigned x2N = 5; // data points per cluster
  421. unsigned x2D = 2; // data dimensions
  422. unsigned x2K = 3; // count of clusters
  423. cmReal_t m2M[] = { 0, 0, 1, 1, 2, 2 }; // cluster means
  424. cmReal_t v2M[] = { .5,.5, .5, .5, .5, .5 }; // cluster variances
  425. cmReal_t x2[ x2D * x2N * x2K ]; // data matrix
  426. char x2c[][6] = { "red","green","blue" };
  427. cmVOR_RandomGaussMM( x2, x2D, x2N*x2K, m2M, v2M, x2K );
  428. cmVOR_Print(rpt,x2D,x2N*x2K,x2);
  429. unsigned k = 0;
  430. //for(k=0; k<x2K; ++k)
  431. // cmPlotGauss2D( rpt, m2M + (k*x2D), v2M + (k*x2D), x2 + (k*x2D*x2N), x2D, x2N, x2c[k] );
  432. unsigned classIdxV[ x2N*x2K ]; // kmeans generated point assignment vector
  433. cmReal_t centroidM[ x2D*x2K ]; // kmeans generated centroids
  434. // classify each point using kmeans
  435. unsigned iterCnt = cmVOR_Kmeans(classIdxV, centroidM, x2K, x2, x2D, x2N*x2K, NULL, 0, false, _cmProcTestKmeansDist, NULL );
  436. cmTestPrint(rpt,"kmeans iterations:%i\n",iterCnt);
  437. // plot kmeans clusters
  438. for(k=0; k<x2K; ++k)
  439. {
  440. cmReal_t pM[ x2D * x2N ];
  441. cmReal_t dV[ x2D ];
  442. cmReal_t vV[ x2D ];
  443. unsigned i,j=0;
  444. cmVOR_Fill(vV,x2D,0);
  445. // for each data point in cluster k
  446. for(i=0; i<x2N*x2K; ++i)
  447. if( classIdxV[i] == k )
  448. {
  449. // store the data point for later plotting
  450. cmVOR_Copy( pM+(j++*x2D), x2D, x2 + (i*x2D) );
  451. // calculate the variance
  452. cmVOR_SubVVV( dV, x2D, centroidM + (k*x2D), x2 + (i*x2D ));
  453. cmVOR_MultVV( dV, x2D, dV);
  454. cmVOR_AddVV( vV, x2D, dV );
  455. }
  456. // normalize the variance
  457. if( j - 1 > 0 )
  458. cmVOR_DivVS( vV, x2D, j-1 );
  459. // plot data points in cluster k
  460. cmPlotGauss2D( rpt, centroidM + (k*x2D), vV, pM, x2D, j, x2c[k] );
  461. }
  462. cmPlotDraw();
  463. cmKeyPress(NULL);
  464. }
  465. void cmFPExceptTest()
  466. {
  467. double n = DBL_EPSILON;
  468. double d = DBL_MAX;
  469. //double r0 = sqrt(-1.0);
  470. double r0 = n/d;
  471. printf("%e\n",r0);
  472. }
  473. void cmLaTest(cmRpt_t* rpt)
  474. {
  475. unsigned a0n = 3;
  476. unsigned a1n = 4;
  477. cmReal_t A0[] = { 1,4,7,2,5,8,3,6,9 };
  478. cmReal_t A1[] = { 1,0,0,0, 0,-2,0,0, 0,0,3,0, 0,0,0,4 };
  479. cmReal_t t[ a1n*a1n ];
  480. cmReal_t det = cmVOR_DetM(A0,a0n);
  481. printf("det:%e\n",det);
  482. cmVOR_InvM(A0,a0n);
  483. cmVOR_PrintE(rpt,a0n,a0n,A0);
  484. det = cmVOR_DetDiagM(A1,a1n);
  485. printf("det:%e\n",det);
  486. cmVOR_InvDiagM(A1,a1n);
  487. cmVOR_PrintE(rpt,a1n,a1n,A1);
  488. cmReal_t A[] = { 8, 3, 4, 1, 5, 9, 6, 7, 2 }; // magic(3)
  489. cmReal_t B[] = { 1, 2, 3, 4, 5, 6};
  490. //cmVOR_SolveLS(A,3,B,2);
  491. cmVOR_Print(rpt,3,2,B);
  492. cmReal_t D[2*3];
  493. cmVOR_MultMMM( D, 2, 3, B, A, 3 );
  494. cmVOR_Print(rpt,2,3,D);
  495. cmVOR_RandSymPosDef( A1, a1n, t );
  496. cmVOR_PrintL("A1:\n",rpt,a1n,a1n,A1);
  497. cmVOR_Chol(A1,a1n);
  498. cmVOR_PrintL("A1:\n",rpt,a1n,a1n,A1);
  499. cmVOR_RandSymPosDef( A1, a1n, t );
  500. cmVOR_PrintL("A1:\n",rpt,a1n,a1n,A1);
  501. cmVOR_CholZ(A1,a1n);
  502. cmVOR_PrintL("A1:\n",rpt,a1n,a1n,A1);
  503. }
  504. typedef struct
  505. {
  506. const cmReal_t* xM; // matrix base
  507. unsigned D; // row cnt
  508. unsigned N; // col cnt
  509. } cmTestReadFuncData_t;
  510. const cmReal_t* cmTestReadFunc( void* userPtr, unsigned colIdx )
  511. {
  512. cmTestReadFuncData_t* p = (cmTestReadFuncData_t*)userPtr;
  513. assert( userPtr != NULL && colIdx < p->N );
  514. return p->xM + (colIdx*p->D);
  515. }
  516. void cmMvnProbTest(cmRpt_t* rpt)
  517. {
  518. unsigned D = 2;
  519. unsigned N = 3;
  520. cmReal_t uV[] = { 0, 0 };
  521. cmReal_t sM[] = { 1, 0, 0, 1 };
  522. cmReal_t xM[] = { 1, 2, 3, 4, 5, 6 };
  523. cmReal_t yV[ N ];
  524. cmTestReadFuncData_t r;
  525. cmVOR_MultVarGaussPDF( yV, xM, uV, sM, D, N, false );
  526. //cmVOR_PrintE(rpt,1,N,yV);
  527. cmVOR_RandSymPosDef(sM, D, NULL );
  528. //cmVOR_PrintL("sM:\n",rpt,D,D,sM);
  529. cmVOR_RandomGaussNonDiagM( xM, D, N, uV, sM, NULL );
  530. //cmVOR_PrintL("yM:\n",rpt,D,N,xM);
  531. cmReal_t S[] ={ 0.67462, 0.49828, 0.49828, 0.36804 };
  532. cmVOR_MultVS(S,N*N,10);
  533. cmReal_t mu[] = {0, 0};
  534. cmReal_t logDet = cmVOR_LogDetM(S,D);
  535. cmReal_t x[] = {-.1, -.1, 0, 0, .1, .1};
  536. bool diagFl = false;
  537. cmVOR_InvM(S,D);
  538. r.xM = x;
  539. r.D = D;
  540. r.N = N;
  541. cmVOR_MultVarGaussPDF3( yV, cmTestReadFunc, &r, mu, S, logDet, D, N, diagFl );
  542. cmVOR_PrintL("pr:\n",rpt,1,N,yV);
  543. cmVOR_MultVarGaussPDF2( yV, x, mu, S, logDet, D, N, diagFl );
  544. cmVOR_PrintL("pr:\n",rpt,1,N,yV);
  545. }
  546. void cmCovarTest(cmRpt_t* rpt)
  547. {
  548. unsigned D = 2;
  549. unsigned N = 1000;
  550. cmReal_t xM[D*N];
  551. cmReal_t uV[D];
  552. cmReal_t sM[D*D];
  553. cmReal_t scM[D*D];
  554. srand((unsigned)time(NULL));
  555. if(1)
  556. {
  557. cmVOR_RandSymPosDef(sM,D,NULL);
  558. cmVOR_Random(uV,D,0,1);
  559. }
  560. if(0)
  561. {
  562. sM[0] = .1;
  563. sM[1] = 0;
  564. sM[2] = 0;
  565. sM[3] = .1;
  566. }
  567. if(0)
  568. {
  569. uV[0] = 0;
  570. uV[1] = 0;
  571. sM[0] = .48533;
  572. sM[1] = .27140;
  573. sM[2] = .27140;
  574. sM[3] = .15191;
  575. }
  576. cmVOR_RandomGaussNonDiagM(xM,D,N,uV,sM,NULL);
  577. cmVOR_GaussCovariance(scM,D,xM,N,NULL,NULL,0);
  578. //cmVOR_PrintL("xM:\n", rpt, D,N,xM);
  579. cmVOR_PrintL("uV: ", rpt, 1, D, uV);
  580. cmVOR_PrintL("sM:\n",rpt, D, D, sM);
  581. cmVOR_PrintL("covar:\n",rpt, D,D,scM);
  582. }
  583. const cmReal_t* cmCovarSrcFunc( void* p, unsigned idx )
  584. { return ((const cmReal_t*)p) + 3*idx; }
  585. void cmCovarTest2(cmRpt_t* rpt )
  586. {
  587. const int D = 3;
  588. const int N = 10;
  589. // each data point is in a column of xM[]
  590. cmReal_t xM[] = {0.18621, 0.39466, 0.29122, 0.49663, 0.58397, 0.98434, 0.26542, 0.88850, 0.10009, 0.18815, 0.42153, 0.30218, 0.56357, 0.55696, 0.50647, 0.64502, 0.78920, 0.70395, 0.88892, 0.26669, 0.27277, 0.74299, 0.32620, 0.89648, 0.99930, 0.78351, 0.35355, 0.86343, 0.87964, 0.21095};
  591. cmReal_t sM[ D * D ];
  592. cmVOR_GaussCovariance2(sM,D,cmCovarSrcFunc,N,xM,NULL,NULL,0);
  593. cmVOR_PrintL("1 covar: ",rpt, D,D,sM);
  594. cmVOR_GaussCovariance(sM,D,xM,N,NULL,NULL,0);
  595. cmVOR_PrintL("2 covar: ",rpt, D,D,sM);
  596. /*
  597. m = [
  598. 0.18621 0.39466 0.29122;
  599. 0.49663 0.58397 0.98434;
  600. 0.26542 0.88850 0.10009;
  601. 0.18815 0.42153 0.30218;
  602. 0.56357 0.55696 0.50647;
  603. 0.64502 0.78920 0.70395;
  604. 0.88892 0.26669 0.27277;
  605. 0.74299 0.32620 0.89648;
  606. 0.99930 0.78351 0.35355;
  607. 0.86343 0.87964 0.21095;
  608. ]
  609. octave> cov(m)
  610. ans =
  611. 0.0885575 0.0092695 0.0123219
  612. 0.0092695 0.0546553 -0.0168115
  613. 0.0123219 -0.0168115 0.0909342
  614. */
  615. }
  616. void cmMahalanobisTest(cmRpt_t* rpt )
  617. {
  618. const int D = 3;
  619. const int N = 10;
  620. // each data point is in a column of xM[]
  621. cmReal_t xM[] = {0.18621, 0.39466, 0.29122, 0.49663, 0.58397, 0.98434, 0.26542, 0.88850, 0.10009, 0.18815, 0.42153, 0.30218, 0.56357, 0.55696, 0.50647, 0.64502, 0.78920, 0.70395, 0.88892, 0.26669, 0.27277, 0.74299, 0.32620, 0.89648, 0.99930, 0.78351, 0.35355, 0.86343, 0.87964, 0.21095};
  622. cmReal_t uV[D];
  623. //cmReal_t xV[] = {0.633826, 0.349463, 0.053582 };
  624. cmReal_t xV[] = {0.72477, 0.98973, 0.11622};
  625. cmReal_t sM[ D * D ];
  626. // find the mean of the data set (mean across the columns)
  627. cmVOR_Mean2( uV, cmCovarSrcFunc, D, N, xM );
  628. cmVOR_PrintL("mean: ",rpt, 1, D, uV );
  629. cmVOR_GaussCovariance(sM,D,xM,N,uV,NULL,0);
  630. cmVOR_PrintL("covar: ",rpt, D,D,sM);
  631. cmReal_t* r = cmVOR_InvM(sM,D);
  632. cmVOR_PrintL("inv covar: ",rpt, D,D,sM);
  633. cmReal_t d = cmVOR_MahalanobisDistance( xV, D, uV, sM );
  634. cmRptPrintf(rpt,"Mahalanobis dist:%f %p\n",d,r);
  635. /*
  636. octave>m =
  637. 0.18621 0.39466 0.29122
  638. 0.49663 0.58397 0.98434
  639. 0.26542 0.88850 0.10009
  640. 0.18815 0.42153 0.30218
  641. 0.56357 0.55696 0.50647
  642. 0.64502 0.78920 0.70395
  643. 0.88892 0.26669 0.27277
  644. 0.74299 0.32620 0.89648
  645. 0.99930 0.78351 0.35355
  646. 0.86343 0.87964 0.21095
  647. octave> u = mean(m)
  648. u = 0.58396 0.58908 0.46220
  649. octave> y
  650. y = 0.633826 0.349463 0.053582
  651. octave> sqrt((y-u)*inv(cov(m))*(y-u)')
  652. ans = 2.0322
  653. */
  654. }
  655. void cmRandIntSeqTest(cmRpt_t* rpt)
  656. {
  657. unsigned vn = 10;
  658. unsigned v[vn];
  659. unsigned i = 0;
  660. for(i=0; i<10; ++i)
  661. {
  662. cmVOU_RandomSeq(v,vn);
  663. cmVOU_PrintL("v: ", rpt, 1, vn, v );
  664. }
  665. }
  666. //------------------------------------------------------------------------------------------------------------
  667. void cmSynthTest()
  668. {
  669. unsigned vn = 128;
  670. unsigned blkN = 2;
  671. cmSample_t v0[ vn*blkN ];
  672. cmSample_t v1[ vn*blkN ];
  673. cmSample_t v2[ vn*blkN ];
  674. cmSample_t v3[ vn*blkN ];
  675. cmSample_t v4[ vn*blkN ];
  676. cmSample_t v5[ vn*blkN ];
  677. cmSample_t v6[ vn*blkN ];
  678. cmSample_t v7[ vn*blkN ];
  679. double srate = vn;
  680. double hz = 1;
  681. unsigned sinPhs = 0,cosPhs = 0, sqrPhs=0, sawPhs=0, triPhs=0, pulPhs=0, impPhs=0, phsPhs=0;
  682. unsigned otCnt = 7;
  683. unsigned i;
  684. for(i=0; i<blkN; ++i)
  685. {
  686. sinPhs = cmVOS_SynthSine( v0+(i*vn), ((i+1)*vn), sinPhs, srate, hz );
  687. cosPhs = cmVOS_SynthCosine( v1+(i*vn), ((i+1)*vn), cosPhs, srate, hz );
  688. sqrPhs = cmVOS_SynthSquare( v2+(i*vn), ((i+1)*vn), sqrPhs, srate, hz, otCnt );
  689. sawPhs = cmVOS_SynthSawtooth( v3+(i*vn), ((i+1)*vn), sawPhs, srate, hz, otCnt );
  690. triPhs = cmVOS_SynthTriangle( v4+(i*vn), ((i+1)*vn), triPhs, srate, hz, otCnt );
  691. pulPhs = cmVOS_SynthPulseCos( v5+(i*vn), ((i+1)*vn), pulPhs, srate, hz, otCnt );
  692. impPhs = cmVOS_SynthImpulse( v6+(i*vn), ((i+1)*vn), impPhs, srate, hz );
  693. phsPhs = cmVOS_SynthPhasor( v7+(i*vn), ((i+1)*vn), phsPhs, srate, hz );
  694. }
  695. cmPlotSetup("Test Proc Impl",2,1);
  696. cmPlotLineS( "cos", NULL, v1, NULL, vn*blkN, NULL, kSolidPlotLineId );
  697. cmPlotLineS( "sqr", NULL, v2, NULL, vn*blkN, NULL, kSolidPlotLineId );
  698. cmPlotLineS( "imp", NULL, v6, NULL, vn*blkN, NULL, kSolidPlotLineId );
  699. cmPlotLineS( "tri", NULL, v4, NULL, vn*blkN, NULL, kSolidPlotLineId );
  700. cmPlotSelectSubPlot( 1, 0 );
  701. cmPlotLineS( "sin", NULL, v0, NULL, vn*blkN, NULL, kSolidPlotLineId );
  702. cmPlotLineS( "saw", NULL, v3, NULL, vn*blkN, NULL, kSolidPlotLineId );
  703. cmPlotLineS( "pul", NULL, v5, NULL, vn*blkN, NULL, kSolidPlotLineId );
  704. cmPlotLineS( "phs", NULL, v7, NULL, vn*blkN, NULL, kSolidPlotLineId );
  705. cmPlotDraw();
  706. //cmPlotPrint(false);
  707. }
  708. void cmMeanVarTest(cmRpt_t* rpt )
  709. {
  710. enum { cnt = 7, dim=2 };
  711. cmReal_t v[cnt*dim] = {0, 1, 3, 4, 6, 7, 9, 10, 12, 13, 15, 16, 18, 19 };
  712. cmReal_t mean;
  713. cmTestPrint(rpt, "sum:%f\n", cmVOR_SumN(v,cnt, dim));
  714. cmTestPrint(rpt, "mean:%f\n", mean = cmVOR_MeanN(v,cnt,dim));
  715. cmTestPrint(rpt, "var:%f\n", cmVOR_VarianceN(v,cnt,dim,NULL));
  716. cmTestPrint(rpt, "var:%f\n", cmVOR_VarianceN(v,cnt,dim,&mean));
  717. unsigned rn = 3;
  718. unsigned cn = 5;
  719. cmReal_t mM[] = {0,1,2,3,4,5,6,7,8,9,10,11,12,13,14};
  720. cmReal_t aV[cn];
  721. cmReal_t vV[cn];
  722. cmVOR_MeanM(aV,mM,rn,cn,0);
  723. cmVOR_PrintL("mean cols: ", rpt,1, cn, aV );
  724. cmVOR_VarianceM(vV,mM,rn,cn,aV,0);
  725. cmVOR_PrintL("var cols: ", rpt,1, cn, vV );
  726. cmVOR_MeanM(aV,mM,rn,cn,1);
  727. cmVOR_PrintL("mean rows: ", rpt,1, rn, aV );
  728. cmVOR_VarianceM(vV,mM,rn,cn,aV,1);
  729. cmVOR_PrintL("var rows: ", rpt,1, rn, vV );
  730. cmVOR_VarianceM(vV,mM,rn,cn,NULL,0);
  731. cmVOR_PrintL("var cols: ", rpt,1, cn, vV );
  732. cmVOR_VarianceM(vV,mM,rn,cn,NULL,1);
  733. cmVOR_PrintL("var rows: ", rpt,1, rn, vV );
  734. cmReal_t mV[] = { 1,2,2,2 };
  735. cmTestPrint(rpt,"Mode:%f\n",cmVOR_Mode(mM,rn*cn));
  736. cmTestPrint(rpt,"Mode:%f\n",cmVOR_Mode(mV,5));
  737. }
  738. void cmMedianFilterTest( cmRpt_t* rpt )
  739. {
  740. enum { xn=5 };
  741. unsigned wn;
  742. cmReal_t x[xn] = { 0, 1, 2, 3, 4 };
  743. cmReal_t y[xn];
  744. for( wn=1; wn<=7; ++wn)
  745. {
  746. cmVOR_MedianFilt( x, xn, wn, y, 1 );
  747. cmTestPrint(NULL,"%i : ",wn);
  748. cmVOR_Print(rpt,1,xn,y);
  749. }
  750. }
  751. void cmConstQTest1( cmConstQ* p, cmRpt_t* rpt )
  752. {
  753. const char* fn = "/home/kevin/src/ac/m2i.txt";
  754. cmCtx* c = cmCtxAlloc(NULL,rpt, cmLHeapNullHandle,cmSymTblNullHandle);
  755. cmMatrixBuf* m0p = cmMatrixBufAllocFile( c, NULL, fn );
  756. printf("target mtx:%i x %i\n", m0p->rn, m0p->cn);
  757. printf("sparse mtx:%i x %i\n", p->wndSmpCnt, p->constQBinCnt);
  758. unsigned ri,ci;
  759. for(ri=0; ri<m0p->rn; ++ri)
  760. {
  761. double sum0 = 0;
  762. double sum1 = 0;
  763. double dsum = 0;
  764. double md = 0;
  765. unsigned mi = -1;
  766. for(ci=0; ci<m0p->cn; ++ci)
  767. {
  768. double v0 = cmMatrixBufColPtr( m0p, ci)[ri];
  769. double v1 = cimag( p->skM[ (ci*p->wndSmpCnt)+ri ] );
  770. double d = fabs(v1-v0);
  771. sum0 += v0;
  772. sum1 += v1;
  773. dsum += d;
  774. if( d > md )
  775. {
  776. mi = ci;
  777. md = d;
  778. }
  779. }
  780. printf("%3i (%4i % 9e) % 9e%% s0:% 9e s1:% 9e\n",ri,mi,md,dsum/sum0,sum0,sum1);
  781. }
  782. cmMatrixBufFree(&m0p);
  783. cmCtxFree(&c);
  784. }
  785. void cmSRCTest( cmRpt_t* rpt )
  786. {
  787. const char* ifn = "/home/kevin/src/ac/test0.aif";
  788. //ifn = "/home/kevin/src/st/fv/audio/00-11-060-I-Shapeshifter-TranquilVapor.aiff";
  789. const char* ofn = "/home/kevin/src/ac/temp1.aif";
  790. unsigned upFact = 1;
  791. unsigned dnFact = 8;
  792. unsigned procSmpCnt = 64;
  793. unsigned iChIdx = 0;
  794. unsigned dnProcSmpCnt = procSmpCnt / dnFact;
  795. unsigned oChCnt = 1;
  796. unsigned oBitsPerSmp = 16;
  797. unsigned oChIdx = 0;
  798. cmCtx* c = cmCtxAlloc(NULL,rpt,cmLHeapNullHandle,cmSymTblNullHandle);
  799. cmAudioFileRd* arp = cmAudioFileRdAlloc( c, NULL, procSmpCnt, ifn, iChIdx, 0, cmInvalidIdx );
  800. if( arp != NULL )
  801. {
  802. cmSRC* srp = cmSRCAlloc( c, NULL, arp->info.srate,procSmpCnt, upFact, dnFact );
  803. cmSRC* urp = cmSRCAlloc( c, NULL, arp->info.srate/dnFact,dnProcSmpCnt, dnFact, upFact );
  804. cmAudioFileWr* awp = cmAudioFileWrAlloc( c, NULL, procSmpCnt, ofn, arp->info.srate, oChCnt, oBitsPerSmp);
  805. //cmSRCShow(c,srp);
  806. printf("frame cnt:%i\n",arp->info.frameCnt);
  807. unsigned i = 0, j=0, cnt=10000;
  808. while( cmAudioFileRdRead( arp ) == cmOkRC )
  809. {
  810. if( arp->eofFl || arp->outN != procSmpCnt )
  811. break;
  812. cmSRCExec(srp,arp->outV, arp->outN);
  813. cmSRCExec(urp,srp->outV, srp->outN);
  814. cmAudioFileWrExec( awp, oChIdx, urp->outV, urp->outN );
  815. i+=arp->outN;
  816. if( i > procSmpCnt * cnt )
  817. {
  818. j += i;
  819. i = 0;
  820. printf("%i ",j);
  821. fflush(stdout);
  822. }
  823. }
  824. printf("%i done\n",j + i);
  825. fflush(stdout);
  826. cmAudioFileWrFree(&awp);
  827. cmSRCFree(&urp);
  828. cmSRCFree(&srp);
  829. }
  830. cmAudioFileRdFree(&arp);
  831. cmCtxFree(&c);
  832. }
  833. void cmBeatHistTest( cmRpt_t* rpt )
  834. {
  835. unsigned i;
  836. cmCtx c;
  837. cmCtxInit(&c, rpt,cmLHeapNullHandle,cmSymTblNullHandle);
  838. cmPlotSetup("Beat Histogram Test",1,1);
  839. enum { wndN = 13 };
  840. unsigned frmCnt = 512; // df[] element count
  841. double srate = 1/.0116; // df[] sample rate
  842. double bpm = 120; // beats per minute
  843. //unsigned spb = floor(60*srate/bpm); // samples per beat
  844. cmSample_t wndV[ wndN ];
  845. cmSample_t df[ frmCnt ];
  846. cmSample_t is[ frmCnt ];
  847. // create a df signal by convolving a impulse train with a hanning window
  848. cmVOS_HannMatlab(wndV,wndN);
  849. cmVOS_SynthImpulse( is, frmCnt, 0, srate, bpm/60 );
  850. cmConvolveSignal( &c, wndV, wndN, is, frmCnt, df, frmCnt );
  851. //cmPlotLineS( "df", NULL, df, NULL, frmCnt, NULL, kSolidPlotLineId );
  852. //cmVOS_Print(rpt, 1, 25, df );
  853. cmBeatHist* p = cmBeatHistAlloc(&c,NULL,frmCnt);
  854. //cmVOR_Print(rpt, p->frmCnt, p->maxLagCnt, p->m );
  855. for(i=0; i<frmCnt; ++i)
  856. cmBeatHistExec(p,df[i]);
  857. cmBeatHistCalc(p);
  858. cmPlotDraw();
  859. cmBeatHistFree(&p);
  860. }
  861. void cmOlaProcTest( cmRpt_t* rpt )
  862. {
  863. cmCtx c;
  864. cmOla ola;
  865. unsigned wndSmpCnt = 32;
  866. unsigned hopSmpCnt = 8;
  867. unsigned procSmpCnt = 4;
  868. unsigned wndCnt = 8;
  869. unsigned i = 0;
  870. cmSample_t wndV[] = { 0,1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,15,14,13,12,11,10,9,8,7,6,5,4,3,2,1 };
  871. assert( sizeof(wndV)/sizeof(wndV[0]) == wndSmpCnt );
  872. cmCtxAlloc(&c, rpt,cmLHeapNullHandle,cmSymTblNullHandle);
  873. cmOlaAlloc(&c,&ola,wndSmpCnt,hopSmpCnt, procSmpCnt, kUnityWndId );
  874. for(i=0; i<wndCnt; ++i)
  875. {
  876. cmOlaExecS(&ola,wndV,wndSmpCnt);
  877. int j;
  878. for(j=0; j<hopSmpCnt/procSmpCnt; ++j)
  879. {
  880. cmVOS_Print( rpt, 1, ola.procSmpCnt, ola.outPtr );
  881. cmOlaExecOut(&ola);
  882. }
  883. }
  884. cmObjFreeStatic( cmOlaFree, cmOla, ola );
  885. cmObjFreeStatic( cmCtxFree, cmCtx, c );
  886. }
  887. void cmTestBarkFiltMask(cmRpt_t* rpt)
  888. {
  889. // 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23
  890. //cmReal_t bl[]= { 0,100,200,300,400,510,630,770, 920,1080,1270,1480,1720,2000,2320,2700,3150,3700,4400,5300,6400,7700, 9500,12000};
  891. //cmReal_t bc[]= { 50,150,250,350,450,570,700,840,1000,1170,1370,1600,1850,2150,2500,2900,3400,4000,4800,5800,7000,8500,10500,13500};
  892. //cmReal_t bh[]= {100,200,300,400,510,630,770,920,1080,1270,1480,1720,2000,2320,2700,3150,3700,4400,5300,6400,7700,9500,12000,15500};
  893. // -1 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 (23+1)
  894. cmReal_t b[]= {0, 50,150,250,350,450,570,700,840,1000,1170,1370,1600,1850,2150,2500,2900,3400,4000,4800,5800,7000,8500,10500,13500, 15500 };
  895. unsigned bandCnt = 24;
  896. unsigned binCnt = 64;
  897. unsigned binHz = 100.0;
  898. cmReal_t maskMtx[ bandCnt * binCnt ];
  899. cmReal_t stSpread = 0;
  900. cmVOR_TriangleMask(maskMtx, bandCnt, binCnt, b+1, binHz, stSpread, b+0, b+2 );
  901. cmVOR_Print( rpt, bandCnt, binCnt, maskMtx);
  902. }
  903. void cmTestFftTemplate(cmRpt_t* rpt)
  904. {
  905. double srate = 16;
  906. double hz = 1;
  907. unsigned n = 16;
  908. unsigned otCnt = 3;
  909. cmCtx* c = cmCtxAlloc( NULL, rpt,cmLHeapNullHandle,cmSymTblNullHandle);
  910. // take the fft of a vector of cmSample_t values
  911. cmFftSR* p0;
  912. cmIFftRS* p2;
  913. cmSample_t vs[n];
  914. cmVOS_SynthSquare( vs, n, 0, srate, hz, otCnt );
  915. cmVOS_PrintL( "\nsig:\n",rpt, 1, n, vs);
  916. p0 = cmFftAllocSR(c,NULL, vs, n, kToPolarFftFl );
  917. p2 = cmIFftAllocRS(c,NULL,p0->binCnt);
  918. cmFftExecSR( p0, NULL, 0 );
  919. cmVOR_DivVS( p0->magV, p0->binCnt, p0->wndSmpCnt );
  920. cmVOR_PrintL( "\nmag:\n",rpt, 1, p0->binCnt, p0->magV );
  921. cmIFftExecPolarRS( p2, p0->magV, p0->phsV );
  922. cmVOS_PrintL( "\nifft:\n",rpt,1, p2->binCnt, p2->outV );
  923. cmIFftFreeRS(&p2);
  924. cmFftFreeSR(&p0);
  925. // take the fft of a vector of cmReal_t values
  926. cmFftRR* p1;
  927. cmIFftRR* p3;
  928. cmReal_t vr[n];
  929. cmVOR_SynthSquare( vr, n, 0, srate, hz, otCnt );
  930. cmVOR_PrintL("\nsig\n", rpt, 1, n, vr);
  931. p1 = cmFftAllocRR(c,NULL, vr, n, kToPolarFftFl );
  932. p3 = cmIFftAllocRR(c,NULL, p1->binCnt);
  933. cmFftExecRR( p1, NULL, 0 );
  934. cmVOR_DivVS( p1->magV, p1->binCnt, p1->wndSmpCnt );
  935. cmVOR_PrintL( "\nmag:\n",rpt, 1, p1->binCnt, p1->magV );
  936. cmIFftExecPolarRR( p3, p1->magV, p1->phsV );
  937. cmVOR_PrintL( "\nifft:\n",rpt,1, p3->binCnt, p3->outV );
  938. cmIFftFreeRR(&p3);
  939. cmFftFreeRR(&p1);
  940. cmCtxFree(&c);
  941. }
  942. cmReal_t cmSelDistFunc( void* userPtr, const cmReal_t* v0, const cmReal_t* v1, unsigned vn )
  943. { return cmVOR_EuclidDistance(v0,v1,vn); }
  944. void cmSelectCols(cmRpt_t* rpt)
  945. {
  946. unsigned mcn = 10;
  947. unsigned mrn = 3;
  948. unsigned selN = 3;
  949. cmReal_t m0[ mrn * mcn ];
  950. cmReal_t m1[ mrn * selN ];
  951. unsigned selIdxV[ selN ];
  952. cmVOR_Random(m0,mrn*mcn,0.0,1.0);
  953. cmVOR_PrintL("\norg\n",rpt, mrn, mcn, m0 );
  954. cmVOR_SelectRandom(m1,selIdxV,selN,m0,mrn,mcn);
  955. cmVOR_PrintL("\nrandom\n",rpt, mrn, selN, m1 );
  956. cmVOU_PrintL("\nindexes\n",rpt, 1, selN, selIdxV );
  957. cmVOR_SelectMaxDist(m1, selIdxV, selN, m0, mrn, mcn, cmSelDistFunc, NULL );
  958. cmVOR_PrintL("\nselect max dist\n",rpt, mrn, selN, m1 );
  959. cmVOU_PrintL("\nindexes\n",rpt, 1, selN, selIdxV );
  960. cmVOR_SelectMaxAvgDist(m1, selIdxV, selN, m0, mrn, mcn, cmSelDistFunc, NULL );
  961. cmVOR_PrintL("\nselect max avg dist\n",rpt, mrn, selN, m1 );
  962. cmVOU_PrintL("\nindexes\n",rpt, 1, selN, selIdxV );
  963. }
  964. void cmNMF_Test( cmRpt_t* rpt )
  965. {
  966. //const char* plotDev = "wxwidgets";
  967. int r = 2;
  968. int n = 100;
  969. int m = 1000;
  970. int vn = 25;
  971. unsigned maxIterCnt = 2000;
  972. unsigned convergeCnt = 40;
  973. cmReal_t* V = cmMemAllocZ( cmReal_t, n*m );
  974. cmReal_t* W = cmMemAllocZ( cmReal_t, n*r );
  975. cmReal_t* H = cmMemAllocZ( cmReal_t, r*m );
  976. cmReal_t* d = cmMemAllocZ( cmReal_t, vn );
  977. cmReal_t* t = cmMemAllocZ( cmReal_t, vn );
  978. unsigned i;
  979. cmCtx* ctx = cmCtxAlloc(NULL,rpt,cmLHeapNullHandle,cmSymTblNullHandle);
  980. cmNmf_t* nmf = cmNmfAlloc(ctx, NULL, n, m, r, maxIterCnt, convergeCnt );
  981. cmVOR_Hann(d,vn);
  982. // prevent zeros in V[]
  983. cmVOR_Random(V,n*m,0.0,0.1);
  984. for(i=0; i<m; ++i)
  985. {
  986. cmVOR_MultVVS( t, vn, d, (cmReal_t)i/m );
  987. cmVOR_AddVV( V + (i*n), vn, t );
  988. }
  989. for(i=0; i<m; ++i)
  990. {
  991. cmVOR_MultVVS( t, vn, d, 1.0 - ((cmReal_t)i/m) );
  992. cmVOR_AddVV( V + (i*n) + n - (vn+1), vn, t );
  993. }
  994. cmNmfExec( nmf, V, m );
  995. //cmPlviewPlotImage(plotDev, vRptFunc, NULL, V, n, m, 0, m, 0, n );
  996. //cmPlviewPlotColY(plotDev, vRptFunc, NULL, nmf->W, n, r, 0, n, 0.0, 1.0 );
  997. //cmPlviewPlotRowY(plotDev, vRptFunc, NULL, nmf->H, r, m, 0, m, 0.0, 1.0 );
  998. cmNmfFree(&nmf);
  999. cmCtxFree(&ctx);
  1000. cmMemPtrFree(&V);
  1001. cmMemPtrFree(&W);
  1002. cmMemPtrFree(&H);
  1003. cmMemPtrFree(&d);
  1004. cmMemPtrFree(&t);
  1005. }
  1006. void cmLinearMapTest( cmRpt_t* rpt )
  1007. {
  1008. int sN = 4;
  1009. int dN = 15;
  1010. cmReal_t sV[sN];
  1011. cmVOR_Seq( sV, sN, 0, 1 );
  1012. cmReal_t dV[dN];
  1013. cmVOR_LinearMap(dV, dN, sV, sN );
  1014. cmVOR_Print(rpt, 1, dN, dV );
  1015. int s2N = 15;
  1016. int d2N = 4;
  1017. cmReal_t s2V[sN];
  1018. cmReal_t d2V[dN];
  1019. cmVOR_Seq( s2V, s2N, 0, 1 );
  1020. cmVOR_LinearMap(d2V, d2N, s2V, s2N );
  1021. cmVOR_Print(rpt, 1, d2N, d2V );
  1022. }
  1023. void cmMemErrCallback( void* userDataPtr, const char* fmt, va_list vl )
  1024. { vprintf(fmt,vl);}
  1025. // Test a cmProc
  1026. void cmProcTestProc(cmCtx_t* ctx )
  1027. {
  1028. unsigned baseSymId = 1000;
  1029. unsigned dfltBlockByteCnt = 1024;
  1030. cmLHeapH_t lhH = cmLHeapCreate(dfltBlockByteCnt,ctx);
  1031. cmSymTblH_t stH = cmSymTblCreate(cmSymTblNullHandle,baseSymId,ctx);
  1032. cmCtx* c = cmCtxAlloc(NULL,&ctx->rpt,lhH,stH);
  1033. assert( cmLHeapIsValid(lhH));
  1034. assert( cmSymTblIsValid(stH));
  1035. //cmShiftBufTest(c);
  1036. cmPvAnlTest(c);
  1037. cmLHeapDestroy(&lhH);
  1038. cmSymTblDestroy(&stH);
  1039. cmCtxFree(&c);
  1040. }
  1041. void cmPuTest(cmCtx_t* ctx);
  1042. void cmAudLabelFileTest(cmCtx_t* ctx);
  1043. void cmProcTestNoInit(cmCtx_t* ctx)
  1044. {
  1045. //cmSelectCols(rpt);
  1046. //cmNMF_Test(rpt);
  1047. //cmGmmTest( rpt );
  1048. //cmMvnProbTest(rpt);
  1049. //cmShiftBufTest(rpt,NULL);
  1050. //cmLinearMapTest(&ctx->rpt);
  1051. //cmCovarTest2(&ctx->rpt);
  1052. //cmMahalanobisTest(&ctx->rpt);
  1053. //cmProcTestProc(ctx); // test a cmProcObj
  1054. //cmAudioFileProcTest(ctx);
  1055. //cmAudioFileReadWriteTest(ctx);
  1056. //cmPuTest(ctx);
  1057. //cmFileGetLineTest(ctx);
  1058. //cmAudLabelFileTest(ctx);
  1059. //cmStackTest(ctx);
  1060. //cmBinMtxFileTest(ctx );
  1061. //cmMtxMultTest(ctx->err.rpt);
  1062. //cmStandardizeTest(ctx);
  1063. cmRbmBinaryTest(ctx);
  1064. }
  1065. void cmProcTestGnuPlot( cmCtx_t* ctx )
  1066. {
  1067. cmPlotInitialize(NULL);
  1068. cmProcTestNoInit(ctx);
  1069. cmTestPrint(&ctx->rpt,"%s\n","press any key");
  1070. cmKeyPress(NULL);
  1071. cmPlotFinalize();
  1072. }
  1073. void cmProcTest(cmCtx_t* ctx)
  1074. {
  1075. cmKbRecd kb;
  1076. cmMdInitialize( ctx->guardByteCnt, ctx->alignByteCnt, ctx->mmFlags, &ctx->rpt );
  1077. cmPlotInitialize(NULL);
  1078. // cmAudioFileTest();
  1079. // this should be added to cmCtxInit()
  1080. //cmSetupFloatPointExceptHandler(NULL);
  1081. //cmDelayTest();
  1082. //cmMelTest();
  1083. //cmDctTest();
  1084. //cmMtxMultTest(rpt);
  1085. //cmSonesTest();
  1086. //cmProcImplTestSynth();
  1087. //cmZeroCrossTest();
  1088. //cmFIRTest();
  1089. //cmFftTest();
  1090. //cmFftTestComplex();
  1091. //cmAudioFileWrTest();
  1092. //cmSRCTest();
  1093. //cmProcImplTest1();
  1094. //cmFuncFilterTest();
  1095. //cmRandomTest(rpt);
  1096. //cmTestGaussWin();
  1097. //cmTestMtxBuf();
  1098. //cmDhmmTest();
  1099. //cmGaussTest(rpt);
  1100. //cmLaTest(rpt);
  1101. //cmFPExceptTest();
  1102. //cmMvnProbTest(rpt);
  1103. //cmShiftRotateTest(rpt);
  1104. //cmMeanVarTest(rpt);
  1105. //cmStatsProcTest(rpt);
  1106. //cmMedianFilterTest(rpt);
  1107. //cmFilterTest(rpt);
  1108. //cmBeatHistTest(rpt);
  1109. //cmIFftTest(rpt);
  1110. //cmConvolveTest(rpt);
  1111. //cmBeatHistTest(rpt);
  1112. //cmWndFuncTest(rpt);a
  1113. //cmGmmTest( rpt );
  1114. //cmChmmTest( rpt );
  1115. //cmCovarTest(rpt);
  1116. //cmRandIntSeqTest(rpt);
  1117. //cmChordTest(rpt);
  1118. //cmConstQTest1()
  1119. //void cmSRCTest( rpt );
  1120. //cmFuncFilterTest();
  1121. //cmOlaProcTest(rpt);
  1122. //cmTestBarkFiltMask(rpt);
  1123. //cmTestFftTemplate(rpt);
  1124. cmProcTestNoInit(ctx);
  1125. cmTestPrint(&ctx->rpt,"%s\n","press any key");
  1126. cmKeyPress(&kb);
  1127. cmPlotFinalize();
  1128. cmMdReport( kIgnoreNormalMmFl );
  1129. cmMdFinalize();
  1130. }