cml/libcm
1
0
Fork 0
Code Issues Pull Requests Releases Wiki Activity

Merge branch 'master' of klarke.webfactional.com:webapps/git/repos/libcm

Browse Source
This commit is contained in:
Kevin Larke 2015-07-24 10:47:01 -07:00
commit 1f74830cb8
15 changed files with 1381 additions and 46 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

2
Makefile.am
View File

@ -4,7 +4,7 @@ cmHDR =
cmSRC =
cmHDR += src/libcm/cmErr.h src/libcm/cmCtx.h src/libcm/cmRpt.h src/libcm/cmGlobal.h src/libcm/cmComplexTypes.h src/libcm/cmFloatTypes.h src/libcm/cmPrefix.h
cmSRC += src/libcm/cmErr.c src/libcm/cmCtx.c src/libcm/cmRpt.c src/libcm/cmGlobal.c
cmSRC += src/libcm/cmErr.c src/libcm/cmCtx.c src/libcm/cmRpt.c src/libcm/cmGlobal.c src/libcm/cmComplexTypes.c
cmHDR += src/libcm/cmSerialize.h src/libcm/cmSymTbl.h src/libcm/cmHashTbl.h src/libcm/cmFileSys.h src/libcm/cmFile.h
cmSRC += src/libcm/cmSerialize.c src/libcm/cmSymTbl.c src/libcm/cmHashTbl.c src/libcm/cmFileSys.c src/libcm/cmFile.c

4
cmAudioPort.c
View File

@ -363,6 +363,8 @@ void cmApReport( cmRpt_t* rpt )
}
}
//cmApAlsaDeviceReport(rpt);
}
/// [cmAudioPortExample]
@ -668,6 +670,7 @@ int cmApPortTest( bool runFl, cmRpt_t* rpt, int argc, const char* argv[] )
runFl = _cmApGetOpt(argc,argv,"-p",!runFl,true)?false:true;
r.srate = _cmApGetOpt(argc,argv,"-r",44100,false);
r.chIdx = _cmApGetOpt(argc,argv,"-a",0,false);
r.chCnt = _cmApGetOpt(argc,argv,"-c",2,false);
r.bufCnt = _cmApGetOpt(argc,argv,"-b",3,false);
@ -685,7 +688,6 @@ int cmApPortTest( bool runFl, cmRpt_t* rpt, int argc, const char* argv[] )
r.outDevIdx = _cmGlobalOutDevIdx = _cmApGetOpt(argc,argv,"-o",2,false);
r.phase = 0;
r.frqHz = 2000;
r.srate = 44100;
r.bufInIdx = 0;
r.bufOutIdx = 0;
r.bufFullCnt = 0;

60
cmComplexTypes.c Normal file
View File

@ -0,0 +1,60 @@
#include "cmPrefix.h"
#include "cmGlobal.h"
#include "cmFloatTypes.h"
#include "cmComplexTypes.h"
void cmVOCR_MultVVV( cmComplexR_t* y, const cmComplexS_t* x0, const cmComplexR_t* x1, unsigned n )
{
unsigned i;
for(i=0; i<n; ++i)
{
y[i] = x0[i] * x1[i];
/*
cmReal_t ab = x0[i].r * x1[i].r;
cmReal_t bd = x0[i].i * x1[i].i;
cmReal_t bc = x0[i].i * x1[i].r;
cmReal_t ad = x0[i].r * x1[i].i;
y[i].r = ab - bd;
y[i].i = bc + ad;
*/
}
}
void cmVOCR_MultVFV( cmComplexR_t* y, const float* x, unsigned n )
{
unsigned i;
for(i=0; i<n; ++i)
{
y[i] *= x[i];
}
}
void cmVOCR_DivVFV( cmComplexR_t* y, const float* x, unsigned n )
{
unsigned i;
for(i=0; i<n; ++i)
{
y[i] /= x[i];
}
}
void cmVOCR_Abs( cmSample_t* y, const cmComplexR_t* x, unsigned n )
{
unsigned i;
for(i=0; i<n; ++i)
y[i] = (cmSample_t)cmCabsR(x[i]);
}
void cmVOCR_DivR_VV( cmComplexR_t* y, const cmReal_t* x, unsigned n )
{
unsigned i;
for(i=0; i<n; ++i)
{
y[i] /= x[i];
//y[i].r /= x[i];
//y[i].i /= x[i];
}
}

12
cmComplexTypes.h
View File

@ -11,6 +11,7 @@
#define cmCrealS crealf
#define cmCimagS cimagf
#define cmCargS cargf
#define cmCconjS conjf
#define cmFftPlanAllocS fftwf_plan_dft_r2c_1d
#define cmFft1dPlanAllocS fftwf_plan_dft_1d
@ -29,6 +30,7 @@
#define cmCrealS creal
#define cmCimagS cimag
#define cmCargS carg
#define cmCconjS conj
#define cmFftPlanAllocS fftw_plan_dft_r2c_1d
#define cmFft1dPlanAllocS fftw_plan_dft_1d
@ -53,6 +55,7 @@
#define cmCrealR crealf
#define cmCimagR cimagf
#define cmCargR cargf
#define cmCconjR conjf
#define cmFftPlanAllocR fftwf_plan_dft_r2c_1d
#define cmFft1dPlanAllocR fftwf_plan_dft_1d
@ -71,6 +74,7 @@
#define cmCrealR creal
#define cmCimagR cimag
#define cmCargR carg
#define cmCconjR conj
#define cmFftPlanAllocR fftw_plan_dft_r2c_1d
#define cmFft1dPlanAllocR fftw_plan_dft_1d
@ -84,4 +88,12 @@
#endif
void cmVOCR_MultVVV( cmComplexR_t* y, const cmComplexS_t* x0, const cmComplexR_t* x1, unsigned n );
void cmVOCR_MultVFV( cmComplexR_t* y, const float* x, unsigned n );
void cmVOCR_DivVFV( cmComplexR_t* y, const float_t* x, unsigned n );
void cmVOCR_Abs( cmSample_t* y, const cmComplexR_t* x, unsigned n );
void cmVOCR_MultVS( cmComplexR_t* y, cmReal_t v, unsigned n );
void cmVOCR_DivVS( cmComplexR_t* y, cmReal_t v, unsigned n );
#endif

145
cmMath.c
View File

@ -1,5 +1,10 @@
#include "cmPrefix.h"
#include "cmGlobal.h"
#include "cmRpt.h"
#include "cmErr.h"
#include "cmCtx.h"
#include "cmMem.h"
#include "cmMallocDebug.h"
#include "cmFloatTypes.h"
#include "cmMath.h"
#include <sys/types.h> // u_char
@ -149,6 +154,20 @@ unsigned cmPrevOddU( unsigned v ) { return cmIsOddU(v) ? v : v-1; }
unsigned cmNextEvenU( unsigned v ) { return cmIsEvenU(v) ? v : v+1; }
unsigned cmPrevEvenU( unsigned v ) { return cmIsEvenU(v) ? v : v-1; }
unsigned cmModIncr(int idx, int delta, int maxN )
{
int sum = idx + delta;
if( sum >= maxN )
return sum - maxN;
if( sum < 0 )
return maxN + sum;
return sum;
}
// modified bessel function of first kind, order 0
// ref: orfandis appendix B io.m
double cmBessel0( double x )
@ -464,6 +483,128 @@ bool cmIsCloseU( unsigned x0, unsigned x1, double eps )
{
if( x0 == x1 )
return true;
return abs(x0-x1)/(x0+x1) < eps;
if( x0 > x1 )
return (x0-x1)/(x0+x1) < eps;
else
return (x1-x0)/(x0+x1) < eps;
}
//=================================================================
// cmLFSR() implementation based on note at bottom of:
// http://www.ece.cmu.edu/~koopman/lfsr/index.html
void cmLFSR( unsigned lfsrN, unsigned tapMask, unsigned seed, unsigned* yV, unsigned yN )
{
assert( 0 < lfsrN && lfsrN < 32 );
unsigned i;
for(i=0; i<yN; ++i)
{
if( (yV[i] = seed & 1)==1 )
seed = (seed >> 1) ^ tapMask;
else
seed = (seed >> 1);
}
}
bool cmMLS_IsBalanced( const unsigned* xV, int xN)
{
int a = 0;
unsigned i;
for(i=0; i<xN; ++i)
if( xV[i] == 1 )
++a;
return abs(a - (xN-a)) == 1;
}
unsigned _cmGenGoldCopy( int* y, unsigned yi, unsigned yN, unsigned* x, unsigned xN)
{
unsigned i;
for(i=0; i<xN; ++i,++yi)
y[yi] = x[i]==1 ? -1 : 1;
assert(yi <= yN);
return yi;
}
bool cmGenGoldCodes( unsigned lfsrN, unsigned poly_coeff0, unsigned poly_coeff1, unsigned goldN, int* yM, unsigned mlsN )
{
bool retFl = true;
unsigned yi = 0;
unsigned yN = goldN * mlsN;
unsigned* mls0V = cmMemAllocZ(unsigned,mlsN);
unsigned* mls1V = cmMemAllocZ(unsigned,mlsN);
unsigned* xorV = cmMemAllocZ(unsigned,mlsN);
unsigned i,j;
cmLFSR(lfsrN, poly_coeff0, 1 << (lfsrN-1), mls0V, mlsN);
cmLFSR(lfsrN, poly_coeff1, 1 << (lfsrN-1), mls1V, mlsN);
if( cmMLS_IsBalanced(mls0V,mlsN) )
yi = _cmGenGoldCopy(yM, yi, yN, mls0V, mlsN);
if( yi<yN && cmMLS_IsBalanced(mls1V,mlsN) )
yi = _cmGenGoldCopy(yM, yi, yN, mls1V, mlsN);
for(i=0; yi < yN && i<mlsN-1; ++i )
{
for(j=0; j<mlsN; ++j)
xorV[j] = (mls0V[j] + mls1V[ (i+j) % mlsN ]) % 2;
if( cmMLS_IsBalanced(xorV,mlsN) )
yi = _cmGenGoldCopy(yM,yi,yN,xorV,mlsN);
}
if(yi < yN )
{
//rc = cmErrMsg(err,kOpFailAtRC,"Gold code generation failed. Insuffient balanced pairs.");
retFl = false;
}
cmMemFree(mls0V);
cmMemFree(mls1V);
cmMemFree(xorV);
return retFl;
}
bool cmLFSR_Test()
{
// lfsrN = 5; % 5 6 7;
// poly_coeff0 = 0x12; % 0x12 0x21 0x41;
// poly_coeff1 = 0x1e; % 0x1e 0x36 0x72;
unsigned lfsrN = 7;
unsigned pc0 = 0x41;
unsigned pc1 = 0x72;
unsigned mlsN = (1 << lfsrN)-1;
unsigned yN = mlsN*2;
unsigned yV[ yN ];
unsigned i;
cmLFSR( lfsrN, pc0, 1 << (lfsrN-1), yV, yN );
for(i=0; i<mlsN; ++i)
if( yV[i] != yV[i+mlsN] )
return false;
//atVOU_PrintL(NULL,"0x12",yV,mlsN,2);
cmLFSR( lfsrN, pc1, 1 << (lfsrN-1), yV, yN );
//atVOU_PrintL(NULL,"0x17",yV,mlsN,2);
for(i=0; i<mlsN; ++i)
if( yV[i] != yV[i+mlsN] )
return false;
return true;
}

