cmProc2.h/c: Many changes to cmVectArray. Fixed high-pass generation in cmFIR.

This commit is contained in:
Kevin Larke 2015-05-22 14:12:04 -07:00
parent 834b6f421f
commit 962ec25635
2 changed files with 534 additions and 387 deletions

673
cmProc2.c
View File

@ -766,22 +766,22 @@ void cmDelayTest( cmRpt_t* rpt, cmLHeapH_t lhH, cmSymTblH_t stH )
//------------------------------------------------------------------------------------------------------------ //------------------------------------------------------------------------------------------------------------
cmFIR* cmFIRAllocKaiser(cmCtx* c, cmFIR* p, unsigned procSmpCnt, double srate, double passHz, double stopHz, double passDb, double stopDb ) cmFIR* cmFIRAllocKaiser(cmCtx* c, cmFIR* p, unsigned procSmpCnt, double srate, double passHz, double stopHz, double passDb, double stopDb, unsigned flags )
{ {
cmFIR* op = cmObjAlloc( cmFIR, c, p ); cmFIR* op = cmObjAlloc( cmFIR, c, p );
if( procSmpCnt > 0 && srate > 0 ) if( procSmpCnt > 0 && srate > 0 )
if( cmFIRInitKaiser(op,procSmpCnt,srate,passHz,stopHz,passDb,stopDb) != cmOkRC ) if( cmFIRInitKaiser(op,procSmpCnt,srate,passHz,stopHz,passDb,stopDb,flags) != cmOkRC )
cmObjFree(&op); cmObjFree(&op);
return op; return op;
} }
cmFIR* cmFIRAllocSinc( cmCtx* c, cmFIR* p, unsigned procSmpCnt, double srate, unsigned sincSmpCnt, double fcHz, unsigned flags ) cmFIR* cmFIRAllocSinc( cmCtx* c, cmFIR* p, unsigned procSmpCnt, double srate, unsigned sincSmpCnt, double fcHz, unsigned flags, const double* wndV )
{ {
cmFIR* op = cmObjAlloc( cmFIR, c, p ); cmFIR* op = cmObjAlloc( cmFIR, c, p );
if( srate > 0 && sincSmpCnt > 0 ) if( srate > 0 && sincSmpCnt > 0 )
if( cmFIRInitSinc(op,procSmpCnt,srate,sincSmpCnt,fcHz,flags) != cmOkRC ) if( cmFIRInitSinc(op,procSmpCnt,srate,sincSmpCnt,fcHz,flags,wndV) != cmOkRC )
cmObjFree(&op); cmObjFree(&op);
return op; return op;
} }
@ -794,7 +794,7 @@ cmRC_t cmFIRFree( cmFIR** pp )
{ {
cmFIR* p = *pp; cmFIR* p = *pp;
if((rc = cmFIRFinal(*pp)) != cmOkRC ) if((rc = cmFIRFinal(*pp)) == cmOkRC )
{ {
cmMemPtrFree(&p->coeffV); cmMemPtrFree(&p->coeffV);
cmMemPtrFree(&p->outV); cmMemPtrFree(&p->outV);
@ -806,7 +806,7 @@ cmRC_t cmFIRFree( cmFIR** pp )
return rc; return rc;
} }
cmRC_t cmFIRInitKaiser( cmFIR* p, unsigned procSmpCnt, double srate, double passHz, double stopHz, double passDb, double stopDb ) cmRC_t cmFIRInitKaiser( cmFIR* p, unsigned procSmpCnt, double srate, double passHz, double stopHz, double passDb, double stopDb, unsigned flags )
{ {
// pass/stop frequencies above nyquist produce incorrect results // pass/stop frequencies above nyquist produce incorrect results
assert( passHz <= srate/2 && stopHz<=srate/2); assert( passHz <= srate/2 && stopHz<=srate/2);
@ -815,13 +815,12 @@ cmRC_t cmFIRInitKaiser( cmFIR* p, unsigned procSmpCnt, double srate, double pass
double fcHz = (passHz + stopHz) / 2; // fc is half way between passHz and stopHz double fcHz = (passHz + stopHz) / 2; // fc is half way between passHz and stopHz
double dHz = fabs(stopHz-passHz); double dHz = fabs(stopHz-passHz);
//double signFcmt = stopHz > passHz ? 1 : -1;
// convert ripple spec from db to linear // convert ripple spec from db to linear
double dPass = (pow(10,passDb/20)-1) / (pow(10,passDb/20)+1); double dPass = (pow(10,passDb/20)-1) / (pow(10,passDb/20)+1);
double dStop = pow(10,-stopDb/20); double dStop = pow(10,-stopDb/20);
// in prcmtice the ripple must be equal in the stop and pass band - so take the minimum between the two // 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); double d = cmMin(dPass,dStop);
// convert the ripple bcmk to db // convert the ripple bcmk to db
@ -835,7 +834,6 @@ cmRC_t cmFIRInitKaiser( cmFIR* p, unsigned procSmpCnt, double srate, double pass
{ {
if( A > 21 ) if( A > 21 )
alpha = (0.5842 * (pow(A-21.0,0.4))) + (0.07886*(A-21)); alpha = (0.5842 * (pow(A-21.0,0.4))) + (0.07886*(A-21));
} }
double D = 0.922; double D = 0.922;
@ -843,7 +841,6 @@ cmRC_t cmFIRInitKaiser( cmFIR* p, unsigned procSmpCnt, double srate, double pass
if( A > 21 ) if( A > 21 )
D = (A - 7.95) / 14.36; D = (A - 7.95) / 14.36;
// compute the filter order // compute the filter order
unsigned N = (unsigned)(floor(D * srate / dHz) + 1); unsigned N = (unsigned)(floor(D * srate / dHz) + 1);
@ -853,18 +850,12 @@ cmRC_t cmFIRInitKaiser( cmFIR* p, unsigned procSmpCnt, double srate, double pass
//printf("fc=%f df=%f dPass=%f dStop=%f d=%f alpha=%f A=%f D=%f N=%i\n",fcHz,dHz,dPass,dStop,d,alpha,A,D,N); //printf("fc=%f df=%f dPass=%f dStop=%f d=%f alpha=%f A=%f D=%f N=%i\n",fcHz,dHz,dPass,dStop,d,alpha,A,D,N);
// form an ideal FIR LP impulse response based on a sinc function
cmFIRInitSinc(p,procSmpCnt,srate,N,fcHz,0);
// compute a kaiser function to truncate the sinc // compute a kaiser function to truncate the sinc
double wnd[ N ]; double wnd[ N ];
cmVOD_Kaiser( wnd, N, alpha ); cmVOD_Kaiser( wnd, N, alpha );
// apply the kaiser window to the sinc function // form an ideal FIR LP impulse response based on a sinc function
cmVOD_MultVV( p->coeffV, p->coeffCnt, wnd ); cmFIRInitSinc(p,procSmpCnt,srate,N,fcHz,flags, wnd);
double sum = cmVOD_Sum(p->coeffV,p->coeffCnt);
cmVOD_DivVS(p->coeffV,p->coeffCnt,sum );
//cmVOD_Print(stdout,1,p->coeffCnt,p->coeffV); //cmVOD_Print(stdout,1,p->coeffCnt,p->coeffV);
@ -872,7 +863,7 @@ cmRC_t cmFIRInitKaiser( cmFIR* p, unsigned procSmpCnt, double srate, double pass
} }
cmRC_t cmFIRInitSinc( cmFIR* p, unsigned procSmpCnt, double srate, unsigned sincSmpCnt, double fcHz, unsigned flags ) cmRC_t cmFIRInitSinc( cmFIR* p, unsigned procSmpCnt, double srate, unsigned sincSmpCnt, double fcHz, unsigned flags, const double* wndV )
{ {
cmRC_t rc; cmRC_t rc;
@ -886,7 +877,11 @@ cmRC_t cmFIRInitSinc( cmFIR* p, unsigned procSmpCnt, double srate, unsigned sinc
p->delayV = cmMemResizeZ( double, p->delayV, p->coeffCnt-1 ); // there is always one less delay than coeff p->delayV = cmMemResizeZ( double, p->delayV, p->coeffCnt-1 ); // there is always one less delay than coeff
p->delayIdx = 0; p->delayIdx = 0;
cmVOD_LP_Sinc(p->coeffV, p->coeffCnt, srate, fcHz, cmIsFlag(flags,kHighPassFIRFl) ? kHighPass_LPSincFl : 0 ); unsigned lp_flags = kNormalize_LPSincFl;
lp_flags |= cmIsFlag(flags,kHighPassFIRFl) ? kHighPass_LPSincFl : 0;
cmVOD_LP_Sinc(p->coeffV, p->coeffCnt, wndV, srate, fcHz, lp_flags );
return cmOkRC; return cmOkRC;
} }
@ -963,7 +958,7 @@ void cmFIRTest( cmRpt_t* rpt, cmLHeapH_t lhH, cmSymTblH_t stH )
cmVOS_Random(in,procSmpCnt, -1.0, 1.0 ); cmVOS_Random(in,procSmpCnt, -1.0, 1.0 );
cmFIR* ffp = cmFIRAllocKaiser( &c, NULL, procSmpCnt,srate, srate*0.025, srate/2, 10, 60 ); cmFIR* ffp = cmFIRAllocKaiser( &c, NULL, procSmpCnt,srate, srate*0.025, srate/2, 10, 60, 0 );
//cmFIR* ffp = cmFIRAllocSinc( &c, NULL, 32, 1000, 0 ); //cmFIR* ffp = cmFIRAllocSinc( &c, NULL, 32, 1000, 0 );
cmFftSR* ftp = cmFftAllocSR( &c, NULL, ffp->outV, ffp->outN, kToPolarFftFl ); cmFftSR* ftp = cmFftAllocSR( &c, NULL, ffp->outV, ffp->outN, kToPolarFftFl );
@ -3661,6 +3656,22 @@ cmRC_t cmNmfExec( cmNmf_t* p, const cmReal_t* vM, unsigned cn )
//------------------------------------------------------------------------------------------------------------ //------------------------------------------------------------------------------------------------------------
unsigned _cmVectArrayTypeByteCnt( cmVectArray_t* p, unsigned flags )
{
switch( flags & kVaMask )
{
case kFloatVaFl: return sizeof(float);
case kDoubleVaFl: return sizeof(double);
case kIntVaFl: return sizeof(int);
case kUIntVaFl: return sizeof(unsigned);
}
if( p != NULL )
cmCtxRtCondition(&p->obj,cmInvalidArgRC,"Unknown data type.");
return 0;
}
cmRC_t _cmVectArrayAppend( cmVectArray_t* p, const void* v, unsigned typeByteCnt, unsigned valCnt ) cmRC_t _cmVectArrayAppend( cmVectArray_t* p, const void* v, unsigned typeByteCnt, unsigned valCnt )
{ {
cmRC_t rc = cmSubSysFailRC; cmRC_t rc = cmSubSysFailRC;
@ -3671,7 +3682,7 @@ cmRC_t _cmVectArrayAppend( cmVectArray_t* p, const void* v, unsigned typeByteCnt
return rc; return rc;
// verify that all vectors written to this vector array contain the same data type. // verify that all vectors written to this vector array contain the same data type.
if( (cmIsFlag(p->flags,kFloatVaFl) && typeByteCnt!=sizeof(float)) || (cmIsFlag(p->flags,kDoubleVaFl) && typeByteCnt!=sizeof(double))) if( typeByteCnt != _cmVectArrayTypeByteCnt(p,p->flags) )
return cmCtxRtCondition(&p->obj,cmInvalidArgRC,"All data stored to a cmVectArray_t must be a consistent type."); return cmCtxRtCondition(&p->obj,cmInvalidArgRC,"All data stored to a cmVectArray_t must be a consistent type.");
// allocate space for the link record // allocate space for the link record
@ -3726,8 +3737,18 @@ cmVectArray_t* cmVectArrayAlloc( cmCtx* ctx, unsigned flags )
assert(p != NULL); assert(p != NULL);
switch( flags & (kFloatVaFl | kDoubleVaFl | kSampleVaFl | kRealVaFl ) ) switch( flags & kVaMask )
{ {
case kIntVaFl:
p->flags |= kIntVaFl;
p->typeByteCnt = sizeof(int);
break;
case kUIntVaFl:
p->flags |= kUIntVaFl;
p->typeByteCnt = sizeof(unsigned);
break;
case kFloatVaFl: case kFloatVaFl:
p->flags |= kFloatVaFl; p->flags |= kFloatVaFl;
p->typeByteCnt = sizeof(float); p->typeByteCnt = sizeof(float);
@ -3738,38 +3759,10 @@ cmVectArray_t* cmVectArrayAlloc( cmCtx* ctx, unsigned flags )
p->typeByteCnt = sizeof(double); p->typeByteCnt = sizeof(double);
break; break;
case kSampleVaFl:
if( sizeof(cmSample_t) == sizeof(float) )
p->flags |= kFloatVaFl;
else
{
if( sizeof(cmSample_t) == sizeof(double))
p->flags |= kDoubleVaFl;
else
rc = cmCtxRtCondition(&p->obj,cmInvalidArgRC,"The size of the cmSample_t is not consistent with float or double.");
}
p->typeByteCnt = sizeof(cmSample_t);
break;
case kRealVaFl:
if( sizeof(cmReal_t) == sizeof(float) )
p->flags |= kFloatVaFl;
else
{
if( sizeof(cmReal_t) == sizeof(double))
p->flags |= kDoubleVaFl;
else
rc = cmCtxRtCondition(&p->obj,cmInvalidArgRC,"The size of the cmReal_t is not consistent with float or double.");
}
p->typeByteCnt = sizeof(cmReal_t);
break;
default: default:
rc = cmCtxRtCondition(&p->obj,cmInvalidArgRC,"The vector array value type flag was not recognized."); rc = cmCtxRtCondition(&p->obj,cmInvalidArgRC,"The vector array value type flag was not recognized.");
} }
if(rc != cmOkRC) if(rc != cmOkRC)
cmVectArrayFree(&p); cmVectArrayFree(&p);
@ -3807,7 +3800,6 @@ cmVectArray_t* cmVectArrayAllocFromFile(cmCtx* ctx, const char* fn )
buf = cmMemAlloc(char,maxEleCnt*typeByteCnt); buf = cmMemAlloc(char,maxEleCnt*typeByteCnt);
if((p = cmVectArrayAlloc(ctx, flags )) == NULL ) if((p = cmVectArrayAlloc(ctx, flags )) == NULL )
goto errLabel; goto errLabel;
@ -3833,10 +3825,8 @@ cmVectArray_t* cmVectArrayAllocFromFile(cmCtx* ctx, const char* fn )
rc = cmCtxRtCondition(&p->obj,rc,"The vector array data store failed on vector index:%i in '%s'.",i,cmStringNullGuard(fn)); rc = cmCtxRtCondition(&p->obj,rc,"The vector array data store failed on vector index:%i in '%s'.",i,cmStringNullGuard(fn));
goto errLabel; goto errLabel;
} }
} }
errLabel: errLabel:
if( fp != NULL ) if( fp != NULL )
@ -3851,7 +3841,6 @@ cmVectArray_t* cmVectArrayAllocFromFile(cmCtx* ctx, const char* fn )
} }
cmRC_t cmVectArrayFree( cmVectArray_t** pp ) cmRC_t cmVectArrayFree( cmVectArray_t** pp )
{ {
cmRC_t rc = cmOkRC; cmRC_t rc = cmOkRC;
@ -3894,19 +3883,22 @@ cmRC_t cmVectArrayClear( cmVectArray_t* p )
unsigned cmVectArrayCount( const cmVectArray_t* p ) unsigned cmVectArrayCount( const cmVectArray_t* p )
{ return p->vectCnt; } { return p->vectCnt; }
/* unsigned cmVectArrayMaxRowCount( const cmVectArray_t* p )
unsigned cmVectArrayEleCount( const cmVectArray_t* p )
{ {
cmVectArrayVect_t* ep = p->bp; const cmVectArrayVect_t* np = p->bp;
unsigned n = 0; unsigned maxN = 0;
for(; ep!=NULL; ep=ep->link )
n += ep->n;
return n; for(; np!=NULL; np=np->link)
if( np->n > maxN )
maxN = np->n;
return maxN;
} }
*/
cmRC_t cmVectArrayAppendV( cmVectArray_t* p, const void* v, unsigned vn )
{ return _cmVectArrayAppend(p,v,_cmVectArrayTypeByteCnt(p,p->flags), vn); }
cmRC_t cmVectArrayAppendS( cmVectArray_t* p, const cmSample_t* v, unsigned vn ) cmRC_t cmVectArrayAppendS( cmVectArray_t* p, const cmSample_t* v, unsigned vn )
{ return _cmVectArrayAppend(p,v,sizeof(v[0]),vn); } { return _cmVectArrayAppend(p,v,sizeof(v[0]),vn); }
@ -3920,27 +3912,10 @@ cmRC_t cmVectArrayAppendD( cmVectArray_t* p, const double* v, unsigned vn )
{ return _cmVectArrayAppend(p,v,sizeof(v[0]),vn); } { return _cmVectArrayAppend(p,v,sizeof(v[0]),vn); }
cmRC_t cmVectArrayAppendI( cmVectArray_t* p, const int* v, unsigned vn ) cmRC_t cmVectArrayAppendI( cmVectArray_t* p, const int* v, unsigned vn )
{ { return _cmVectArrayAppend(p,v,sizeof(v[0]),vn); }
p->tempV = cmMemResize(double,p->tempV,vn);
unsigned i;
for(i=0; i<vn; ++i)
p->tempV[i] = v[i];
cmRC_t rc = cmVectArrayAppendD(p,p->tempV,vn);
return rc;
}
cmRC_t cmVectArrayAppendU( cmVectArray_t* p, const unsigned* v, unsigned vn ) cmRC_t cmVectArrayAppendU( cmVectArray_t* p, const unsigned* v, unsigned vn )
{ { return _cmVectArrayAppend(p,v,sizeof(v[0]),vn); }
p->tempV = cmMemResize(double,p->tempV,vn);
unsigned i;
for(i=0; i<vn; ++i)
p->tempV[i] = v[i];
cmRC_t rc = cmVectArrayAppendD(p,p->tempV,vn);
return rc;
}
cmRC_t cmVectArrayWrite( cmVectArray_t* p, const char* fn ) cmRC_t cmVectArrayWrite( cmVectArray_t* p, const char* fn )
{ {
@ -3985,13 +3960,47 @@ cmRC_t cmVectArrayWrite( cmVectArray_t* p, const char* fn )
} }
} }
errLabel: errLabel:
if( fp != NULL ) if( fp != NULL )
fclose(fp); fclose(fp);
return rc; return rc;
} }
cmRC_t cmVectArrayPrint( cmVectArray_t* p, cmRpt_t* rpt )
{
cmRC_t rc = cmOkRC;
cmVectArrayVect_t* rp = p->bp;
for(; rp!=NULL; rp=rp->link)
{
switch( p->flags & kVaMask )
{
case kFloatVaFl:
cmVOF_Print(rpt,1,rp->n,rp->u.fV);
break;
case kDoubleVaFl:
cmVOD_Print(rpt,1,rp->n,rp->u.dV);
break;
case kIntVaFl:
cmVOI_Print(rpt,1,rp->n,rp->u.iV);
break;
case kUIntVaFl:
cmVOU_Print(rpt,1,rp->n,rp->u.uV);
break;
default:
rc = cmCtxRtCondition(&p->obj,cmInvalidArgRC,"The vector array value type flag was not recognized.");
break;
}
}
return rc;
}
unsigned cmVectArrayForEachS( cmVectArray_t* p, unsigned idx, unsigned cnt, cmVectArrayForEachFuncS_t func, void* arg ) unsigned cmVectArrayForEachS( cmVectArray_t* p, unsigned idx, unsigned cnt, cmVectArrayForEachFuncS_t func, void* arg )
{ {
cmVectArrayVect_t* ep = p->bp; cmVectArrayVect_t* ep = p->bp;
@ -4022,145 +4031,180 @@ unsigned cmVectArrayForEachS( cmVectArray_t* p, unsigned idx, unsigned cnt, cmVe
return n; return n;
} }
cmRC_t _cmVectArrayWriteMatrix( cmCtx* ctx, const char* fn, unsigned flags, const void* m, unsigned rn, unsigned cn )
{
cmRC_t rc = cmOkRC;
cmVectArray_t* p;
const char* b = (const char*)m;
unsigned tbc = _cmVectArrayTypeByteCnt( NULL, flags );
unsigned ri = 0;
char* vv = cmMemAlloc(char,cn*tbc);
if((p = cmVectArrayAlloc(ctx,flags)) == NULL )
return cmCtxRtCondition(&ctx->obj,cmSubSysFailRC,"Unable to allocate a cmVectArray_t in %s().",__FUNCTION__);
for(ri=0; ri<rn; ++ri)
{
// get ptr to first element in row 'ri' or m[]
const char* v = b + ri*tbc;
unsigned ci;
// for each column in m[ri,:]
for(ci=0; ci<cn; ++ci)
memcpy(vv + ci*tbc, v + ci*rn*tbc, tbc );
// append the row to the VectArray
if((rc = cmVectArrayAppendV(p,v,cn)) != cmOkRC )
{
rc = cmCtxRtCondition(&p->obj,rc,"Vector append failed in %s().",__FUNCTION__);
goto errLabel;
}
}
if((rc = cmVectArrayWrite(p,fn)) != cmOkRC )
rc = cmCtxRtCondition(&p->obj,rc,"Vector array write failed in %s().",__FUNCTION__);
errLabel:
if((rc = cmVectArrayFree(&p)) != cmOkRC )
rc = cmCtxRtCondition(&ctx->obj,rc,"Vector array free failed in %s().",__FUNCTION__);
cmMemFree(vv);
return rc;
}
cmRC_t cmVectArrayWriteVectorV( cmCtx* ctx, const char* fn, const void* v, unsigned vn, unsigned flags )
{ return _cmVectArrayWriteMatrix( ctx, fn, flags, v, 1, vn ); }
cmRC_t cmVectArrayWriteVectorS( cmCtx* ctx, const char* fn, const cmSample_t* v, unsigned vn ) cmRC_t cmVectArrayWriteVectorS( cmCtx* ctx, const char* fn, const cmSample_t* v, unsigned vn )
{ { return _cmVectArrayWriteMatrix( ctx, fn, kSampleVaFl, v, 1, vn ); }
cmRC_t rc = cmOkRC;
cmVectArray_t* p;
if((p = cmVectArrayAlloc(ctx,kSampleVaFl)) == NULL ) cmRC_t cmVectArrayWriteVectorR( cmCtx* ctx, const char* fn, const cmReal_t* v, unsigned vn )
return cmCtxRtCondition(&ctx->obj,cmSubSysFailRC,"Unable to allocate an cmVectArray_t in cmVectArrayWriteVectorS()."); { return _cmVectArrayWriteMatrix( ctx, fn, kRealVaFl, v, 1, vn ); }
if((rc = cmVectArrayAppendS(p,v,vn)) != cmOkRC ) cmRC_t cmVectArrayWriteVectorD( cmCtx* ctx, const char* fn, const double* v, unsigned vn )
{ { return _cmVectArrayWriteMatrix( ctx, fn, kDoubleVaFl, v, 1, vn ); }
rc = cmCtxRtCondition(&p->obj,rc,"Vector append failed in cmVectArrayWriteVectorS().");
goto errLabel;
}
if((rc = cmVectArrayWrite(p,fn)) != cmOkRC ) cmRC_t cmVectArrayWriteVectorF( cmCtx* ctx, const char* fn, const float* v, unsigned vn )
rc = cmCtxRtCondition(&p->obj,rc,"Vector array write failed in cmVectArrayWriteVectorS()."); { return _cmVectArrayWriteMatrix( ctx, fn, kFloatVaFl, v, 1, vn ); }
errLabel:
if((rc = cmVectArrayFree(&p)) != cmOkRC )
rc = cmCtxRtCondition(&ctx->obj,rc,"Free failed on cmVectArrayFree() in cmVectArrayWriteVectorS().");
return rc;
}
cmRC_t cmVectArrayWriteMatrixS( cmCtx* ctx, const char* fn, const cmSample_t* m, unsigned rn, unsigned cn )
{
cmRC_t rc = cmOkRC;
cmVectArray_t* p;
unsigned i;
if((p = cmVectArrayAlloc(ctx,kSampleVaFl)) == NULL )
return cmCtxRtCondition(&ctx->obj,cmSubSysFailRC,"Unable to allocate an cmVectArray_t in cmVectArrayWriteVectorS().");
for(i=0; i<cn; ++i)
{
if((rc = cmVectArrayAppendS(p,m + (i*rn), rn)) != cmOkRC )
{
rc = cmCtxRtCondition(&p->obj,rc,"Vector append failed in cmVectArrayWriteVectorS().");
goto errLabel;
}
}
if((rc = cmVectArrayWrite(p,fn)) != cmOkRC )
rc = cmCtxRtCondition(&p->obj,rc,"Vector array write failed in cmVectArrayWriteVectorS().");
errLabel:
if((rc = cmVectArrayFree(&p)) != cmOkRC )
rc = cmCtxRtCondition(&ctx->obj,rc,"Free failed on cmVectArrayFree() in cmVectArrayWriteVectorS().");
return rc;
}
cmRC_t cmVectArrayWriteMatrixR( cmCtx* ctx, const char* fn, const cmReal_t* m, unsigned rn, unsigned cn )
{
cmRC_t rc = cmOkRC;
cmVectArray_t* p;
unsigned i;
if((p = cmVectArrayAlloc(ctx,kRealVaFl)) == NULL )
return cmCtxRtCondition(&ctx->obj,cmSubSysFailRC,"Unable to allocate an cmVectArray_t in cmVectArrayWriteVectorR().");
for(i=0; i<cn; ++i)
{
if((rc = cmVectArrayAppendR(p,m + (i*rn), rn)) != cmOkRC )
{
rc = cmCtxRtCondition(&p->obj,rc,"Vector append failed in cmVectArrayWriteVectorR().");
goto errLabel;
}
}
if((rc = cmVectArrayWrite(p,fn)) != cmOkRC )
rc = cmCtxRtCondition(&p->obj,rc,"Vector array write failed in cmVectArrayWriteVectorR().");
errLabel:
if((rc = cmVectArrayFree(&p)) != cmOkRC )
rc = cmCtxRtCondition(&ctx->obj,rc,"Free failed on cmVectArrayFree() in cmVectArrayWriteVectorR().");
return rc;
}
cmRC_t cmVectArrayWriteVectorI( cmCtx* ctx, const char* fn, const int* v, unsigned vn ) cmRC_t cmVectArrayWriteVectorI( cmCtx* ctx, const char* fn, const int* v, unsigned vn )
{ return _cmVectArrayWriteMatrix( ctx, fn, kIntVaFl, v, 1, vn ); }
cmRC_t cmVectArrayWriteVectorU( cmCtx* ctx, const char* fn, const unsigned* v, unsigned vn )
{ return _cmVectArrayWriteMatrix( ctx, fn, kUIntVaFl, v, 1, vn ); }
cmRC_t cmVectArrayWriteMatrixV( cmCtx* ctx, const char* fn, const void* v, unsigned rn, unsigned cn, unsigned flags )
{ return _cmVectArrayWriteMatrix( ctx, fn, flags, v, rn, cn); }
cmRC_t cmVectArrayWriteMatrixS( cmCtx* ctx, const char* fn, const cmSample_t* v, unsigned rn, unsigned cn )
{ return _cmVectArrayWriteMatrix( ctx, fn, kSampleVaFl, v, rn, cn); }
cmRC_t cmVectArrayWriteMatrixR( cmCtx* ctx, const char* fn, const cmReal_t* v, unsigned rn, unsigned cn )
{ return _cmVectArrayWriteMatrix( ctx, fn, kRealVaFl, v, rn, cn); }
cmRC_t cmVectArrayWriteMatrixD( cmCtx* ctx, const char* fn, const double* v, unsigned rn, unsigned cn )
{ return _cmVectArrayWriteMatrix( ctx, fn, kDoubleVaFl, v, rn, cn); }
cmRC_t cmVectArrayWriteMatrixF( cmCtx* ctx, const char* fn, const float* v, unsigned rn, unsigned cn )
{ return _cmVectArrayWriteMatrix( ctx, fn, kFloatVaFl, v, rn, cn); }
cmRC_t cmVectArrayWriteMatrixI( cmCtx* ctx, const char* fn, const int* v, unsigned rn, unsigned cn )
{ return _cmVectArrayWriteMatrix( ctx, fn, kIntVaFl, v, rn, cn); }
cmRC_t cmVectArrayWriteMatrixU( cmCtx* ctx, const char* fn, const unsigned* v, unsigned rn, unsigned cn )
{ return _cmVectArrayWriteMatrix( ctx, fn, kUIntVaFl, v, rn, cn); }
// Fill v[(*vnRef)*tbc] with the data from the current row of p.
// Return the count of elements copied to v[] in *vnRef.
cmRC_t _cmVectArrayGetV( cmVectArray_t* p, void* v, unsigned* vnRef, unsigned tbc )
{ {
assert( tbc == p->typeByteCnt );
if( cmVectArrayIsEOL(p) )
return cmCtxRtCondition(&p->obj,cmSubSysFailRC,"%s failed because the state is EOL.",__FUNCTION__);
unsigned n = cmMin((*vnRef)*tbc, p->cur->n * p->typeByteCnt );
memcpy(v, p->cur->u.v, n );
*vnRef = n/tbc;
return cmOkRC;
}
cmRC_t _cmVectArrayReadMatrixV( cmCtx* ctx, const char* fn, void** mRef, unsigned* rnRef, unsigned* cnRef )
{
assert( mRef != NULL );
assert( cnRef != NULL );
assert( rnRef != NULL );
*mRef = NULL;
*cnRef = 0;
*rnRef = 0;
cmRC_t rc = cmOkRC; cmRC_t rc = cmOkRC;
cmVectArray_t* p; cmVectArray_t* va;
if((p = cmVectArrayAlloc(ctx,kIntVaFl)) == NULL ) if((va = cmVectArrayAllocFromFile(ctx, fn )) == NULL )
return cmCtxRtCondition(&ctx->obj,cmSubSysFailRC,"Unable to allocate an cmVectArray_t in cmVectArrayWriteVectorS()."); rc = cmCtxRtCondition(&ctx->obj,cmSubSysFailRC,"Unable to read the vectarray from the file '%s'.", cmStringNullGuard(fn));
else
if((rc = cmVectArrayAppendI(p,v,vn)) != cmOkRC )
{ {
rc = cmCtxRtCondition(&p->obj,rc,"Vector append failed in cmVectArrayWriteVectorS()."); unsigned rn = cmVectArrayCount(va); // count of rows
unsigned cn = cmVectArrayMaxRowCount(va); // max count of ele's among all rows
char* m = cmMemAllocZ(char,va->typeByteCnt*rn*cn); // allocate the matrix
unsigned ci = 0;
cmVectArrayRewind(va);
// read each vector into a column of m[]
for(; !cmVectArrayIsEOL(va); ++ci)
{
unsigned n = cmVectArrayEleCount(va);
assert( m+(ci*rn+n)*va->typeByteCnt <= m + rn*cn*va->typeByteCnt );
if( _cmVectArrayGetV(va, m + ci*rn*va->typeByteCnt, &n, va->typeByteCnt) != cmOkRC )
goto errLabel; goto errLabel;
cmVectArrayAdvance(va,1);
} }
if((rc = cmVectArrayWrite(p,fn)) != cmOkRC ) *mRef = m;
rc = cmCtxRtCondition(&p->obj,rc,"Vector array write failed in cmVectArrayWriteVectorS()."); *cnRef = cn;
*rnRef = rn;
}
errLabel: errLabel:
if((rc = cmVectArrayFree(&p)) != cmOkRC ) if( va != NULL )
rc = cmCtxRtCondition(&ctx->obj,rc,"Free failed on cmVectArrayFree() in cmVectArrayWriteVectorS()."); cmVectArrayFree(&va);
return rc; return rc;
} }
cmRC_t cmVectArrayReadMatrixV( cmCtx* ctx, const char* fn, void** mRef, unsigned* rnRef, unsigned* cnRef )
{ return _cmVectArrayReadMatrixV(ctx, fn, mRef, rnRef, cnRef ); }
cmRC_t cmVectArrayWriteMatrixI( cmCtx* ctx, const char* fn, const int* m, unsigned rn, unsigned cn ) cmRC_t cmVectArrayReadMatrixS( cmCtx* ctx, const char* fn, cmSample_t** mRef, unsigned* rnRef, unsigned* cnRef )
{ { return _cmVectArrayReadMatrixV(ctx, fn, (void**)mRef, rnRef, cnRef ); }
cmRC_t rc = cmOkRC;
cmVectArray_t* p;
unsigned i;
if((p = cmVectArrayAlloc(ctx,kIntVaFl)) == NULL ) cmRC_t cmVectArrayReadMatrixR( cmCtx* ctx, const char* fn, cmReal_t** mRef, unsigned* rnRef, unsigned* cnRef )
return cmCtxRtCondition(&ctx->obj,cmSubSysFailRC,"Unable to allocate an cmVectArray_t in cmVectArrayWriteVectorI()."); { return _cmVectArrayReadMatrixV(ctx, fn, (void**)mRef, rnRef, cnRef ); }
for(i=0; i<cn; ++i) cmRC_t cmVectArrayReadMatrixD( cmCtx* ctx, const char* fn, double** mRef, unsigned* rnRef, unsigned* cnRef )
{ { return _cmVectArrayReadMatrixV(ctx, fn, (void**)mRef, rnRef, cnRef ); }
if((rc = cmVectArrayAppendI(p,m + (i*rn), rn)) != cmOkRC )
{
rc = cmCtxRtCondition(&p->obj,rc,"Vector append failed in cmVectArrayWriteVectorI().");
goto errLabel;
}
} cmRC_t cmVectArrayReadMatrixF( cmCtx* ctx, const char* fn, float** mRef, unsigned* rnRef, unsigned* cnRef )
{ return _cmVectArrayReadMatrixV(ctx, fn, (void**)mRef, rnRef, cnRef ); }
if((rc = cmVectArrayWrite(p,fn)) != cmOkRC ) cmRC_t cmVectArrayReadMatrixI( cmCtx* ctx, const char* fn, int** mRef, unsigned* rnRef, unsigned* cnRef )
rc = cmCtxRtCondition(&p->obj,rc,"Vector array write failed in cmVectArrayWriteVectorI()."); { return _cmVectArrayReadMatrixV(ctx, fn, (void**)mRef, rnRef, cnRef ); }
errLabel:
if((rc = cmVectArrayFree(&p)) != cmOkRC )
rc = cmCtxRtCondition(&ctx->obj,rc,"Free failed on cmVectArrayFree() in cmVectArrayWriteVectorI().");
return rc;
}
cmRC_t cmVectArrayReadMatrixU( cmCtx* ctx, const char* fn, unsigned** mRef, unsigned* rnRef, unsigned* cnRef )
{ return _cmVectArrayReadMatrixV(ctx, fn, (void**)mRef, rnRef, cnRef ); }
cmRC_t cmVectArrayForEachTextFuncS( void* arg, unsigned idx, const cmSample_t* xV, unsigned xN ) cmRC_t cmVectArrayForEachTextFuncS( void* arg, unsigned idx, const cmSample_t* xV, unsigned xN )
{ {
@ -4187,13 +4231,11 @@ cmRC_t cmVectArrayAdvance( cmVectArray_t* p, unsigned n )
} }
return cmOkRC; return cmOkRC;
} }
bool cmVectArrayIsEOL( const cmVectArray_t* p ) bool cmVectArrayIsEOL( const cmVectArray_t* p )
{ return p->cur == NULL; } { return p->cur == NULL; }
unsigned cmVectArrayEleCount( const cmVectArray_t* p ) unsigned cmVectArrayEleCount( const cmVectArray_t* p )
{ {
if( p->cur == NULL ) if( p->cur == NULL )
@ -4201,102 +4243,137 @@ unsigned cmVectArrayEleCount( const cmVectArray_t* p )
return p->cur->n; return p->cur->n;
} }
cmRC_t cmVectArrayGetV( cmVectArray_t* p, void* v, unsigned* vnRef )
{ return _cmVectArrayGetV(p,v,vnRef,_cmVectArrayTypeByteCnt(p,p->flags)); }
cmRC_t cmVectArrayGetS( cmVectArray_t* p, cmSample_t* v, unsigned* vnRef )
cmRC_t cmVectArrayGetF( cmVectArray_t* p, float* v, unsigned* vnRef )
{ {
assert( cmIsFlag(p->flags,kFloatVaFl) ); assert( cmIsFlag(p->flags,kSampleVaFl) );
return _cmVectArrayGetV(p,v,vnRef,sizeof(cmSample_t));
}
if( cmVectArrayIsEOL(p) ) cmRC_t cmVectArrayGetR( cmVectArray_t* p, cmReal_t* v, unsigned* vnRef )
return cmCtxRtCondition(&p->obj,cmSubSysFailRC,"VectArray get failed because the state is EOL."); {
assert( cmIsFlag(p->flags,kRealVaFl) );
unsigned n = cmMin(*vnRef,p->cur->n); return _cmVectArrayGetV(p,v,vnRef,sizeof(cmReal_t));
cmVOF_Copy(v,n,p->cur->u.fV);
*vnRef = n;
return cmOkRC;
} }
cmRC_t cmVectArrayGetD( cmVectArray_t* p, double* v, unsigned* vnRef ) cmRC_t cmVectArrayGetD( cmVectArray_t* p, double* v, unsigned* vnRef )
{ {
assert( cmIsFlag(p->flags,kDoubleVaFl) ); assert( cmIsFlag(p->flags,kDoubleVaFl) );
return _cmVectArrayGetV(p,v,vnRef,sizeof(double));
}
if( cmVectArrayIsEOL(p) ) cmRC_t cmVectArrayGetF( cmVectArray_t* p, float* v, unsigned* vnRef )
return cmCtxRtCondition(&p->obj,cmSubSysFailRC,"VectArray get failed because the state is EOL."); {
assert( cmIsFlag(p->flags,kFloatVaFl) );
unsigned n = cmMin(*vnRef,p->cur->n); return _cmVectArrayGetV(p,v,vnRef,sizeof(float));
cmVOD_Copy(v,n,p->cur->u.dV);
*vnRef = n;
return cmOkRC;
} }
cmRC_t cmVectArrayGetI( cmVectArray_t* p, int* v, unsigned* vnRef ) cmRC_t cmVectArrayGetI( cmVectArray_t* p, int* v, unsigned* vnRef )
{ {
if( cmVectArrayIsEOL(p) ) assert( cmIsFlag(p->flags,kIntVaFl) );
return cmCtxRtCondition(&p->obj,cmSubSysFailRC,"VectArray get failed because the state is EOL."); return _cmVectArrayGetV(p,v,vnRef,sizeof(int));
unsigned n = cmMin(*vnRef,p->cur->n);
unsigned i;
if( cmIsFlag(p->flags,kDoubleVaFl) )
{
for(i=0; i<n; ++i)
v[i] = (int)(p->cur->u.dV[i]);
}
else
{
if( cmIsFlag(p->flags,kFloatVaFl) )
{
for(i=0; i<n; ++i)
v[i] = (int)(p->cur->u.fV[i]);
}
else
{
assert(0);
}
}
*vnRef = n;
return cmOkRC;
} }
cmRC_t cmVectArrayGetU( cmVectArray_t* p, unsigned* v, unsigned* vnRef ) cmRC_t cmVectArrayGetU( cmVectArray_t* p, unsigned* v, unsigned* vnRef )
{ {
if( cmVectArrayIsEOL(p) ) assert( cmIsFlag(p->flags,kUIntVaFl) );
return cmCtxRtCondition(&p->obj,cmSubSysFailRC,"VectArray get failed because the state is EOL."); return _cmVectArrayGetV(p,v,vnRef,sizeof(unsigned));
unsigned n = cmMin(*vnRef,p->cur->n);
unsigned i;
if( cmIsFlag(p->flags,kDoubleVaFl) )
{
for(i=0; i<n; ++i)
v[i] = (unsigned)(p->cur->u.dV[i]);
}
else
{
if( cmIsFlag(p->flags,kFloatVaFl) )
{
for(i=0; i<n; ++i)
v[i] = (unsigned)(p->cur->u.fV[i]);
}
else
{
assert(0);
}
}
*vnRef = n;
return cmOkRC;
} }
cmRC_t _cmVectArrayMatrixIsEqual( cmCtx* ctx, const char* fn, const void* mm, unsigned rn, unsigned cn, unsigned flags, bool* resultFlRef )
{
assert( resultFlRef != NULL );
cmRC_t rc = cmOkRC;
cmVectArray_t* p = NULL;
const char* m = (const char*)mm;
unsigned tbc = _cmVectArrayTypeByteCnt(NULL,flags);
unsigned ri = 0;
char* vv = cmMemAlloc(char,cn*tbc);
*resultFlRef = false;
// read the vector array
if((p = cmVectArrayAllocFromFile(ctx, fn )) == NULL)
{
rc = cmCtxRtCondition(&ctx->obj,cmSubSysFailRC,"Unable to read the VectArray from the file '%s'.", cmStringNullGuard(fn));
goto errLabel;
}
// verify that the matrix type matches the vector array type
if( (p->flags & kVaMask) != (flags & kVaMask) )
{
rc = cmCtxRtCondition(&ctx->obj,cmInvalidArgRC,"Invalid type conversion in '%s'.",__FUNCTION__);
goto errLabel;
}
// the row count of the VectArray and m[] must be the same
if( cmVectArrayCount(p) != rn )
{
*resultFlRef = false;
goto errLabel;
}
// for each row in VectArray
for(; !cmVectArrayIsEOL(p); ++ri )
{
unsigned vn = cmVectArrayEleCount(p);
char v[ vn*p->typeByteCnt ];
unsigned ci;
// get the current row from the VectArray into v[vn]
if( _cmVectArrayGetV(p,v,&vn,p->typeByteCnt) != cmOkRC )
goto errLabel;
// if the size of the current row does not match the row element count of the matrix
if( vn != cn )
goto errLabel;
for(ci=0; ci<cn; ++ci)
memcpy(vv + ci*tbc, m + (ri*tbc) + (ci*rn*tbc), tbc );
// the current row does not match the matrix column vector
if( memcmp(v, vv, vn*p->typeByteCnt ) != 0 )
goto errLabel;
cmVectArrayAdvance(p,1);
}
*resultFlRef = true;
errLabel:
if( p != NULL )
cmVectArrayFree(&p);
cmMemFree(vv);
return rc;
}
cmRC_t cmVectArrayMatrixIsEqualV( cmCtx* ctx, const char* fn, const void* m, unsigned rn, unsigned cn, bool* resultFlRef, unsigned flags )
{ return _cmVectArrayMatrixIsEqual(ctx, fn, m, rn, cn, flags, resultFlRef); }
cmRC_t cmVectArrayMatrixIsEqualS( cmCtx* ctx, const char* fn, const cmSample_t* m, unsigned rn, unsigned cn, bool* resultFlRef )
{ return _cmVectArrayMatrixIsEqual(ctx, fn, m, rn, cn, kSampleVaFl, resultFlRef ); }
cmRC_t cmVectArrayMatrixIsEqualR( cmCtx* ctx, const char* fn, const cmReal_t* m, unsigned rn, unsigned cn, bool* resultFlRef )
{ return _cmVectArrayMatrixIsEqual(ctx, fn, m, rn, cn, kRealVaFl, resultFlRef ); }
cmRC_t cmVectArrayMatrixIsEqualD( cmCtx* ctx, const char* fn, const double* m, unsigned rn, unsigned cn, bool* resultFlRef )
{ return _cmVectArrayMatrixIsEqual(ctx, fn, m, rn, cn, kDoubleVaFl, resultFlRef ); }
cmRC_t cmVectArrayMatrixIsEqualF( cmCtx* ctx, const char* fn, const float* m, unsigned rn, unsigned cn, bool* resultFlRef )
{ return _cmVectArrayMatrixIsEqual(ctx, fn, m, rn, cn, kFloatVaFl, resultFlRef ); }
cmRC_t cmVectArrayMatrixIsEqualI( cmCtx* ctx, const char* fn, const int* m, unsigned rn, unsigned cn, bool* resultFlRef )
{ return _cmVectArrayMatrixIsEqual(ctx, fn, m, rn, cn, kIntVaFl, resultFlRef ); }
cmRC_t cmVectArrayMatrixIsEqualU( cmCtx* ctx, const char* fn, const unsigned* m, unsigned rn, unsigned cn, bool* resultFlRef )
{ return _cmVectArrayMatrixIsEqual(ctx, fn, m, rn, cn, kUIntVaFl, resultFlRef ); }
unsigned cmVectArrayVectEleCount( cmVectArray_t* p, unsigned groupIdx, unsigned groupCnt ) unsigned cmVectArrayVectEleCount( cmVectArray_t* p, unsigned groupIdx, unsigned groupCnt )
{ {
unsigned n = 0; unsigned n = 0;
@ -4321,9 +4398,9 @@ unsigned cmVectArrayVectEleCount( cmVectArray_t* p, unsigned groupIdx, unsigned
p->cur = pos; p->cur = pos;
return n; return n;
} }
cmRC_t cmVectArrayFormVectF( cmVectArray_t* p, unsigned groupIdx, unsigned groupCnt, float** vRef, unsigned* vnRef ) cmRC_t cmVectArrayFormVectF( cmVectArray_t* p, unsigned groupIdx, unsigned groupCnt, float** vRef, unsigned* vnRef )
{ {
cmRC_t rc = cmOkRC; cmRC_t rc = cmOkRC;
@ -4421,7 +4498,6 @@ cmRC_t cmVectArrayFormVectColF( cmVectArray_t* p, unsigned groupIdx, unsigned
return rc; return rc;
} }
cmRC_t cmVectArrayFormVectColU( cmVectArray_t* p, unsigned groupIdx, unsigned groupCnt, unsigned colIdx, unsigned** vRef, unsigned* vnRef ) cmRC_t cmVectArrayFormVectColU( cmVectArray_t* p, unsigned groupIdx, unsigned groupCnt, unsigned colIdx, unsigned** vRef, unsigned* vnRef )
{ {
cmRC_t rc = cmOkRC; cmRC_t rc = cmOkRC;
@ -4472,16 +4548,19 @@ cmRC_t cmVectArrayFormVectColU( cmVectArray_t* p, unsigned groupIdx, unsigned
return rc; return rc;
} }
cmRC_t cmVectArrayTest( cmCtx* ctx, const char* fn ) cmRC_t cmVectArrayTest( cmCtx* ctx, const char* fn, bool genFl )
{ {
cmRC_t rc = cmOkRC; cmRC_t rc = cmOkRC;
cmVectArray_t* p = NULL; cmVectArray_t* p = NULL;
unsigned flags = kSampleVaFl;
cmSample_t v[] = { 0, 1, 2, 3, 4, 5 };
if( fn == NULL || strlen(fn)==0 ) if( fn == NULL || strlen(fn)==0 )
return cmCtxRtCondition(&p->obj,cmSubSysFailRC,"Invalid test output file name."); return cmCtxRtCondition(&p->obj,cmSubSysFailRC,"Invalid test output file name.");
if( genFl )
{
unsigned flags = kSampleVaFl;
cmSample_t v[] = { 0, 1, 2, 3, 4, 5 };
if( (p = cmVectArrayAlloc(ctx,flags)) == NULL ) if( (p = cmVectArrayAlloc(ctx,flags)) == NULL )
return cmCtxRtCondition(&p->obj,cmSubSysFailRC,"The vectory array object allocation failed."); return cmCtxRtCondition(&p->obj,cmSubSysFailRC,"The vectory array object allocation failed.");
@ -4511,20 +4590,20 @@ cmRC_t cmVectArrayTest( cmCtx* ctx, const char* fn )
//cmVectArrayForEachS(p,0,cmVectArrayEleCount(p),cmVectArrayForEachTextFuncS,&ctx->printRpt); //cmVectArrayForEachS(p,0,cmVectArrayEleCount(p),cmVectArrayForEachTextFuncS,&ctx->printRpt);
if( cmVectArrayFree(&p) != cmOkRC ) //if( cmVectArrayFree(&p) != cmOkRC )
{ //{
rc = cmCtxRtCondition(&p->obj,cmSubSysFailRC,"The vectory array release failed."); // rc = cmCtxRtCondition(&p->obj,cmSubSysFailRC,"The vectory array release failed.");
goto errLabel; // goto errLabel;
//}
} }
else
{
if((p = cmVectArrayAllocFromFile(ctx, fn )) == NULL ) if((p = cmVectArrayAllocFromFile(ctx, fn )) == NULL )
{ {
rc = cmCtxRtCondition(&p->obj,cmSubSysFailRC,"VectArray alloc from file failed."); rc = cmCtxRtCondition(&p->obj,cmSubSysFailRC,"VectArray alloc from file failed.");
goto errLabel; goto errLabel;
} }
while(!cmVectArrayIsEOL(p)) while(!cmVectArrayIsEOL(p))
{ {
unsigned n = cmVectArrayEleCount(p); unsigned n = cmVectArrayEleCount(p);
@ -4537,15 +4616,28 @@ cmRC_t cmVectArrayTest( cmCtx* ctx, const char* fn )
goto errLabel; goto errLabel;
} }
cmVOS_PrintL("v:",NULL,1,n,v); //cmVOS_PrintL("v:",NULL,1,n,v);
cmVectArrayAdvance(p,1); cmVectArrayAdvance(p,1);
} }
// Test matrix reading
cmSample_t* m;
unsigned rn,cn;
if( cmVectArrayReadMatrixS(ctx, fn, &m, &rn, &cn ) != cmOkRC )
goto errLabel;
else
{
//cmVOS_PrintL("v:",NULL,rn,cn,m);
cmMemFree(m);
}
}
errLabel: errLabel:
if( cmVectArrayFree(&p) != cmOkRC ) if( cmVectArrayFree(&p) != cmOkRC )
rc = cmCtxRtCondition(&p->obj,cmSubSysFailRC,"The vectory array release failed."); rc = cmCtxRtCondition(&p->obj,cmSubSysFailRC,"The vector array release failed.");
return rc; return rc;
} }
@ -4753,6 +4845,7 @@ cmRC_t cmFrqTrkFree( cmFrqTrk** pp )
cmMemFree(p->logFn); cmMemFree(p->logFn);
cmMemFree(p->levelFn); cmMemFree(p->levelFn);
cmMemFree(p->specFn); cmMemFree(p->specFn);
cmMemFree(p->attenFn);
cmObjFree(pp); cmObjFree(pp);
return rc; return rc;

102
cmProc2.h
View File

@ -176,12 +176,18 @@ extern "C" {
enum { kHighPassFIRFl = 0x01 }; enum { kHighPassFIRFl = 0x01 };
/// Set p to NULL to dynamically allocate the object. // Note that the relative values of passHz and stopHz do not matter
cmFIR* cmFIRAllocKaiser(cmCtx* c, cmFIR* p, unsigned procSmpCnt, double srate, double passHz, double stopHz, double passDb, double stopDb ); // for low-pass vs high-pass filters. In practice passHz and
cmFIR* cmFIRAllocSinc( cmCtx* c, cmFIR* p, unsigned procSmpCnt, double srate, unsigned sincSmpCnt, double fcHz, unsigned flags ); // stopHz can be swapped with no effect on the filter in either
// case. Set p to NULL to dynamically allocate the object.
cmFIR* cmFIRAllocKaiser(cmCtx* c, cmFIR* p, unsigned procSmpCnt, double srate, double passHz, double stopHz, double passDb, double stopDb, unsigned flags );
// Set wndV[sincSmpCnt] to NULL to use a unity window otherwise set it to a window
// function of length sincSmpCnt.
cmFIR* cmFIRAllocSinc( cmCtx* c, cmFIR* p, unsigned procSmpCnt, double srate, unsigned sincSmpCnt, double fcHz, unsigned flags, const double* wndV );
cmRC_t cmFIRFree( cmFIR** pp ); cmRC_t cmFIRFree( cmFIR** pp );
cmRC_t cmFIRInitKaiser( cmFIR* p, unsigned procSmpCnt, double srate, double passHz, double stopHz, double passDb, double stopDb ); cmRC_t cmFIRInitKaiser( cmFIR* p, unsigned procSmpCnt, double srate, double passHz, double stopHz, double passDb, double stopDb, unsigned flags );
cmRC_t cmFIRInitSinc( cmFIR* p, unsigned procSmpCnt, double srate, unsigned sincSmpCnt, double fcHz, unsigned flags ); 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 cmFIRFinal( cmFIR* p );
cmRC_t cmFIRExec( cmFIR* p, const cmSample_t* sp, unsigned sn ); cmRC_t cmFIRExec( cmFIR* p, const cmSample_t* sp, unsigned sn );
void cmFIRTest(); void cmFIRTest();
@ -783,6 +789,8 @@ extern "C" {
double* dV; // dV[n] vector of doubles double* dV; // dV[n] vector of doubles
float* fV; // fV[n] vecotr of floats float* fV; // fV[n] vecotr of floats
cmSample_t* sV; // sV[n] vector of cmSample_t cmSample_t* sV; // sV[n] vector of cmSample_t
int* iV;
unsigned* uV;
} u; } u;
struct cmVectArrayVect_str* link; // link to next element record struct cmVectArrayVect_str* link; // link to next element record
@ -791,12 +799,13 @@ extern "C" {
enum enum
{ {
kFloatVaFl = 0x01, kDoubleVaFl = 0x01,
kDoubleVaFl = 0x02, kRealVaFl = 0x01,
kSampleVaFl = 0x04, kFloatVaFl = 0x02,
kRealVaFl = 0x08, kSampleVaFl = 0x02,
kIntVaFl = 0x02, // int and uint is converted to double kIntVaFl = 0x04,
kUIntVaFl = 0x02 // kUIntVaFl = 0x08,
kVaMask = 0x0f
}; };
typedef struct typedef struct
@ -824,7 +833,13 @@ extern "C" {
// Return the count of vectors contained in the vector array. // Return the count of vectors contained in the vector array.
cmRC_t cmVectArrayCount( const cmVectArray_t* p ); cmRC_t cmVectArrayCount( const cmVectArray_t* p );
// Return the maximum element count among all rows.
unsigned cmVectArrayMaxRowCount( const cmVectArray_t* p );
// Store a new vector by appending it to the end of the internal vector list. // 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
// the data type set in p->flags.
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 cmVectArrayAppendS( cmVectArray_t* p, const cmSample_t* v, unsigned vn );
cmRC_t cmVectArrayAppendR( cmVectArray_t* p, const cmReal_t* v, unsigned vn ); cmRC_t cmVectArrayAppendR( cmVectArray_t* p, const cmReal_t* v, unsigned vn );
cmRC_t cmVectArrayAppendF( cmVectArray_t* p, const float* v, unsigned vn ); cmRC_t cmVectArrayAppendF( cmVectArray_t* p, const float* v, unsigned vn );
@ -835,45 +850,84 @@ extern "C" {
// Write a vector array in a format that can be read by readVectArray.m. // 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 cmVectArrayWrite( cmVectArray_t* p, const char* fn );
// Print the vector array to rpt.
cmRC_t cmVectArrayPrint( cmVectArray_t* p, cmRpt_t* rpt );
typedef cmRC_t (*cmVectArrayForEachFuncS_t)( void* arg, unsigned idx, const cmSample_t* xV, unsigned xN ); typedef cmRC_t (*cmVectArrayForEachFuncS_t)( void* arg, unsigned idx, const cmSample_t* xV, unsigned xN );
unsigned cmVectArrayForEachS( cmVectArray_t* p, unsigned idx, unsigned cnt, cmVectArrayForEachFuncS_t func, void* arg ); 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
// data type set in the 'flags' parameter.
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 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 );
cmRC_t cmVectArrayWriteVectorD( cmCtx* ctx, const char* fn, const double* v, unsigned vn );
cmRC_t cmVectArrayWriteVectorF( cmCtx* ctx, const char* fn, const float* v, unsigned vn );
cmRC_t cmVectArrayWriteVectorI( cmCtx* ctx, const char* fn, const int* v, unsigned vn ); cmRC_t cmVectArrayWriteVectorI( cmCtx* ctx, const char* fn, const int* v, unsigned vn );
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 matrix is transposed as it is // Write the column-major matrix m[rn,cn] to the file 'fn'.
// written and therefore will be read back as a 'cn' by 'rn' matrix. // Note that the true type of m[] in cmVectArrayWriteMatrixV() must match the
// data type set in the 'flags' parameter.
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 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 ); cmRC_t cmVectArrayWriteMatrixR( cmCtx* ctx, const char* fn, const cmReal_t* m, unsigned rn, unsigned cn );
cmRC_t cmVectArrayWriteMatrixD( cmCtx* ctx, const char* fn, const double* m, unsigned rn, unsigned cn );
cmRC_t cmVectArrayWriteMatrixF( cmCtx* ctx, const char* fn, const float* m, unsigned rn, unsigned cn );
cmRC_t cmVectArrayWriteMatrixI( cmCtx* ctx, const char* fn, const int* m, unsigned rn, unsigned cn ); cmRC_t cmVectArrayWriteMatrixI( cmCtx* ctx, const char* fn, const int* m, unsigned rn, unsigned cn );
cmRC_t cmVectArrayWriteMatrixU( cmCtx* ctx, const char* fn, const unsigned* m, unsigned rn, unsigned cn );
// Read a VectArray file and return it as a matrix.
// The returned memory must be released with a subsequent call to cmMemFree().
// Note that the true type of the pointer address 'mRef' in the call to
// cmVectArrayReadMatrixV() must match the data type of the cmVectArray_t
// specified by 'fn'.
cmRC_t cmVectArrayReadMatrixV( cmCtx* ctx, const char* fn, void** mRef, unsigned* rnRef, unsigned* cnRef );
cmRC_t cmVectArrayReadMatrixS( cmCtx* ctx, const char* fn, cmSample_t** mRef, unsigned* rnRef, unsigned* cnRef );
cmRC_t cmVectArrayReadMatrixR( cmCtx* ctx, const char* fn, cmReal_t** mRef, unsigned* rnRef, unsigned* cnRef );
cmRC_t cmVectArrayReadMatrixD( cmCtx* ctx, const char* fn, double** mRef, unsigned* rnRef, unsigned* cnRef );
cmRC_t cmVectArrayReadMatrixF( cmCtx* ctx, const char* fn, float** mRef, unsigned* rnRef, unsigned* cnRef );
cmRC_t cmVectArrayReadMatrixI( cmCtx* ctx, const char* fn, int** mRef, unsigned* rnRef, unsigned* cnRef );
cmRC_t cmVectArrayReadMatrixU( cmCtx* ctx, const char* fn, unsigned** mRef, unsigned* rnRef, unsigned* cnRef );
// Row iteration control functions.
cmRC_t cmVectArrayRewind( cmVectArray_t* p ); cmRC_t cmVectArrayRewind( cmVectArray_t* p );
cmRC_t cmVectArrayAdvance( cmVectArray_t* p, unsigned n ); cmRC_t cmVectArrayAdvance( cmVectArray_t* p, unsigned n );
bool cmVectArrayIsEOL( const cmVectArray_t* p ); bool cmVectArrayIsEOL( const cmVectArray_t* p );
unsigned cmVectArrayEleCount( const cmVectArray_t* p ); unsigned cmVectArrayEleCount( const cmVectArray_t* p );
cmRC_t cmVectArrayGetF( cmVectArray_t* p, float* v, unsigned* vnRef );
// Copy the current row vector to v[].
// Note that the true type of v[] in cmVectArrayGetV() must match the data type of 'p'.
cmRC_t cmVectArrayGetV( cmVectArray_t* p, void* v, unsigned* vnRef );
cmRC_t cmVectArrayGetS( cmVectArray_t* p, cmSample_t* v, unsigned* vnRef );
cmRC_t cmVectArrayGetR( cmVectArray_t* p, cmReal_t* v, unsigned* vnRef );
cmRC_t cmVectArrayGetD( cmVectArray_t* p, double* v, unsigned* vnRef ); cmRC_t cmVectArrayGetD( cmVectArray_t* p, double* v, unsigned* vnRef );
cmRC_t cmVectArrayGetF( cmVectArray_t* p, float* v, unsigned* vnRef );
cmRC_t cmVectArrayGetI( cmVectArray_t* p, int* v, unsigned* vnRef ); cmRC_t cmVectArrayGetI( cmVectArray_t* p, int* v, unsigned* vnRef );
cmRC_t cmVectArrayGetU( cmVectArray_t* p, unsigned* v, unsigned* vnRef ); cmRC_t cmVectArrayGetU( cmVectArray_t* p, unsigned* v, unsigned* vnRef );
// Set *resultFlRef to true if m[rn,cn] is equal to the cmVectArray_t specified by 'fn'.
// Note that the true type of 'm[]' in the call to cmVectArrayMatrixIsEqualV()
// must match the data type set in 'flags'.
cmRC_t cmVectArrayMatrixIsEqualV( cmCtx* ctx, const char* fn, const void* m, unsigned rn, unsigned cn, bool* resultFlRef, unsigned flags );
cmRC_t cmVectArrayMatrixIsEqualS( cmCtx* ctx, const char* fn, const cmSample_t* m, unsigned rn, unsigned cn, bool* resultFlRef );
cmRC_t cmVectArrayMatrixIsEqualR( cmCtx* ctx, const char* fn, const cmReal_t* m, unsigned rn, unsigned cn, bool* resultFlRef );
cmRC_t cmVectArrayMatrixIsEqualD( cmCtx* ctx, const char* fn, const double* m, unsigned rn, unsigned cn, bool* resultFlRef );
cmRC_t cmVectArrayMatrixIsEqualF( cmCtx* ctx, const char* fn, const float* m, unsigned rn, unsigned cn, bool* resultFlRef );
cmRC_t cmVectArrayMatrixIsEqualI( cmCtx* ctx, const char* fn, const int* m, unsigned rn, unsigned cn, bool* resultFlRef );
cmRC_t cmVectArrayMatrixIsEqualU( cmCtx* ctx, const char* fn, const unsigned* m, unsigned rn, unsigned cn, bool* resultFlRef );
// If a vector array is composed of repeating blocks of 'groupCnt' sub-vectors // If a vector array is composed of repeating blocks of 'groupCnt' sub-vectors
// where the concatenated ith sub-vectors in each group form a single super-vector then // where the concatenated ith sub-vectors in each group form a single super-vector then
// this function will return the super-vector. Use cmMemFree(*vRef) to release // this function will return the super-vector. Use cmMemFree(*vRef) to release
// the returned super-vector. // the returned super-vector.
cmRC_t cmVectArrayFormVectR( cmVectArray_t* p, unsigned groupIdx, unsigned groupCnt, cmReal_t** vRef, unsigned* vnRef );
cmRC_t cmVectArrayFormVectF( cmVectArray_t* p, unsigned groupIdx, unsigned groupCnt, float** vRef, unsigned* vnRef ); cmRC_t cmVectArrayFormVectF( cmVectArray_t* p, unsigned groupIdx, unsigned groupCnt, float** vRef, unsigned* vnRef );
cmRC_t cmVectArrayFormVectColF( cmVectArray_t* p, unsigned groupIdx, unsigned groupCnt, unsigned colIdx, float** vRef, unsigned* vnRef ); cmRC_t cmVectArrayFormVectColF( cmVectArray_t* p, unsigned groupIdx, unsigned groupCnt, unsigned colIdx, float** vRef, unsigned* vnRef );
cmRC_t cmVectArrayFormVectColU( cmVectArray_t* p, unsigned groupIdx, unsigned groupCnt, unsigned colIdx, unsigned** vRef, unsigned* vnRef ); cmRC_t cmVectArrayFormVectColU( cmVectArray_t* p, unsigned groupIdx, unsigned groupCnt, unsigned colIdx, unsigned** vRef, unsigned* vnRef );
cmRC_t cmVectArrayTest( cmCtx* ctx, const char* fn ); cmRC_t cmVectArrayTest( cmCtx* ctx, const char* fn, bool genFl );
#if CM_FLOAT_SMP == 1
#define cmVectArrayGetS cmVectArrayGetF
#define cmVectArrayFormVectS cmVectArrayFormVectF
#else
#define cmVectArrayGetS cmVectArrayGetD
#define cmVectArrayFormVectS cmVectArrayFormVectD
#endif
//----------------------------------------------------------------------------------------------------------------------- //-----------------------------------------------------------------------------------------------------------------------
// Spectral whitening filter. // Spectral whitening filter.