43
cmMath.h
View File

@ -15,6 +15,12 @@ unsigned cmPrevOddU( unsigned v );
unsigned cmNextEvenU( unsigned v );
unsigned cmPrevEvenU( unsigned v );
/// Increment or decrement 'idx' by 'delta' always wrapping the result into the range
/// 0 to (maxN-1).
/// 'idx': initial value
/// 'delta': incremental amount
/// 'maxN' - 1 : maximum return value.
unsigned cmModIncr(int idx, int delta, int maxN );
// modified bessel function of first kind, order 0
// ref: orfandis appendix B io.m
@ -74,4 +80,41 @@ bool cmIsCloseF( float x0, float x1, double eps );
bool cmIsCloseI( int x0, int x1, double eps );
bool cmIsCloseU( unsigned x0, unsigned x1, double eps );
//=================================================================
// Run a length 'lfsrN' linear feedback shift register (LFSR) for 'yN' iterations to
// produce a length 'yN' bit string in yV[yN].
// 'lfsrN' count of bits in the shift register range: 2<= lfsrN <= 32.
// 'tapMask' is a bit mask which gives the tap indexes positions for the LFSR.
// The least significant bit corresponds to the maximum delay tap position.
// The min tap position is therefore denoted by the tap mask bit location 1 << (lfsrN-1).
// A minimum of two taps must exist.
// 'seed' sets the initial delay state.
// 'yV[yN]' is the the output vector
// 'yN' is count of elements in yV.
// The function resturn kOkAtRC on success or kInvalidArgsRCRC if any arguments are invalid.
// /sa cmLFSR_Test.
void cmLFSR( unsigned lfsrN, unsigned tapMask, unsigned seed, unsigned* yV, unsigned yN );
// Example and test code for cmLFSR()
bool cmLFSR_Test();
// Generate a set of 'goldN' Gold codes using the Maximum Length Sequences (MLS) generated
// by a length 'lfsrN' linear feedback shift register.
// 'err' is an error object to be set if the the function fails.
// 'lfsrN' is the length of the Linear Feedback Shift Registers (LFSR) used to generate the MLS.
// 'poly_coeff0' tap mask for the first LFSR.
// 'coeff1' tap mask the the second LFSR.
// 'goldN' is the count of Gold codes to generate.
// 'yM[mlsN', goldN] is a column major output matrix where each column contains a Gold code.
// 'mlsN' is the length of the maximum length sequence for each Gold code which can be
// calculated as mlsN = (1 << a->lfsrN) - 1.
// Note that values of 'lfsrN' and the 'poly_coeffx' must be carefully selected such that
// they will produce a MLS. For example to generate a MLS with length 31 set 'lfsrN' to 5 and
// then select poly_coeff from two different elements of the set {0x12 0x14 0x17 0x1B 0x1D 0x1E}.
// See http://www.ece.cmu.edu/~koopman/lfsr/index.html for a complete set of MSL polynomial
// coefficients for given LFSR lengths.
// Returns false if insufficient balanced pairs exist.
bool cmGenGoldCodes( unsigned lfsrN, unsigned poly_coeff0, unsigned poly_coeff1, unsigned goldN, int* yM, unsigned mlsN );
#endif

2
cmProc.c
View File

@ -4393,7 +4393,7 @@ cmRC_t cmChmmTrain( cmChmm_t* p, const cmReal_t* oM, unsigned T, unsigned ite
cmReal_t logProb0 = cmChmmForward( p, oM, T, alphaM, logPrV );
// check for convergence
cmReal_t dLogProb = labs(logProb0-logProb) / ((labs(logProb0)+labs(logProb)+cmReal_EPSILON)/2);
cmReal_t dLogProb = fabs(logProb0-logProb) / ((fabs(logProb0)+fabs(logProb)+cmReal_EPSILON)/2);
if( dLogProb < thresh )
break;

57
cmProc2.c
View File

@ -823,7 +823,7 @@ cmRC_t cmFIRInitKaiser( cmFIR* p, unsigned procSmpCnt, double srate, double pass
// in practice the ripple must be equal in the stop and pass band - so take the minimum between the two
double d = cmMin(dPass,dStop);
// convert the ripple bcmk to db
// convert the ripple back to db
double A = -20 * log10(d);
// compute the kaiser alpha coeff
@ -914,7 +914,7 @@ cmRC_t cmFIRExec( cmFIR* p, const cmSample_t* sbp, unsigned sn )
// calc the output sample
while( cbp<cep)
{
// note that the delay is being iterated bcmkwards
// note that the delay is being iterated backwards
if( di == -1 )
di=delayCnt-1;
@ -936,7 +936,7 @@ cmRC_t cmFIRExec( cmFIR* p, const cmSample_t* sbp, unsigned sn )
return cmOkRC;
}
void cmFIRTest( cmRpt_t* rpt, cmLHeapH_t lhH, cmSymTblH_t stH )
void cmFIRTest0( cmRpt_t* rpt, cmLHeapH_t lhH, cmSymTblH_t stH )
{
unsigned N = 512;
cmKbRecd kb;
@ -978,6 +978,37 @@ void cmFIRTest( cmRpt_t* rpt, cmLHeapH_t lhH, cmSymTblH_t stH )
cmFIRFree(&ffp);
}
void cmFIRTest1( cmCtx* ctx )
{
const char* sfn = "/home/kevin/temp/sig.va";
const char* ffn = "/home/kevin/temp/fir.va";
unsigned N = 44100;
unsigned srate = N;
unsigned procSmpCnt = N;
double passHz = 15000;
double stopHz = 14000;
double passDb = 1.0;
double stopDb = 60.0;
unsigned flags = kHighPassFIRFl;
cmSample_t x[ procSmpCnt ];
cmVOS_Fill(x,procSmpCnt,0);
x[0] = 1;
cmVOS_Random(x,procSmpCnt, -1.0, 1.0 );
cmFIR* f = cmFIRAllocKaiser( ctx, NULL, procSmpCnt, srate, passHz, stopHz, stopDb, passDb, flags );
cmFIRExec( f, x, procSmpCnt );
cmVectArrayWriteMatrixS(ctx, ffn, f->outV, 1, f->outN );
cmVectArrayWriteMatrixS(ctx, sfn, x, 1, N );
cmFIRFree(&f);
}
//------------------------------------------------------------------------------------------------------------
@ -3966,6 +3997,14 @@ cmRC_t cmVectArrayWrite( cmVectArray_t* p, const char* fn )
return rc;
}
cmRC_t cmVectArrayWriteDirFn(cmVectArray_t* p, const char* dir, const char* fn )
{
assert( dir!=NULL && fn!=NULL );
const cmChar_t* path = cmFsMakeFn( dir, fn, NULL, NULL );
cmRC_t rc = cmVectArrayWrite(p,path);
cmFsFreeFn(path);
return rc;
}
cmRC_t cmVectArrayPrint( cmVectArray_t* p, cmRpt_t* rpt )
{
@ -4054,7 +4093,7 @@ cmRC_t _cmVectArrayWriteMatrix( cmCtx* ctx, const char* fn, unsigned flags, cons
memcpy(vv + ci*tbc, v + ci*rn*tbc, tbc );
// append the row to the VectArray
if((rc = cmVectArrayAppendV(p,v,cn)) != cmOkRC )
if((rc = cmVectArrayAppendV(p,vv,cn)) != cmOkRC )
{
rc = cmCtxRtCondition(&p->obj,rc,"Vector append failed in %s().",__FUNCTION__);
goto errLabel;
@ -5660,7 +5699,7 @@ cmRC_t cmExpanderInit( cmExpander* p,
p->envV[atkN+i] = p->rlsLvl + (G*i/rlsN);
}
printf("rmsN:%i atkN:%i rlsN:%i thr:%f %f rls:%f %f\n",p->rmsN,atkN,rlsN,threshDb,p->threshLvl,rlsDb,p->rlsLvl);
//printf("rmsN:%i atkN:%i rlsN:%i thr:%f %f rls:%f %f\n",p->rmsN,atkN,rlsN,threshDb,p->threshLvl,rlsDb,p->rlsLvl);
//for(i=0; i<p->envN; ++i)
// printf("%i %f\n",i,p->envV[i]);
@ -5680,7 +5719,7 @@ cmRC_t cmExpanderExec( cmExpander* p, cmSample_t* x, cmSample_t* y, unsigne
for(i=0; i<xyN; ++i)
{
// NOTE: using abs() instead of pow(x,2)
p->rmsV[p->rmsIdx] = abs(x[i]);
p->rmsV[p->rmsIdx] = fabsf(x[i]);
if( ++p->rmsIdx >= p->rmsN )
p->rmsIdx = 0;
@ -5910,7 +5949,7 @@ cmSpecDist_t* cmSpecDistAlloc( cmCtx* ctx,cmSpecDist_t* ap, unsigned procSmpCnt,
cmSpecDist_t* p = cmObjAlloc( cmSpecDist_t, ctx, ap );
//p->iSpecVa = cmVectArrayAlloc(ctx,kRealVaFl);
p->oSpecVa = cmVectArrayAlloc(ctx,kRealVaFl);
//p->oSpecVa = cmVectArrayAlloc(ctx,kRealVaFl);
if( procSmpCnt != 0 )
{
@ -5931,7 +5970,7 @@ cmRC_t cmSpecDistFree( cmSpecDist_t** pp )
cmSpecDistFinal(p);
//cmVectArrayFree(&p->iSpecVa);
cmVectArrayFree(&p->oSpecVa);
//cmVectArrayFree(&p->oSpecVa);
cmMemPtrFree(&p->hzV);
cmMemPtrFree(&p->iSpecM);
cmMemPtrFree(&p->oSpecM);
@ -6055,7 +6094,7 @@ cmRC_t cmSpecDistFinal(cmSpecDist_t* p )
cmRC_t rc = cmOkRC;
//cmVectArrayWrite(p->iSpecVa, "/home/kevin/temp/frqtrk/iSpec.va");
cmVectArrayWrite(p->oSpecVa, "/home/kevin/temp/expand/oSpec.va");
//cmVectArrayWrite(p->oSpecVa, "/home/kevin/temp/expand/oSpec.va");
cmPvAnlFree(&p->pva);
cmPvSynFree(&p->pvs);

24
cmProc2.h
View File

@ -190,7 +190,8 @@ extern "C" {
cmRC_t cmFIRInitSinc( cmFIR* p, unsigned procSmpCnt, double srate, unsigned sincSmpCnt, double fcHz, unsigned flags, const double* wndV );
cmRC_t cmFIRFinal( cmFIR* p );
cmRC_t cmFIRExec( cmFIR* p, const cmSample_t* sp, unsigned sn );
void cmFIRTest();
void cmFIRTest0( cmRpt_t* rpt, cmLHeapH_t lhH, cmSymTblH_t stH );
void cmFIRTest1( cmCtx* ctx );
//------------------------------------------------------------------------------------------------------------
// Apply a generic function to a windowed signal with a one sample hop size.
@ -771,7 +772,7 @@ extern "C" {
//------------------------------------------------------------------------------------------------------------
// cmVectArray buffers row vectors of arbitrary lenght in memory.
// cmVectArray buffers row vectors of arbitrary length in memory.
// The buffers may then be access using the cmVectArrayGetXXX() functions.
// The entire contents of the file may be written to a file using atVectArrayWrite().
// The file may then be read in back into memory using cmVectArrayAllocFromFile()
@ -837,8 +838,12 @@ extern "C" {
unsigned cmVectArrayMaxRowCount( const cmVectArray_t* p );
// Store a new vector by appending it to the end of the internal vector list.
// Note that the true type of v[] in the call to cmVectArrayAppendV() must match
// Note:
// 1. The true type of v[] in the call to cmVectArrayAppendV() must match
// the data type set in p->flags.
// 2. The 'vn' argument to atVectArrayAppendV() is an element count not
// a byte count. The size of each element is determined by the data type
// as set by atVectArrayAlloc().
cmRC_t cmVectArrayAppendV( cmVectArray_t* p, const void* v, unsigned vn );
cmRC_t cmVectArrayAppendS( cmVectArray_t* p, const cmSample_t* v, unsigned vn );
cmRC_t cmVectArrayAppendR( cmVectArray_t* p, const cmReal_t* v, unsigned vn );
@ -849,6 +854,7 @@ extern "C" {
// Write a vector array in a format that can be read by readVectArray.m.
cmRC_t cmVectArrayWrite( cmVectArray_t* p, const char* fn );
cmRC_t cmVectArrayWriteDirFn(cmVectArray_t* p, const char* dir, const char* fn );
// Print the vector array to rpt.
cmRC_t cmVectArrayPrint( cmVectArray_t* p, cmRpt_t* rpt );
@ -857,8 +863,12 @@ extern "C" {
unsigned cmVectArrayForEachS( cmVectArray_t* p, unsigned idx, unsigned cnt, cmVectArrayForEachFuncS_t func, void* arg );
// Write the vector v[vn] in the VectArray file format.
// Note that the true type of v[] in cmVectArrayWriteVectoV() must match the
// Note:
// 1. The true type of v[] in cmVectArrayWriteVectoV() must match the
// data type set in the 'flags' parameter.
// 2. The 'vn' argument to atVectArrayWriteVectorV() is an element count not
// a byte count. The size of each element is determined by the data type
// as set by atVectArrayAlloc().
cmRC_t cmVectArrayWriteVectorV( cmCtx* ctx, const char* fn, const void* v, unsigned vn, unsigned flags );
cmRC_t cmVectArrayWriteVectorS( cmCtx* ctx, const char* fn, const cmSample_t* v, unsigned vn );
cmRC_t cmVectArrayWriteVectorR( cmCtx* ctx, const char* fn, const cmReal_t* v, unsigned vn );
@ -868,8 +878,12 @@ extern "C" {
cmRC_t cmVectArrayWriteVectorU( cmCtx* ctx, const char* fn, const unsigned* v, unsigned vn );
// Write the column-major matrix m[rn,cn] to the file 'fn'.
// Note that the true type of m[] in cmVectArrayWriteMatrixV() must match the
// Notes:
// 1. The true type of m[] in cmVectArrayWriteMatrixV() must match the
// data type set in the 'flags' parameter.
// 2. The 'rn','cn' arguments to atVectWriteMatrixV() is are element counts not
// byte counts. The size of each element is determined by the data type
// as set by atVectArrayAlloc().
cmRC_t cmVectArrayWriteMatrixV( cmCtx* ctx, const char* fn, const void* m, unsigned rn, unsigned cn, unsigned flags );
cmRC_t cmVectArrayWriteMatrixS( cmCtx* ctx, const char* fn, const cmSample_t* m, unsigned rn, unsigned cn );
cmRC_t cmVectArrayWriteMatrixR( cmCtx* ctx, const char* fn, const cmReal_t* m, unsigned rn, unsigned cn );

743
cmProc5.c
View File

@ -16,7 +16,11 @@
#include "cmProcObj.h"
#include "cmProcTemplate.h"
#include "cmMath.h"
#include "cmFile.h"
#include "cmTime.h"
#include "cmMidi.h"
#include "cmProc.h"
#include "cmProc2.h"
#include "cmProc5.h"
#include "cmVectOps.h"
@ -123,3 +127,742 @@ cmRC_t cmGoertzelExec( cmGoertzel* p, const cmSample_t* inpV, unsigned procSmpCn
}
//=======================================================================================================================
double _cmGoldSigSinc( double t, double T )
{
double x = t/T;
return x == 0 ? 1.0 : sin(M_PI*x)/(M_PI*x);
}
void _cmGoldSigRaisedCos( cmSample_t* yV, int yN, double sPc, double beta )
{
int i;
for(i=0; i<yN; ++i)
{
double t = i - yN/2;
double den = 1 - (4*(beta*beta)*(t*t) / (sPc*sPc));
double a;
if(fabs(den) < 0.00001 )
a = 1;
else
a = cos(M_PI * beta * t/ sPc ) / den;
yV[i] = _cmGoldSigSinc(t,sPc) * a;
}
}
void _cmGoldSigConv( cmGoldSig_t* p, unsigned chIdx )
{
int i;
int sPc = p->a.samplesPerChip;
int osf = p->a.rcosOSFact;
// for each bit in the spreading-code
for(i=0; i<p->mlsN; ++i)
{
int j = (i*sPc) + sPc/2; // index into bbV[] of center of impulse response
int k = j - (sPc*osf)/2; // index into bbV[] of start of impulse response
int h;
// for each sample in the impulse response
for(h=0; h<p->rcosN; ++h,++k)
{
while( k<0 )
k += p->sigN;
while( k>=p->sigN )
k -= p->sigN;
p->ch[chIdx].bbV[k] += p->ch[chIdx].pnV[i] * p->rcosV[h];
}
}
}
void _cmGoldSigModulate( cmGoldSig_t* p, unsigned chIdx )
{
unsigned i;
double rps = 2.0 * M_PI * p->a.carrierHz / p->a.srate;
cmSample_t* yV = p->ch[chIdx].mdV;
cmSample_t* bbV = p->ch[chIdx].bbV;
for(i=0; i<p->sigN; ++i)
yV[ i ] = bbV[i]*cos(rps*i) + bbV[i]*sin(rps*i);
// apply a half Hann envelope to the onset/offset of the id signal
if( p->a.envMs > 0 )
{
unsigned wndMs = p->a.envMs * 2;
unsigned wndN = wndMs * p->a.srate / 1000;
wndN += wndN % 2 ? 0 : 1; // force the window length to be odd
unsigned wNo2 = wndN/2 + 1;
cmSample_t wndV[ wndN ];
cmVOS_Hann(wndV,wndN);
cmVOS_MultVV(yV,wNo2,wndV);
cmVOS_MultVV(yV + p->sigN - wNo2, wNo2, wndV + wNo2 - 1);
}
}
cmGoldSig_t* cmGoldSigAlloc( cmCtx* ctx, cmGoldSig_t* p, const cmGoldSigArg_t* a )
{
cmGoldSig_t* op = cmObjAlloc(cmGoldSig_t,ctx,p);
if( a != NULL )
if( cmGoldSigInit(op,a) != cmOkRC )
cmGoldSigFree(&op);
return op;
}
cmRC_t cmGoldSigFree( cmGoldSig_t** pp )
{
cmRC_t rc = cmOkRC;
if( pp == NULL || *pp == NULL )
return rc;
cmGoldSig_t* p = *pp;
if((rc = cmGoldSigFinal(p)) != cmOkRC )
return rc;
unsigned i;
for(i=0; i<p->a.chN; ++i)
{
cmMemFree(p->ch[i].bbV);
cmMemFree(p->ch[i].mdV);
}
cmMemFree(p->ch);
cmMemFree(p->rcosV);
cmMemFree(p->pnM);
cmMemFree(p);
*pp = NULL;
return rc;
}
cmRC_t cmGoldSigInit( cmGoldSig_t* p, const cmGoldSigArg_t* a )
{
cmRC_t rc = cmOkRC;
unsigned i;
p->a = *a; // store arg recd
p->ch = cmMemResizeZ(cmGoldSigCh_t,p->ch,a->chN); // alloc channel array
p->mlsN = (1 << a->lfsrN) - 1; // calc spreading code length
p->rcosN = a->samplesPerChip * a->rcosOSFact; // calc rcos imp. resp. length
p->rcosN += (p->rcosN % 2)==0; // force rcos imp. length odd
p->rcosV = cmMemResizeZ(cmSample_t,p->rcosV,p->rcosN); // alloc rcos imp. resp. vector
p->pnM = cmMemResizeZ(int,p->pnM,p->mlsN*a->chN); // alloc spreading-code mtx
p->sigN = p->mlsN * a->samplesPerChip; // calc audio signal length
// generate spreading codes
if( cmGenGoldCodes(a->lfsrN, a->mlsCoeff0, a->mlsCoeff1, a->chN, p->pnM, p->mlsN ) == false )
{
rc = cmCtxRtCondition(&p->obj,cmSubSysFailRC,"Unable to generate sufficient balanced Gold codes.");
goto errLabel;
}
// generate the rcos impulse response
_cmGoldSigRaisedCos(p->rcosV,p->rcosN,a->samplesPerChip,a->rcosBeta);
// for each channel
for(i=0; i<a->chN; ++i)
{
// Note: if (i*p->mlsN) is set to 0 in the following line then all channels
// will use the same spreading code.
p->ch[i].pnV = p->pnM + (i*p->mlsN); // get ch. spreading code
p->ch[i].bbV = cmMemResizeZ(cmSample_t,p->ch[i].bbV,p->sigN); // alloc baseband signal vector
p->ch[i].mdV = cmMemResizeZ(cmSample_t,p->ch[i].mdV,p->sigN); // alloc output audio vector
// Convolve spreading code with rcos impulse reponse to form baseband signal.
_cmGoldSigConv(p, i );
// Modulate baseband signal to carrier frq. and apply attack/decay envelope.
_cmGoldSigModulate(p, i );
}
errLabel:
if((rc = cmErrLastRC(&p->obj.err)) != cmOkRC )
cmGoldSigFree(&p);
return rc;
}
cmRC_t cmGoldSigFinal( cmGoldSig_t* p )
{ return cmOkRC; }
cmRC_t cmGoldSigWrite( cmCtx* ctx, cmGoldSig_t* p, const char* fn )
{
cmVectArray_t* vap = NULL;
unsigned i;
vap = cmVectArrayAlloc(ctx,kSampleVaFl);
for(i=0; i<p->a.chN; ++i)
{
cmVectArrayAppendS(vap,p->ch[i].bbV,p->sigN);
cmVectArrayAppendS(vap,p->ch[i].mdV,p->sigN);
}
cmVectArrayWrite(vap,fn);
cmVectArrayFree(&vap);
return cmOkRC;
}
cmRC_t cmGoldSigGen( cmGoldSig_t* p, unsigned chIdx, unsigned prefixN, unsigned dsN, unsigned *bsiV, unsigned bsiN, double noiseGain, cmSample_t** yVRef, unsigned* yNRef )
{
unsigned yN = prefixN + bsiN * (p->sigN + dsN);
cmSample_t* yV = cmMemAllocZ(cmSample_t,yN);
unsigned i;
cmVOS_Random(yV, yN, -noiseGain, noiseGain );
for(i=0; i<bsiN; ++i)
{
bsiV[i] = prefixN + i*(p->sigN + dsN);
cmVOS_AddVV(yV + bsiV[i], p->sigN, p->ch[chIdx].mdV );
}
if( yVRef != NULL )
*yVRef = yV;
if( yNRef != NULL )
*yNRef = yN;
return cmOkRC;
}
//=======================================================================================================================
cmPhat_t* cmPhatAlloc( cmCtx* ctx, cmPhat_t* ap, unsigned chN, unsigned hN, float alpha, unsigned mult, unsigned flags )
{
cmPhat_t* p = cmObjAlloc(cmPhat_t,ctx,ap);
// The FFT buffer and the delay line is at least twice the size of the
// id signal. This will guarantee that at least one complete id signal
// is inside the buffer. In practice it means that it is possible
// that there will be two id's in the buffer therefore if there are
// two correlation spikes it is important that we take the second.
unsigned fhN = cmNextPowerOfTwo(mult*hN);
// allocate the FFT object
cmFftAllocSR(ctx,&p->fft,NULL,fhN,kToPolarFftFl);
cmIFftAllocRS(ctx,&p->ifft,fhN/2 + 1 );
if( chN != 0 )
if( cmPhatInit(p,chN,hN,alpha,mult,flags) != cmOkRC )
cmPhatFree(&p);
return p;
}
cmRC_t cmPhatFree( cmPhat_t** pp )
{
cmRC_t rc = cmOkRC;
if( pp == NULL || *pp == NULL )
return rc;
cmPhat_t* p = *pp;
if((rc = cmPhatFinal(p)) != cmOkRC )
return rc;
cmMemFree(p->t0V);
cmMemFree(p->t1V);
cmMemFree(p->dV);
cmMemFree(p->xV);
cmMemFree(p->fhM);
cmMemFree(p->mhM);
cmMemFree(p->wndV);
cmObjFreeStatic(cmFftFreeSR, cmFftSR, p->fft);
cmObjFreeStatic(cmIFftFreeRS, cmIFftRS, p->ifft);
cmVectArrayFree(&p->ftVa);
cmObjFree(pp);
return rc;
}
cmRC_t cmPhatInit( cmPhat_t* p, unsigned chN, unsigned hN, float alpha, unsigned mult, unsigned flags )
{
cmRC_t rc = cmOkRC;
if((rc = cmPhatFinal(cmOkRC)) != cmOkRC )
return rc;
p->fhN = cmNextPowerOfTwo(mult*hN);
if((cmFftInitSR(&p->fft, NULL, p->fhN, kToPolarFftFl)) != cmOkRC )
return rc;
if((cmFftInitRS(&p->ifft, NULL, p->fft->binCnt )) != cmOkRC )
return rc;
p->alpha = alpha;
p->flags = flags;
// allocate the delay line
p->dV = cmMemResizeZ(cmSample_t,p->dV,p->fhN);
p->di = 0;
// allocate the linear buffer
p->xV = cmMemResizeZ(cmSample_t,p->xV,p->fhN);
p->t0V = cmMemResizeZ(cmComplexR_t,p->t0V,p->fhN);
p->t1V = cmMemResizeZ(cmComplexR_t,p->t1V,p->fhN);
// allocate the window function
p->wndV = cmMemResizeZ(cmSample_t,p->wndV,p->fhN);
cmVOS_Hann(p->wndV,p->fhN);
// allocate the signal id matrix
p->chN = chN;
p->hN = hN;
p->binN = p->fft.binCnt; //atFftRealBinCount(p->fftH);
p->fhM = cmMemResizeZ(cmComplexR_t, p->fhM, p->fhN * chN);
p->mhM = cmMemResizeZ(float, p->mhM, p->binN * chN);
cmPhatReset(p);
//if( cmIsFlag(p->flags,kDebugAtPhatFl))
// cmVectArrayAlloc(ctx, &p->ftVa, kSampleVaFl );
//else
// p->ftVa = NULL;
return rc;
}
cmRC_t cmPhatFinal( cmPhat_t* p )
{ return cmOkRC; }
cmRC_t cmPhatReset( cmPhat_t* p )
{
p->di = 0;
p->absIdx = 0;
cmVOS_Zero(p->dV,p->fhN);
return cmOkRC;
}
cmRC_t cmPhatSetId( cmPhat_t* p, unsigned chIdx, const cmSample_t* hV, unsigned hN )
{
unsigned i;
assert( chIdx < p->chN );
assert( hN == p->hN );
// Allocate a window vector
cmSample_t* wndV = cmMemAllocZ(cmSample_t,hN);
cmVOS_Hann(wndV,hN);
// get ptr to output column in p->fhM[].
cmComplexR_t* yV = p->fhM + (chIdx*p->fhN);
// Zero pad hV[hN] to p->fhN;
assert( hN <= p->fhN );
cmVOS_Zero(p->xV,p->fhN);
cmVOS_Copy(p->xV,hN,hV);
// Apply the window function to the id signal
if(cmIsFlag(p->flags,kHannAtPhatFl) )
cmVOS_MultVVV(p->xV,hN,hV,wndV);
// take FFT of id signal. The result is in fft->complexV and fft->magV,phsV
cmFftExecSR(&p->fft, p->xV, p->fhN );
// Store the magnitude of the id signal
//atFftComplexAbs(p->mhM + (chIdx*p->binN), yV, p->binN);
cmVOF_CopyR(p->mhM + (chIdx*p->binN), p->binN, p->fft.magV );
// Scale the magnitude
cmVOS_MultVS( p->mhM + (chIdx*p->binN), p->binN, p->alpha);
// store the complex conjugate of the FFT result in yV[]
//atFftComplexConj(yV,p->binN);
for(i=0; i<p->binN; ++i)
yV[i] = cmCconjR(p->fft.complexV[i]);
cmMemFree(wndV);
return cmOkRC;
}
cmSample_t* _cmPhatReadVector( cmCtx* ctx, cmPhat_t* p, const char* fn, unsigned* vnRef )
{
cmVectArray_t* vap = NULL;
cmSample_t* v = NULL;
cmRC_t rc = cmOkRC;
// instantiate a VectArray from a file
if( (vap = cmVectArrayAllocFromFile(ctx, fn )) == NULL )
{
rc = cmCtxRtCondition(&p->obj,cmSubSysFailRC,"Id component vector file read failed '%s'.",fn);
goto errLabel;
}
// get the count of elements in the vector
*vnRef = cmVectArrayEleCount(vap);
// allocate memory to hold the vector
v = cmMemAlloc(cmSample_t,*vnRef);
// copy the vector from the vector array object into v[]
if((rc = cmVectArrayGetF(vap,v,vnRef)) != cmOkRC )
{
cmMemFree(v);
v = NULL;
rc = cmCtxRtCondition(&p->obj,cmSubSysFailRC,"Id component vector copy out failed '%s'.",fn);
goto errLabel;
}
cmRptPrintf(p->obj.err.rpt,"%i : %s",*vnRef,fn);
errLabel:
cmVectArrayFree(&vap);
return v;
}
cmRC_t cmPhatExec( cmPhat_t* p, const cmSample_t* xV, unsigned xN )
{
unsigned n = cmMin(xN,p->fhN-p->di);
// update the delay line
cmVOS_Copy(p->dV+p->di,n,xV);
if( n < xN )
cmVOS_Copy(p->dV,xN-n,xV+n);
p->di = cmModIncr(p->di,xN,p->fhN);
// p->absIdx is the absolute sample index associated with di
p->absIdx += xN;
return cmOkRC;
}
void cmPhatChExec(
cmPhat_t* p,
unsigned chIdx,
unsigned sessionId,
unsigned roleId)
{
unsigned n0 = p->fhN - p->di;
unsigned n1 = p->fhN - n0;
// Linearize the delay line into xV[]
cmVOS_Copy(p->xV, n0, p->dV + p->di );
cmVOS_Copy(p->xV+n0, n1, p->dV );
if( cmIsFlag(p->flags,kDebugAtPhatFl))
cmVectArrayAppendS(p->ftVa, p->xV, p->fhN );
// apply a window function to the incoming signal
if( cmIsFlag(p->flags,kHannAtPhatFl) )
cmVOS_MultVV(p->xV,p->fhN,p->wndV);
// Take the FFT of the delay line.
// p->t0V[p->binN] = fft(p->xV)
//atFftRealForward(p->fftH, p->xV, p->fhN, p->t0V, p->binN );
cmFftExecSR(&p->fft, p->xV, p->fhN );
// Calc. the Cross Power Spectrum (aka cross spectral density) of the
// input signal with the id signal.
// Note that the CPS is equivalent to the Fourier Transform of the
// cross-correlation of the two signals.
// t0V[] *= p->fhM[:,chIdx]
//atFftComplexMult( p->t0V, p->fhM + (chIdx * p->fhN), p->binN );
cmVOCR_MultVVV( p->t0V, p->fft.complexV, p->fhM + (chIdx * p->fhN), p->binN);
// Calculate the magnitude of the CPS.
// xV[] = | t0V[] |
cmVOCR_Abs( p->xV, p->t0V, p->binN );
// Weight the CPS by the scaled magnitude of the id signal
// (we want to emphasize the limited frequencies where the
// id signal contains energy)
// t0V[] *= p->mhM[:,chIdx]
if( p->alpha > 0 )
cmVOCR_MultVFV( p->t0V, p->mhM + (chIdx*p->binN), p->binN);
// Divide through by the magnitude of the CPS
// This has the effect of whitening the spectram and thereby
// minimizing the effect of the magnitude correlation
// while maximimizing the effect of the phase correlation.
//
// t0V[] /= xV[]
cmVOCR_DivVFV( p->t0V, p->xV, p->binN );
// Take the IFFT of the weighted CPS to recover the cross correlation.
// xV[] = IFFT(t0V[])
cmIFftExecRS( p->ifft, );
//// ***** atFftRealInverse( p->fftH, p->t0V, p->xV, p->fhN );
// Shift the correlation spike to mark the end of the id
cmVOS_Rotate( p->xV, p->fhN, -((int)p->hN) );
// normalize by the length of the correlation
cmVOS_DivVS(p->xV,p->fhN,p->fhN);
if( cmIsFlag(p->flags,kDebugAtPhatFl))
{
cmVectArrayAppendS(p->ftVa, p->xV, p->fhN );
cmSample_t v[] = { sessionId, roleId };
cmVectArrayAppendS(p->ftVa, v, sizeof(v)/sizeof(v[0]));
}
}
cmRC_t cmPhatWrite( cmPhat_t* p, const char* dirStr )
{
cmRC_t rc = cmOkRC;
if( cmIsFlag(p->flags, kDebugAtPhatFl))
{
const char* path = NULL;
if( p->ftVa != NULL )
if((rc = cmVectArrayWrite(p->ftVa, path = cmFsMakeFn(path,"cmPhatFT","va",dirStr,NULL) )) != cmOkRC )
rc = cmCtxRtCondition(&p->obj,cmSubSysFailRC,"PHAT debug file write failed.");
cmFsFreeFn(path);
}
return rc;
}
#ifdef NOTDEF
cmRC_t cmPhatTest1( cmCtx* ctx, const char* dirStr )
{
cmRC_t rc = cmOkRC;
cmGoldSigArg_t sa;
cmGoldSig_t* s = NULL;
cmPhat_t* p = NULL;
char* path = NULL;
unsigned dspFrmCnt = 256;
unsigned listenDelaySmp = 8196;
double noiseGain = 0.05;
unsigned chIdx = 0;
cmSample_t* yV = NULL;
unsigned yN = 0;
double phatAlpha = 0.5;
unsigned phatMult = 4.0;
double nonLinExpo = 4.0;
cmVectArray_t* outVA = NULL;
cmVectArray_t* inVA = NULL;
cmVectArray_t* statusVA = NULL;
unsigned bsiN = 4;
unsigned bsiV[bsiN]; // known signal onset in absolute samples
unsigned esiV[bsiN]; // known signal offset
unsigned lsiV[bsiN]; // end of listen time (when cmPhatChExec()) is run.
unsigned dsiV[bsiN]; // detection time
unsigned i,j;
sa.chN = 1;
sa.srate = 44100.0;
sa.lfsrN = 8;
sa.mlsCoeff0 = 0x8e;
sa.mlsCoeff1 = 0x96;
sa.samplesPerChip = 64;
sa.rcosBeta = 0.5;
sa.rcosOSFact = 4;
sa.carrierHz = 17000.0;
sa.envMs = 50.0;
// allocate the the id signals
if( (s = cmGoldSigAlloc( ctx, NULL, &sa ) == NULL )
return cmErrMsg(&ctx->err, cmSubSysFailRC, "Signal allocate failed.");
// set the post signal listen delay to half the signal length
listenDelaySmp = s->sigN/2;
// allocate a PHAT detector
if( (p = cmPhatAlloc(ctx,NULL,sa.chN,s->sigN, phatAlpha, phatMult, kDebugAtPhatFl ) == NULL )
{
rc = cmErrMsg(&ctx->err, cmSubSysFailRC, "PHAT allocate failed.");
goto errLabel;
}
// register an id signal with the PHAT detector
if( cmPhatSetId(p, chIdx, s->ch[chIdx].mdV, s->sigN ) != cmOkRC )
{
rc = cmErrMsg(&ctx->err, cmSubSysFailRC, "PHAT setId failed.");
goto errLabel;
}
// generate an input test signal containing bsiN id signals
if( atSignalGen(s,chIdx,p->fhN,s->sigN,bsiV,bsiN,noiseGain,&yV,&yN) != cmOkRC )
{
rc = cmErrMsg(&ctx->err,cmSubSysFailRC,"Signal generation failed.");
goto errLabel;
}
// bsiV[] now holds signal onsets. Set esiV[] to signal offsets.
atVOU_AddVVS(esiV,bsiV,bsiN,s->sigN );
// set lsiV[] to end-of-listen location
atVOU_AddVVS(lsiV,esiV,bsiN,listenDelaySmp);
// zero the detection vector
atVOU_Zero(dsiV,bsiN);
// allocate a vector array to record the PHAT input signals
if( cmVectArrayAlloc(ctx,&inVA,kSampleVaFl) != cmOkRC )
{
rc = cmErrMsg(&ctx->err, cmSubSysFailRC, "vectArray inVA alloc failed.");
goto errLabel;
}
// allocate a vector array to record the PHAT correlation output signals
if( cmVectArrayAlloc(ctx,&outVA,kSampleVaFl) != cmOkRC )
{
rc = cmErrMsg(&ctx->err, cmSubSysFailRC, "vectArray outVA alloc failed.");
goto errLabel;
}
// allocate a vector array to record the PHAT status
if( cmVectArrayAlloc(ctx,&statusVA,kSampleVaFl) != cmOkRC )
{
rc = cmErrMsg(&ctx->err, cmSubSysFailRC, "vectArray statusVA alloc failed.");
goto errLabel;
}
// for each 'dspFrmCnt' samples in the input signal
for(i=0,j=0; j<bsiN && i<=yN-dspFrmCnt; i+=dspFrmCnt)
{
// store a copy of the input signal
cmVectArrayAppendS(inVA,yV+i,dspFrmCnt);
// feed the next dspFrmCnt samples to the PHAT detector
cmPhatExec(p,yV+i,dspFrmCnt);
// if the approximate end of an id signal is encountered
if( lsiV[j] <= i && i < lsiV[j] + dspFrmCnt )
{
// execute the PHAT correlator
cmPhatChExec( p, chIdx, -1, -1 );
// apply non-linear exponent to the correlation vector
cmVOS_PowV(p->xV,p->fhN,nonLinExpo);
// locate the corr. peak inside the listening window
// (the detection window is last 'detectWndSmp' samples in the corr. vector )
unsigned detectWndSmp = 2*listenDelaySmp;
dsiV[j] = cmVOS_ArgMax( p->xV + p->fhN - detectWndSmp, detectWndSmp);
// convert the pk index to absolute time
dsiV[j] = i + dspFrmCnt - detectWndSmp + dsiV[j];
// sig beg sig end detect begin dtct end detect
cmSample_t v[] = { bsiV[j], esiV[j], lsiV[j]-detectWndSmp, lsiV[j], dsiV[j] };
// store the detection information
cmVectArrayAppendS(statusVA,v,sizeof(v)/sizeof(v[0]));
// store the correlation output vector
cmVectArrayAppendS(outVA,p->xV,p->fhN);
j += 1;
}
}
// write inVA
if( cmVectArrayWrite(inVA,path = atMakePath(&ctx->err,path,"phatIn","va",dirStr,NULL)) != cmOkRC )
{
rc = cmErrMsg(&ctx->err, cmSubSysFailRC, "vectArray outVA write failed.");
goto errLabel;
}
// write outVA
if( cmVectArrayWrite(outVA,path = atMakePath(&ctx->err,path,"phatOut","va",dirStr,NULL)) != cmOkRC )
{
rc = cmErrMsg(&ctx->err, cmSubSysFailRC, "vectArray outVA write failed.");
goto errLabel;
}
// write statusVA
if( cmVectArrayWrite(statusVA,path = atMakePath(&ctx->err,path,"phatStatus","va",dirStr,NULL)) != cmOkRC )
{
rc = cmErrMsg(&ctx->err, cmSubSysFailRC, "vectArray statusVA write failed.");
goto errLabel;
}
errLabel:
cmVectArrayFree(&outVA);
cmVectArrayFree(&inVA);
if( cmPhatFree(&p) != cmOkRC )
cmErrMsg(&ctx->err,cmSubSysFailRC,"PHAT free failed.");
if( atSignalFree(&s) != cmOkRC )
cmErrMsg(&ctx->err,cmSubSysFailRC,"Signal free failed.");
return rc;
}
cmRC_t cmPhatTest2( cmCtx* ctx )
{
cmRC_t rc = cmOkRC;
cmPhat_t* p = NULL;
unsigned hN = 16;
float alpha = 1.0;
unsigned mult = 4;
cmSample_t hV[] = { 4,3,2,1, 0,0,0,0, 0,0,0,0, 0,0,0,0 };
cmSample_t x0V[] = { 4,3,2,1, 0,0,0,0, 0,0,0,0, 0,0,0,0 };
cmSample_t x1V[] = { 0,0,0,0, 4,3,2,1, 0,0,0,0, 0,0,0,0 };
cmSample_t x2V[] = { 0,0,0,0, 0,0,0,0, 4,3,2,1, 0,0,0,0 };
cmSample_t x3V[] = { 0,0,0,0, 0,0,0,0, 0,0,0,0, 4,3,2,1 };
cmSample_t* xV[] = { x0V, x1V, x2V, x3V };
unsigned chN = sizeof(xV)/sizeof(xV[0]);
unsigned i;
if(cmPhatAlloc(ctx,&p,chN,hN,alpha,mult,kNoFlagsAtPhatFl) != cmOkRC )
{
rc = cmErrMsg(&ctx->err,cmSubSysFailRC,"cmPhatAlloc() failed.");
goto errLabel;
}
for(i=0; i<chN; ++i)
if( cmPhatSetId(p,i,hV,hN) != cmOkRC )
rc = cmErrMsg(&ctx->err,cmSubSysFailRC,"cmPhatSetId() failed.");
for(i=0; i<chN; ++i)
{
cmPhatReset(p);
if( cmPhatExec(p,xV[i],hN) != cmOkRC )
{
rc = cmErrMsg(&ctx->err,cmSubSysFailRC,"cmPhatExec() failed.");
goto errLabel;
}
cmPhatChExec(p, i, -1, -1);
cmVOS_PrintL(&ctx->printRpt,"x:",p->xV,1,p->fhN);
}
errLabel:
cmPhatFree(&p);
return rc;
}
#endif

155
cmProc5.h
View File

@ -39,6 +39,161 @@ extern "C" {
cmRC_t cmGoertzelExec( cmGoertzel* p, const cmSample_t* in, unsigned procSmpCnt, double* outV, unsigned chCnt );
//=======================================================================================================================
// Gold Code Signal Generator
//
typedef struct
{
unsigned chN; // count of channels (each channel has a unique id)
double srate; // system sample rate (samples/second)
unsigned lfsrN; // linear feedback shift register (LFSR) length used to form Gold codes
unsigned mlsCoeff0; // LFSR coeff. 0
unsigned mlsCoeff1; // LFSR coeff. 1
unsigned samplesPerChip; // samples per spreading code bit
double rcosBeta; // raised cosine impulse response beta coeff.
unsigned rcosOSFact; // raised cosine impulse response oversample factor
double carrierHz; // carrier frequency
double envMs; // attack/decay envelope duration
} cmGoldSigArg_t;
typedef struct
{
int* pnV; // pnV[ mlsN ] spread code (aliased from pnM[:,i])
cmSample_t* bbV; // bbV[ sigN ] baseband signal at audio rate
cmSample_t* mdV; // mdV[ sigN ] modulated signal at audio rate
} cmGoldSigCh_t;
typedef struct
{
cmObj obj; //
cmGoldSigArg_t a; // argument record
cmGoldSigCh_t* ch; // ch[ chN ] channel array
int* pnM; // pnM[mlsN,chN] (aliased to ch[].pnV)
cmSample_t* rcosV; // rcosV[rcosN] raised cosine impulse response
unsigned rcosN; // length of raised cosine impulse response
unsigned mlsN; // length of Gold codes (Maximum length sequence length)
unsigned sigN; // length of channel signals bbV[] and mdV[]
} cmGoldSig_t;
cmGoldSig_t* cmGoldSigAlloc( cmCtx* ctx, cmGoldSig_t* p, const cmGoldSigArg_t* a );
cmRC_t cmGoldSigFree( cmGoldSig_t** pp );
cmRC_t cmGoldSigInit( cmGoldSig_t* p, const cmGoldSigArg_t* a );
cmRC_t cmGoldSigFinal( cmGoldSig_t* p );
cmRC_t cmGoldSigWrite( cmCtx* ctx, cmGoldSig_t* p, const char* fn );
// Generate a signal consisting of underlying white noise with
// bsiN repeated copies of the id signal associated with
// channel 'chIdx'. Each discrete id signal copy is separated by 'dsN' samples.
// The signal will be prefixed with 'prefixN' samples of silence (noise).
// On return sets 'yVRef' to point to the generated signal and 'yNRef'
// to the count of samples in 'yVRef'.
// On error sets yVRef to NULL and yNRef to zero.
// The vector returned in 'yVRef' should be freed via atMemFree().
// On return sets bsiV[bsiN] to the onset sample index of each id signal copy.
// The background noise signal is limited to the range -noiseGain to noiseGain.
cmRC_t cmGoldSigGen(
cmGoldSig_t* p,
unsigned chIdx,
unsigned prefixN,
unsigned dsN,
unsigned *bsiV,
unsigned bsiN,
double noiseGain,
cmSample_t** yVRef,
unsigned* yNRef );
cmRC_t cmGoldSigTest( cmCtx* ctx );
//=======================================================================================================================
// Phase aligned transform generalized cross correlator
//
// Flags for use with the 'flags' argument to cmPhatAlloc()
enum
{
kNoFlagsAtPhatFl= 0x00,
kDebugAtPhatFl = 0x01, // generate debugging file
kHannAtPhatFl = 0x02 // apply a hann window function to the id/audio signals prior to correlation.
};
typedef struct
{
cmObj obj;
cmFftSR fft;
cmIFftRS ifft;
float alpha;
unsigned flags;
cmComplexR_t* fhM; // fhM[fhN,chN] FT of each id signal stored in complex form
float* mhM; // mhM[binN,chN] magnitude of each fhM column
unsigned chN; // count of id signals
unsigned fhN; // length of each FT id signal (fft->xN)
unsigned binN; // length of each mhM column (fft->xN/2);
unsigned hN; // length of each time domain id signal (hN<=fhN/2)
unsigned absIdx; // abs. sample index of p->di
cmSample_t* dV; // dV[fhN] delay line
unsigned di; // next input into delay line
cmSample_t* xV; // xV[fhN] linear delay buffer
cmComplexR_t* t0V; // t0V[fhN]
cmComplexR_t* t1V; // t1V[fhN]
cmSample_t* wndV;
cmVectArray_t* ftVa;
} cmPhat_t;
// Allocate a PHAT based multi-channel correlator.
// 'chN' is the maximum count of id signals to be set via cmPhatSetId().
// 'hN' is the the length of the id signal in samples.
// 'alpha' weight used to emphasize the frequencies where the
// id signal contains energy.
// 'mult' * 'hN' is the correlation length (fhN)
// 'flags' See kDebugAtPhatFl and kWndAtPhatFl.
cmPhat_t* cmPhatAlloc( cmCtx* ctx, cmPhat_t* p, unsigned chN, unsigned hN, float alpha, unsigned mult, unsigned flags );
cmRC_t cmPhatFree( cmPhat_t** pp );
cmRC_t cmPhatInit( cmPhat_t* p, unsigned chN, unsigned hN, float alpha, unsigned mult, unsigned flags );
cmRC_t cmPhatFinal( cmPhat_t* p );
// Zero the audio delay line and reset the current input sample (di)
// and absolute time index (absIdx) to 0.
cmRC_t cmPhatReset( cmPhat_t* p );
// Register an id signal with the correlator.
cmRC_t cmPhatSetId( cmPhat_t* p, unsigned chIdx, const cmSample_t* hV, unsigned hN );
// Update the correlators internal delay buffer.
cmRC_t cmPhatExec( cmPhat_t* p, const cmSample_t* xV, unsigned xN );
// Set p->xV[0:fhN-1] to the correlation function based on
// correlation between the current audio delay line d[] and
// the id signal in fhM[:,chIdx].
// 'sessionId' and 'roleId' are only used to label the
// data stored in the debug file and may be set to any
// arbitrary value if the debug files are not being generated.
void cmPhatChExec(
cmPhat_t* p,
unsigned chIdx,
unsigned sessionId,
unsigned roleId);
cmRC_t cmPhatWrite( cmPhat_t* p, const char* dirStr );
cmRC_t cmPhatTest1( cmCtx* ctx, const char* dirFn );
cmRC_t cmPhatTest2( cmCtx* ctx );
#ifdef __cplusplus
}

2
dsp/cmDspKr.c
View File

@ -2534,7 +2534,7 @@ cmDspRC_t _cmDspRecdPlayParseRsrc( cmDspCtx_t* ctx, cmDspInst_t* inst, cmRecdPla
if( jnp == NULL || cmJsonIsArray(jnp)==false )
{
// this is really a warning - the object does not require preloaded segments.
cmDspInstErr(ctx,inst,kRsrcNotFoundDspRC,"The 'recdPlay' resource used to define pre-loaded segments was not found.");
cmDspInstErr(ctx,inst,kRsrcNotFoundDspRC,"Warning: The 'recdPlay' resource used to define pre-loaded segments was not found.");
return kOkDspRC;
}

4
dsp/cmDspPgm.c
View File

@ -334,6 +334,7 @@ cmDspRC_t _cmDspSysPgm_PlaySine( cmDspSysH_t h, void** userPtrPtr )
cmDspInst_t* ao1p = cmDspSysAllocInst(h,"AudioOut", NULL, 1, useBuiltInFl ? 1 : 3 );
cmDspInst_t* om0p = cmDspSysAllocInst(h,"AMeter","Out", 0);
cmDspInst_t* gain= cmDspSysAllocInst( h,"Scalar", "Gain", 5, kNumberDuiId, 0.0, 10.0, 0.01, 0.0);
cmDspSysConnectAudio(h, php, "out", wtp, "phs" ); // phasor -> wave table
cmDspSysConnectAudio(h, wtp, "out", ao0p, "in" ); // wave table -> audio out
@ -341,6 +342,9 @@ cmDspRC_t _cmDspSysPgm_PlaySine( cmDspSysH_t h, void** userPtrPtr )
cmDspSysConnectAudio(h, wtp, "out", om0p, "in" );
cmDspSysInstallCb( h, chp, "val", ao0p, "ch", NULL);
cmDspSysInstallCb( h, gain, "val", ao0p, "gain", NULL);
cmDspSysInstallCb( h, gain, "val", ao1p, "gain", NULL);
return kOkDspRC;
}

4
dsp/cmDspPgmKr.c
View File

@ -382,7 +382,7 @@ cmDspRC_t _cmDspSysPgm_TimeLine(cmDspSysH_t h, void** userPtrPtr )
cmErr_t err;
krRsrc_t r;
bool fragFl = false;
bool useWtFl = false;
bool useWtFl = true;
bool useChain1Fl = true;
bool useInputEqFl = false;
bool useInCompFl = true;
@ -1048,7 +1048,7 @@ cmDspRC_t _cmDspSysPgm_TimeLine(cmDspSysH_t h, void** userPtrPtr )
cmDspSysInstallCb(h, siRt, "f-out-1", sfp, "smpidx",NULL );
// leave siRt.f-out-1 unconnected because it should be ignored in 'simulate mode'
cmDspSysInstallCb(h, mfp, "mu", muRt, "f-in", NULL );
cmDspSysInstallCb(h, mfp, "id", muRt, "f-in", NULL );
cmDspSysInstallCb(h, muRt, "f-out-1", sfp, "muid", NULL );
// leave muRt.f-out-1 unconnected because it should be ignored in 'simulate mode'

158
linux/cmAudioPortAlsa.c
View File

@ -46,6 +46,12 @@ typedef struct devRecd_str
unsigned iBits; // bits per sample
unsigned oBits;
bool iSignFl; // sample type is signed
bool oSignFl;
bool iSwapFl; // swap the sample bytes
bool oSwapFl;
unsigned iSigBits; // significant bits in each sample beginning
unsigned oSigBits; // with the most sig. bit.
@ -362,6 +368,77 @@ void _cmApDevRtReport( cmRpt_t* rpt, cmApDevRecd_t* drp )
}
void _cmApDevReportFormats( cmRpt_t* rpt, snd_pcm_hw_params_t* hwParams )
{
snd_pcm_format_mask_t* mask;
snd_pcm_format_t fmt[] =
{
SND_PCM_FORMAT_S8,
SND_PCM_FORMAT_U8,
SND_PCM_FORMAT_S16_LE,
SND_PCM_FORMAT_S16_BE,
SND_PCM_FORMAT_U16_LE,
SND_PCM_FORMAT_U16_BE,
SND_PCM_FORMAT_S24_LE,
SND_PCM_FORMAT_S24_BE,
SND_PCM_FORMAT_U24_LE,
SND_PCM_FORMAT_U24_BE,
SND_PCM_FORMAT_S32_LE,
SND_PCM_FORMAT_S32_BE,
SND_PCM_FORMAT_U32_LE,
SND_PCM_FORMAT_U32_BE,
SND_PCM_FORMAT_FLOAT_LE,
SND_PCM_FORMAT_FLOAT_BE,
SND_PCM_FORMAT_FLOAT64_LE,
SND_PCM_FORMAT_FLOAT64_BE,
SND_PCM_FORMAT_IEC958_SUBFRAME_LE,
SND_PCM_FORMAT_IEC958_SUBFRAME_BE,
SND_PCM_FORMAT_MU_LAW,
SND_PCM_FORMAT_A_LAW,
SND_PCM_FORMAT_IMA_ADPCM,
SND_PCM_FORMAT_MPEG,
SND_PCM_FORMAT_GSM,
SND_PCM_FORMAT_SPECIAL,
SND_PCM_FORMAT_S24_3LE,
SND_PCM_FORMAT_S24_3BE,
SND_PCM_FORMAT_U24_3LE,
SND_PCM_FORMAT_U24_3BE,
SND_PCM_FORMAT_S20_3LE,
SND_PCM_FORMAT_S20_3BE,
SND_PCM_FORMAT_U20_3LE,
SND_PCM_FORMAT_U20_3BE,
SND_PCM_FORMAT_S18_3LE,
SND_PCM_FORMAT_S18_3BE,
SND_PCM_FORMAT_U18_3LE,
SND_PCM_FORMAT_U18_3BE,
SND_PCM_FORMAT_G723_24,
SND_PCM_FORMAT_G723_24_1B,
SND_PCM_FORMAT_G723_40,
SND_PCM_FORMAT_G723_40_1B,
SND_PCM_FORMAT_DSD_U8,
//SND_PCM_FORMAT_DSD_U16_LE,
//SND_PCM_FORMAT_DSD_U32_LE,
//SND_PCM_FORMAT_DSD_U16_BE,
//SND_PCM_FORMAT_DSD_U32_BE,
SND_PCM_FORMAT_UNKNOWN
};
snd_pcm_format_mask_alloca(&mask);
snd_pcm_hw_params_get_format_mask(hwParams,mask);
cmRptPrintf(rpt,"Formats: " );
int i;
for(i=0; fmt[i]!=SND_PCM_FORMAT_UNKNOWN; ++i)
if( snd_pcm_format_mask_test(mask, fmt[i] ))
cmRptPrintf(rpt,"%s%s",snd_pcm_format_name(fmt[i]), snd_pcm_format_cpu_endian(fmt[i]) ? " " : " (swap) ");
cmRptPrintf(rpt,"\n");
}
void _cmApDevReport( cmRpt_t* rpt, cmApDevRecd_t* drp )
{
bool inputFl = true;
@ -435,6 +512,8 @@ void _cmApDevReport( cmRpt_t* rpt, cmApDevRecd_t* drp )
ioLabel,minChCnt,maxChCnt,minSrate,maxSrate,minPeriodFrmCnt,maxPeriodFrmCnt,minBufFrmCnt,maxBufFrmCnt,
(snd_pcm_hw_params_is_half_duplex(hwParams) ? "yes" : "no"),
(snd_pcm_hw_params_is_joint_duplex(hwParams) ? "yes" : "no"));
_cmApDevReportFormats( rpt, hwParams );
}
if((err = snd_pcm_close(pcmH)) < 0)
@ -610,6 +689,28 @@ void _cmApStateRecover( snd_pcm_t* pcmH, cmApDevRecd_t* drp, bool inputFl )
}
void _cmApS24_3BE_to_Float( const char* x, cmApSample_t* y, unsigned n )
{
unsigned i;
for(i=0; i<n; ++i,x+=3)
{
int s = (((int)x[0])<<16) + (((int)x[1])<<8) + (((int)x[2]));
y[i] = ((cmApSample_t)s)/0x7fffff;
}
}
void _cmApS24_3BE_from_Float( const cmApSample_t* x, char* y, unsigned n )
{
unsigned i;
for(i=0; i<n; ++i)
{
int s = x[i] * 0x7fffff;
y[i*3+2] = (char)((s & 0x7f0000) >> 16);
y[i*3+1] = (char)((s & 0x00ff00) >> 8);
y[i*3+0] = (char)((s & 0x0000ff) >> 0);
}
}
// Returns count of frames written on success or < 0 on error;
// set smpPtr to NULL to write a buffer of silence
@ -648,9 +749,13 @@ int _cmApWriteBuf( cmApDevRecd_t* drp, snd_pcm_t* pcmH, const cmApSample_t* sp,
case 24:
{
// for use w/ MBox
//_cmApS24_3BE_from_Float(sp, obuf, ep-sp );
int* dp = (int*)obuf;
while( sp < ep )
*dp++ = (int)(*sp++ * 0x7fffff);
}
break;
@ -696,7 +801,6 @@ int _cmApWriteBuf( cmApDevRecd_t* drp, snd_pcm_t* pcmH, const cmApSample_t* sp,
}
// Returns frames read on success or < 0 on error.
// Set smpPtr to NULL to read the incoming buffer and discard it
int _cmApReadBuf( cmApDevRecd_t* drp, snd_pcm_t* pcmH, cmApSample_t* smpPtr, unsigned chCnt, unsigned frmCnt, unsigned bits, unsigned sigBits )
@ -729,7 +833,6 @@ int _cmApReadBuf( cmApDevRecd_t* drp, snd_pcm_t* pcmH, cmApSample_t* smpPtr, uns
// setup the return buffer
cmApSample_t* dp = smpPtr;
cmApSample_t* ep = dp + cmMin(smpCnt,err*chCnt);
switch(bits)
@ -752,6 +855,8 @@ int _cmApReadBuf( cmApDevRecd_t* drp, snd_pcm_t* pcmH, cmApSample_t* smpPtr, uns
case 24:
{
// For use with MBox
//_cmApS24_3BE_to_Float(buf, dp, ep-dp );
int* sp = (int*)buf;
while(dp < ep)
*dp++ = ((cmApSample_t)*sp++) / 0x7fffff;
@ -819,7 +924,7 @@ void _cmApStaticAsyncHandler( snd_async_handler_t* ahandler )
while( (avail = snd_pcm_avail_update(pcmH)) >= (snd_pcm_sframes_t)frmCnt )
{
// Handle inpuut
// Handle input
if( inputFl )
{
// read samples from the device
@ -1024,8 +1129,22 @@ bool _cmApDevSetup( cmApDevRecd_t *drp, unsigned srate, unsigned framesPerCycle,
snd_pcm_uframes_t bufferFrameCnt;
unsigned bits = 0;
int sig_bits = 0;
bool signFl = true;
bool swapFl = false;
cmApRoot_t* p = drp->rootPtr;
snd_pcm_format_t fmt[] =
{
SND_PCM_FORMAT_S32_LE,
SND_PCM_FORMAT_S32_BE,
SND_PCM_FORMAT_S24_LE,
SND_PCM_FORMAT_S24_BE,
SND_PCM_FORMAT_S24_3LE,
SND_PCM_FORMAT_S24_3BE,
SND_PCM_FORMAT_S16_LE,
SND_PCM_FORMAT_S16_BE,
};
// setup input, then output device
for(i=0; i<2; i++,inputFl=!inputFl)
@ -1041,7 +1160,6 @@ bool _cmApDevSetup( cmApDevRecd_t *drp, unsigned srate, unsigned framesPerCycle,
if( _cmApDevShutdown(p, drp, inputFl ) != kOkApRC )
retFl = false;
// attempt to open the sub-device
if((err = snd_pcm_open(&pcmH,drp->nameStr, inputFl ? SND_PCM_STREAM_CAPTURE : SND_PCM_STREAM_PLAYBACK, 0)) < 0 )
retFl = _cmApDevSetupError(p,err,inputFl,drp,"Unable to open the PCM handle");
@ -1076,20 +1194,20 @@ bool _cmApDevSetup( cmApDevRecd_t *drp, unsigned srate, unsigned framesPerCycle,
if((err = snd_pcm_hw_params_set_access(pcmH,hwParams,SND_PCM_ACCESS_RW_INTERLEAVED )) < 0 )
retFl = _cmApDevSetupError(p,err,inputFl, drp, "Unable to set access to: RW Interleaved");
// select the widest possible sample width
if((err = snd_pcm_hw_params_set_format(pcmH,hwParams,SND_PCM_FORMAT_S32)) >= 0 )
bits = 32;
else
{
if((err = snd_pcm_hw_params_set_format(pcmH,hwParams,SND_PCM_FORMAT_S24)) >= 0 )
bits = 24;
else
{
if((err = snd_pcm_hw_params_set_format(pcmH,hwParams,SND_PCM_FORMAT_S16)) >= 0 )
bits = 16;
else
// select the format width
int j;
int fmtN = sizeof(fmt)/sizeof(fmt[0]);
for(j=0; j<fmtN; ++j)
if((err = snd_pcm_hw_params_set_format(pcmH,hwParams,fmt[j])) >= 0 )
break;
if( j == fmtN )
retFl = _cmApDevSetupError(p,err,inputFl, drp, "Unable to set format to: S16");
}
else
{
bits = snd_pcm_format_width(fmt[j]); // bits per sample
signFl = snd_pcm_format_signed(fmt[j]);
swapFl = !snd_pcm_format_cpu_endian(fmt[j]);
}
sig_bits = snd_pcm_hw_params_get_sbits(hwParams);
@ -1167,6 +1285,8 @@ bool _cmApDevSetup( cmApDevRecd_t *drp, unsigned srate, unsigned framesPerCycle,
{
drp->iBits = bits;
drp->iSigBits = sig_bits;
drp->iSignFl = signFl;
drp->iSwapFl = swapFl;
drp->iPcmH = pcmH;
drp->iBuf = cmMemResizeZ( cmApSample_t, drp->iBuf, actFpC * drp->iChCnt );
drp->iFpC = actFpC;
@ -1175,6 +1295,8 @@ bool _cmApDevSetup( cmApDevRecd_t *drp, unsigned srate, unsigned framesPerCycle,
{
drp->oBits = bits;
drp->oSigBits = sig_bits;
drp->oSignFl = signFl;
drp->oSwapFl = swapFl;
drp->oPcmH = pcmH;
drp->oBuf = cmMemResizeZ( cmApSample_t, drp->oBuf, actFpC * drp->oChCnt );
drp->oFpC = actFpC;
@ -1448,7 +1570,7 @@ cmApRC_t cmApAlsaInitialize( cmRpt_t* rpt, unsigned baseApDevIdx )
// this device uses this subdevice in the current direction
dr.flags += inputFl ? kInFl : kOutFl;
printf("%s in:%i chs:%i rate:%i\n",dr.nameStr,inputFl,*chCntPtr,rate);
//printf("%s in:%i chs:%i rate:%i\n",dr.nameStr,inputFl,*chCntPtr,rate);
}