|
@@ -0,0 +1,949 @@
|
|
|
|
|
+/*******************************************************************************
|
|
|
|
|
+版权所有: 珠海欧力配网自动化股份有限公司
|
|
|
|
|
+文件名称: wavelib.c
|
|
|
|
|
+文件版本: 01.00
|
|
|
|
|
+创建作者: ygl
|
|
|
|
|
+创建日期: 2022-03-09
|
|
|
|
|
+功能说明: modwt小波计算
|
|
|
|
|
+其它说明: 支持:DB,sym、coif、haar
|
|
|
|
|
+注意事项:只求近似值,不考虑数据重构
|
|
|
|
|
+*/
|
|
|
|
|
+#ifdef __KERNEL__
|
|
|
|
|
+
|
|
|
|
|
+#include "head.h"
|
|
|
|
|
+#define M_SQRT2 1.41421356237309504880
|
|
|
|
|
+
|
|
|
|
|
+#else
|
|
|
|
|
+#include "head.h"
|
|
|
|
|
+#include <stdio.h>
|
|
|
|
|
+#include <stdlib.h>
|
|
|
|
|
+#include <string.h>
|
|
|
|
|
+#define _USE_MATH_DEFINES
|
|
|
|
|
+#include <math.h>
|
|
|
|
|
+#include "wavelib.h"
|
|
|
|
|
+
|
|
|
|
|
+#define rt_malloc malloc
|
|
|
|
|
+// #define rt_printf printf
|
|
|
|
|
+#define rt_free free
|
|
|
|
|
+
|
|
|
|
|
+#endif
|
|
|
|
|
+
|
|
|
|
|
+hilbert_para coef1; //四组系数
|
|
|
|
|
+hilbert_para coef2;
|
|
|
|
|
+hilbert_para coef3;
|
|
|
|
|
+hilbert_para coef4;
|
|
|
|
|
+struct file *g_fpuo = NULL,*g_fpio = NULL;
|
|
|
|
|
+
|
|
|
|
|
+static const float db1[2] = {
|
|
|
|
|
+ 7.071067811865475244008443621048490392848359376884740365883398e-01,
|
|
|
|
|
+ 7.071067811865475244008443621048490392848359376884740365883398e-01
|
|
|
|
|
+};
|
|
|
|
|
+
|
|
|
|
|
+static const float db2[4] = {
|
|
|
|
|
+ 4.829629131445341433748715998644486838169524195042022752011715e-01,
|
|
|
|
|
+ 8.365163037378079055752937809168732034593703883484392934953414e-01,
|
|
|
|
|
+ 2.241438680420133810259727622404003554678835181842717613871683e-01,
|
|
|
|
|
+ -1.294095225512603811744494188120241641745344506599652569070016e-01
|
|
|
|
|
+};
|
|
|
|
|
+
|
|
|
|
|
+
|
|
|
|
|
+static const float db3[6] = {
|
|
|
|
|
+ 3.326705529500826159985115891390056300129233992450683597084705e-01,
|
|
|
|
|
+ 8.068915093110925764944936040887134905192973949948236181650920e-01,
|
|
|
|
|
+ 4.598775021184915700951519421476167208081101774314923066433867e-01,
|
|
|
|
|
+ -1.350110200102545886963899066993744805622198452237811919756862e-01,
|
|
|
|
|
+ -8.544127388202666169281916918177331153619763898808662976351748e-02,
|
|
|
|
|
+ 3.522629188570953660274066471551002932775838791743161039893406e-02
|
|
|
|
|
+};
|
|
|
|
|
+
|
|
|
|
|
+
|
|
|
|
|
+static const float db4[8] = {
|
|
|
|
|
+ 2.303778133088965008632911830440708500016152482483092977910968e-01,
|
|
|
|
|
+ 7.148465705529156470899219552739926037076084010993081758450110e-01,
|
|
|
|
|
+ 6.308807679298589078817163383006152202032229226771951174057473e-01,
|
|
|
|
|
+ -2.798376941685985421141374718007538541198732022449175284003358e-02,
|
|
|
|
|
+ -1.870348117190930840795706727890814195845441743745800912057770e-01,
|
|
|
|
|
+ 3.084138183556076362721936253495905017031482172003403341821219e-02,
|
|
|
|
|
+ 3.288301166688519973540751354924438866454194113754971259727278e-02,
|
|
|
|
|
+ -1.059740178506903210488320852402722918109996490637641983484974e-02
|
|
|
|
|
+};
|
|
|
|
|
+
|
|
|
|
|
+
|
|
|
|
|
+// static const float db5[10] = {
|
|
|
|
|
+// 1.601023979741929144807237480204207336505441246250578327725699e-01,
|
|
|
|
|
+// 6.038292697971896705401193065250621075074221631016986987969283e-01,
|
|
|
|
|
+// 7.243085284377729277280712441022186407687562182320073725767335e-01,
|
|
|
|
|
+// 1.384281459013207315053971463390246973141057911739561022694652e-01,
|
|
|
|
|
+// -2.422948870663820318625713794746163619914908080626185983913726e-01,
|
|
|
|
|
+// -3.224486958463837464847975506213492831356498416379847225434268e-02,
|
|
|
|
|
+// 7.757149384004571352313048938860181980623099452012527983210146e-02,
|
|
|
|
|
+// -6.241490212798274274190519112920192970763557165687607323417435e-03,
|
|
|
|
|
+// -1.258075199908199946850973993177579294920459162609785020169232e-02,
|
|
|
|
|
+// 3.335725285473771277998183415817355747636524742305315099706428e-03
|
|
|
|
|
+// };
|
|
|
|
|
+
|
|
|
|
|
+
|
|
|
|
|
+// All coif coefficents have to be multiplied by sqrt(2)
|
|
|
|
|
+static const float coif1[6] = {
|
|
|
|
|
+ -5.142972847076845595317549230122688830344559947132656813651045e-02,
|
|
|
|
|
+ 2.389297284707684559531754923012268883034455994713265681365104e-01,
|
|
|
|
|
+ 6.028594569415369119063509846024537766068911989426531362730209e-01,
|
|
|
|
|
+ 2.721405430584630880936490153975462233931088010573468637269790e-01,
|
|
|
|
|
+ -5.142972847076845595317549230122688830344559947132656813651045e-02,
|
|
|
|
|
+ -1.107027152923154404682450769877311169655440052867343186348954e-02
|
|
|
|
|
+};
|
|
|
|
|
+
|
|
|
|
|
+
|
|
|
|
|
+static const float coif2[12] = {
|
|
|
|
|
+ 1.158759673871686817889714882853120395708315073355502818875931e-02,
|
|
|
|
|
+ -2.932013798346856448679594524397843054053420947418409889774786e-02,
|
|
|
|
|
+ -4.763959031100813225872995081511549408622753909592460525840745e-02,
|
|
|
|
|
+ 2.730210465347666137982239328923516270034828327990699588033501e-01,
|
|
|
|
|
+ 5.746823938568638472459483149751499367740786490481481391460366e-01,
|
|
|
|
|
+ 2.948671936956191896750637208703777973914107635455611537640778e-01,
|
|
|
|
|
+ -5.408560709171142997443672832006888537570221990444706777525838e-02,
|
|
|
|
|
+ -4.202648046077160694657530752545884878978719268926222513485613e-02,
|
|
|
|
|
+ 1.674441016327950635146257083249391698866289538037299820224006e-02,
|
|
|
|
|
+ 3.967883612962012109043447090269950094081810916481648252817197e-03,
|
|
|
|
|
+ -1.289203356140659543141355500990678257894936161704492503370186e-03,
|
|
|
|
|
+ -5.095053991076441489598480835620951586540050976664367876412655e-04
|
|
|
|
|
+};
|
|
|
|
|
+
|
|
|
|
|
+static const float coif3[18] = {
|
|
|
|
|
+ -2.682418670922068664584689955153722375535836177157637134187840e-03,
|
|
|
|
|
+ 5.503126707831385107969640263617469178794666057252906037981936e-03,
|
|
|
|
|
+ 1.658356047917034608134280439996549525220639437145367606178002e-02,
|
|
|
|
|
+ -4.650776447872697640390293095170192691113917841041002855534619e-02,
|
|
|
|
|
+ -4.322076356021191118175840907244577856782537221435748296465882e-02,
|
|
|
|
|
+ 2.865033352736474630249006862976158896891076238443844211133873e-01,
|
|
|
|
|
+ 5.612852568703300445990941995240077241406247774064453800050914e-01,
|
|
|
|
|
+ 3.029835717728241602862575774374668529867757043461413348549577e-01,
|
|
|
|
|
+ -5.077014075488886159516471867138370972545857441670871832472707e-02,
|
|
|
|
|
+ -5.819625076158553022607041679522801089624825903982541419721721e-02,
|
|
|
|
|
+ 2.443409432116695639462954438418928805487699080947974989338820e-02,
|
|
|
|
|
+ 1.122924096203786563399489540091488781245346096838814728167341e-02,
|
|
|
|
|
+ -6.369601011048822977293753932627342482077585617391852852955559e-03,
|
|
|
|
|
+ -1.820458915566242322836631665832145136570132777862391313328351e-03,
|
|
|
|
|
+ 7.902051009575939937150950543290226440287715441826917281929124e-04,
|
|
|
|
|
+ 3.296651737931830308416338897758022998655744276957481989605186e-04,
|
|
|
|
|
+ -5.019277455327664998007173088097694083956570594580641192332170e-05,
|
|
|
|
|
+ -2.446573425530813115445387662881902303945941576472342106918209e-05
|
|
|
|
|
+};
|
|
|
|
|
+
|
|
|
|
|
+
|
|
|
|
|
+static const float coif4[24] = {
|
|
|
|
|
+ 6.309612114309468490753696608619526520153127603444406835368201e-04,
|
|
|
|
|
+ -1.152225143769973488683007937016166047881572156705066038094891e-03,
|
|
|
|
|
+ -5.194525163470323267558201363327294331811309729430512113592118e-03,
|
|
|
|
|
+ 1.136246148326482276463392678363118465908960082105224676102131e-02,
|
|
|
|
|
+ 1.886723856956305960822813160712701905823879297781452350370094e-02,
|
|
|
|
|
+ -5.746424190192718517290527411385172124443396690932404284859269e-02,
|
|
|
|
|
+ -3.965265296244913762718094206756579981738035770770645437919302e-02,
|
|
|
|
|
+ 2.936674050161006858761278962798582650835466243678172528509866e-01,
|
|
|
|
|
+ 5.531264550395492870333469741987846570947502710783248169642137e-01,
|
|
|
|
|
+ 3.071573096678856987248881030393884808414165269795297009902001e-01,
|
|
|
|
|
+ -4.711273752389572084912399351781012121935994396763702238263689e-02,
|
|
|
|
|
+ -6.803811467802056988332974920928626798429778679560269769187728e-02,
|
|
|
|
|
+ 2.781363695846951303169163645831936314699164412528991864702607e-02,
|
|
|
|
|
+ 1.773583142270308388403079552822372238681544967313003044695583e-02,
|
|
|
|
|
+ -1.075631615508724933047071603601897536695959225169888787867102e-02,
|
|
|
|
|
+ -4.001010844950535391911552472397083276670126595827549403173754e-03,
|
|
|
|
|
+ 2.652664913530499860820143301690017184933302935238430721089152e-03,
|
|
|
|
|
+ 8.955939276952843603555618778866181384528643960440369133096025e-04,
|
|
|
|
|
+ -4.165001950941708741516836418852536615951250588002878691463468e-04,
|
|
|
|
|
+ -1.838296167136253805617482342622910940008368723403836355183423e-04,
|
|
|
|
|
+ 4.408022661597206973006038672236031501663774161685451815597956e-05,
|
|
|
|
|
+ 2.208284691230832960893331999804142845136324572860276715790883e-05,
|
|
|
|
|
+ -2.304919162676504406778986897925054839632903355820414483306851e-06,
|
|
|
|
|
+ -1.262179179994622253884862172782890488140153502131112374520603e-06
|
|
|
|
|
+};
|
|
|
|
|
+
|
|
|
|
|
+
|
|
|
|
|
+
|
|
|
|
|
+// static const float coif5[30] = {
|
|
|
|
|
+// -1.499645228345950331670593167919531667975440598691604525531231e-04,
|
|
|
|
|
+// 2.535527523580334712936363872191554706055603482812691726895588e-04,
|
|
|
|
|
+// 1.540286725995222360335148244676269541414659303531250711822333e-03,
|
|
|
|
|
+// -2.941078164035693185044038586065593320891475311414770624555173e-03,
|
|
|
|
|
+// -7.164112349410053294382279572472252500899544810929605832362178e-03,
|
|
|
|
|
+// 1.655218330649288840540841623080651353621667424921282557975513e-02,
|
|
|
|
|
+// 1.991901719798432056056857854066125809443504706772520641876273e-02,
|
|
|
|
|
+// -6.499837825472324963374262221660858232544804226063450042795603e-02,
|
|
|
|
|
+// -3.680255347446873527191823500872992242220223547780834450868002e-02,
|
|
|
|
|
+// 2.980959014587191795511466861338063554509597132272839414668911e-01,
|
|
|
|
|
+// 5.475082713540367154128337935687830970431964302909253422329131e-01,
|
|
|
|
|
+// 3.097002590784203529311533316221254677074498876376965941549923e-01,
|
|
|
|
|
+// -4.386731482823615640442730013366750193381707273908757638050452e-02,
|
|
|
|
|
+// -7.464442013283971243472663968192859562973186442054433655531762e-02,
|
|
|
|
|
+// 2.919469277528073666095772398605275751022315529465178441510318e-02,
|
|
|
|
|
+// 2.310457227706684192610065243663928370022983285246219996141160e-02,
|
|
|
|
|
+// -1.397129268638200558584119246355879336305763752871371182932059e-02,
|
|
|
|
|
+// -6.476749751505861835547590642967453082384538848552165075614441e-03,
|
|
|
|
|
+// 4.781116799130657606400088024549264921093190305150784065791191e-03,
|
|
|
|
|
+// 1.719383484385504023022397097446276782318002683055773803854075e-03,
|
|
|
|
|
+// -1.174947934413537690027670037110105795928147523549002426409332e-03,
|
|
|
|
|
+// -4.508222400696236312231932151038336110220594834213702970043431e-04,
|
|
|
|
|
+// 2.134457975086291667348984871136041914777578046177470626552867e-04,
|
|
|
|
|
+// 9.924691139873533169989496559631669037970741600337089424730635e-05,
|
|
|
|
|
+// -2.914684388622130824599478843558087403539428940986384077972155e-05,
|
|
|
|
|
+// -1.504031798197685905639227292876711236513927746903476131955063e-05,
|
|
|
|
|
+// 2.616809660013118152124234488302931243021794024318439103773996e-06,
|
|
|
|
|
+// 1.457502921355163070577152619048168436286350537937563166257584e-06,
|
|
|
|
|
+// -1.148199649902979726237655584441763456854312591680755421569962e-07,
|
|
|
|
|
+// -6.791060677322355511541065559242475254516249773485524025251102e-08
|
|
|
|
|
+// };
|
|
|
|
|
+
|
|
|
|
|
+static const float sym2[4] = {
|
|
|
|
|
+ 0.48296291314469025, 0.83651630373746899,
|
|
|
|
|
+ 0.22414386804185735, -0.12940952255092145 };
|
|
|
|
|
+
|
|
|
|
|
+
|
|
|
|
|
+static const float sym3[6] = {
|
|
|
|
|
+ 0.33267055295095688, 0.80689150931333875,
|
|
|
|
|
+ 0.45987750211933132, -0.13501102001039084,
|
|
|
|
|
+ -0.085441273882241486, 0.035226291882100656 };
|
|
|
|
|
+
|
|
|
|
|
+static const float sym4[8] = {
|
|
|
|
|
+ 0.032223100604042702, -0.012603967262037833,
|
|
|
|
|
+ -0.099219543576847216, 0.29785779560527736,
|
|
|
|
|
+ 0.80373875180591614, 0.49761866763201545,
|
|
|
|
|
+ -0.02963552764599851, -0.075765714789273325 };
|
|
|
|
|
+
|
|
|
|
|
+// static const float sym5[10] = {
|
|
|
|
|
+// 0.019538882735286728, -0.021101834024758855,
|
|
|
|
|
+// -0.17532808990845047, 0.016602105764522319,
|
|
|
|
|
+// 0.63397896345821192, 0.72340769040242059,
|
|
|
|
|
+// 0.1993975339773936, -0.039134249302383094,
|
|
|
|
|
+// 0.029519490925774643, 0.027333068345077982 };
|
|
|
|
|
+
|
|
|
|
|
+/* 要求:
|
|
|
|
|
+* (1)C0^2 + C1^2 + C2^2 + C3^2 = 1
|
|
|
|
|
+* (2)C1 * C3 + C0 * C2 = 0
|
|
|
|
|
+*/
|
|
|
|
|
+
|
|
|
|
|
+/* 矩阵变换需要用到的参数 */
|
|
|
|
|
+static Filter DBn = {
|
|
|
|
|
+ .Len = 4,
|
|
|
|
|
+ .C0 = (const float) 0.482962913144690,
|
|
|
|
|
+ .C1 = (const float) 0.836516303737469,
|
|
|
|
|
+ .C2 = (const float) 0.224143868041857,
|
|
|
|
|
+ .C3 = (const float) -0.129409522550921,
|
|
|
|
|
+};
|
|
|
|
|
+
|
|
|
|
|
+
|
|
|
|
|
+/* 按分解矩阵进行分解 */
|
|
|
|
|
+static void transform(float *input_Data, float *output_Data, int length)
|
|
|
|
|
+{
|
|
|
|
|
+ int i=0;
|
|
|
|
|
+ int count=0;
|
|
|
|
|
+
|
|
|
|
|
+ for(i=0; i<=(length-3); i=i+1)
|
|
|
|
|
+ {
|
|
|
|
|
+ output_Data[count++] = DBn.C0 * input_Data[i] + DBn.C1 * input_Data[i+1] + DBn.C2 * input_Data[i+2] + DBn.C3 * input_Data[i+3];
|
|
|
|
|
+ }
|
|
|
|
|
+ output_Data[count++] = DBn.C2 * input_Data[0] + DBn.C3 * input_Data[1] + DBn.C0 * input_Data[length-2] + DBn.C1 * input_Data[length-1];
|
|
|
|
|
+ output_Data[count] = DBn.C1 * input_Data[0] - DBn.C0 * input_Data[1] + DBn.C3 * input_Data[length-2] - DBn.C2 * input_Data[length-1];
|
|
|
|
|
+}
|
|
|
|
|
+
|
|
|
|
|
+
|
|
|
|
|
+/******************************************************************************
|
|
|
|
|
+* 功能:利用矩阵进行分解
|
|
|
|
|
+* @srcData: 源信号;
|
|
|
|
|
+* @resData: 保存结果;
|
|
|
|
|
+* @srcLen: 源信号长度;
|
|
|
|
|
+*******************************************************************************/
|
|
|
|
|
+void DWT_Matrix_Transform(float *srcData, s16 *resData, int srcLen)
|
|
|
|
|
+{
|
|
|
|
|
+
|
|
|
|
|
+ int i=0;
|
|
|
|
|
+ float out1Data[srcLen];
|
|
|
|
|
+
|
|
|
|
|
+ memset(out1Data, 0, srcLen);
|
|
|
|
|
+ transform(srcData, out1Data, srcLen);
|
|
|
|
|
+
|
|
|
|
|
+ for (i = 0; i < srcLen; i +=1)
|
|
|
|
|
+ {
|
|
|
|
|
+ resData[i] = (s16)rt_round(out1Data[i]); //此处有风险,小心处理
|
|
|
|
|
+ // *((s16 *)&resData[i]) = (s16)(out1Data[i]); //此处有风险,小心处理
|
|
|
|
|
+ }
|
|
|
|
|
+
|
|
|
|
|
+ return ;
|
|
|
|
|
+}
|
|
|
|
|
+
|
|
|
|
|
+/******************************************************************************
|
|
|
|
|
+函数名称: hilbert_para_init
|
|
|
|
|
+函数版本: 01.01
|
|
|
|
|
+创建作者: ygl
|
|
|
|
|
+创建日期: 2021.11.17
|
|
|
|
|
+函数说明: hilbert初始化。
|
|
|
|
|
+参数说明: 无
|
|
|
|
|
+返回值: 长度
|
|
|
|
|
+修改记录:
|
|
|
|
|
+*******************************************************************************/
|
|
|
|
|
+void _hilbert_para_init(void) //初始化系数
|
|
|
|
|
+{
|
|
|
|
|
+ memset(&coef1, 0, sizeof(coef1));
|
|
|
|
|
+ memset(&coef2, 0, sizeof(coef2));
|
|
|
|
|
+ memset(&coef3, 0, sizeof(coef3));
|
|
|
|
|
+ memset(&coef4, 0, sizeof(coef4));
|
|
|
|
|
+}
|
|
|
|
|
+/******************************************************************************
|
|
|
|
|
+函数名称: filter_Init
|
|
|
|
|
+函数版本: 01.01
|
|
|
|
|
+创建作者: ygl
|
|
|
|
|
+创建日期: 2021.11.17
|
|
|
|
|
+函数说明: 卷积核初始化。
|
|
|
|
|
+参数说明: 无
|
|
|
|
|
+返回值: 长度
|
|
|
|
|
+修改记录:
|
|
|
|
|
+*******************************************************************************/
|
|
|
|
|
+void filter_Init(void)
|
|
|
|
|
+{
|
|
|
|
|
+ #if 0
|
|
|
|
|
+ g_fpuo = rt_file_open("/app/uo.txt", O_CREAT|O_RDWR,0); //近似系数存储文件
|
|
|
|
|
+ if(IS_ERR(g_fpuo))
|
|
|
|
|
+ {
|
|
|
|
|
+ rt_printf("LINE %d write file error!\n",__LINE__);
|
|
|
|
|
+ rt_printf("write file error!\n");
|
|
|
|
|
+ // return -1;
|
|
|
|
|
+ }
|
|
|
|
|
+ g_fpio = rt_file_open("/app/io.txt", O_CREAT|O_RDWR,0); //近似系数存储文件
|
|
|
|
|
+ if(IS_ERR(g_fpio))
|
|
|
|
|
+ {
|
|
|
|
|
+ rt_printf("LINE %d write file error!\n",__LINE__);
|
|
|
|
|
+ rt_printf("write file error!\n");
|
|
|
|
|
+ // return -1;
|
|
|
|
|
+ }
|
|
|
|
|
+ #endif
|
|
|
|
|
+ _hilbert_para_init();
|
|
|
|
|
+}
|
|
|
|
|
+
|
|
|
|
|
+float hilbert_filter(float *channel1, float *channel2, int len)
|
|
|
|
|
+{
|
|
|
|
|
+ int i=0;
|
|
|
|
|
+ const float c1=2.0216;
|
|
|
|
|
+ const float c2=1.2562;
|
|
|
|
|
+ const float c3=0.2276;
|
|
|
|
|
+ const float c4=1.5219;
|
|
|
|
|
+ const float c5=0.6205;
|
|
|
|
|
+ const float c6=0.0463;
|
|
|
|
|
+
|
|
|
|
|
+ // const float c1=2.1441;
|
|
|
|
|
+ // const float c2=1.4380;
|
|
|
|
|
+ // const float c3=0.2904;
|
|
|
|
|
+ // const float c4=1.6444;
|
|
|
|
|
+ // const float c5=0.7421;
|
|
|
|
|
+ // const float c6=0.0658;
|
|
|
|
|
+
|
|
|
|
|
+ float channel_f1 = 0;
|
|
|
|
|
+ float channel_f2 = 0;
|
|
|
|
|
+ float channel_f3 = 0;
|
|
|
|
|
+ float channel_f4 = 0;
|
|
|
|
|
+ float Qval = 0;
|
|
|
|
|
+
|
|
|
|
|
+ memset(&coef1,0,sizeof(coef1));
|
|
|
|
|
+ memset(&coef2,0,sizeof(coef2));
|
|
|
|
|
+ memset(&coef3,0,sizeof(coef3));
|
|
|
|
|
+ memset(&coef4,0,sizeof(coef4));
|
|
|
|
|
+
|
|
|
|
|
+ for(i=0;i<len;i++)
|
|
|
|
|
+ {
|
|
|
|
|
+ channel_f1 = c1*(coef1.y2)-c2*(coef1.y4)+c3*(coef1.y6)+c3*(coef1.x1)-c2*(coef1.x3)+c1*(coef1.x5)-coef1.x7;
|
|
|
|
|
+
|
|
|
|
|
+ coef1.y6 = coef1.y5;
|
|
|
|
|
+ coef1.y5 = coef1.y4;
|
|
|
|
|
+ coef1.y4 = coef1.y3;
|
|
|
|
|
+ coef1.y3 = coef1.y2;
|
|
|
|
|
+ coef1.y2 = coef1.y1;
|
|
|
|
|
+ coef1.y1 = channel_f1;
|
|
|
|
|
+ coef1.x7 = coef1.x6;
|
|
|
|
|
+ coef1.x6 = coef1.x5;
|
|
|
|
|
+ coef1.x5 = coef1.x4;
|
|
|
|
|
+ coef1.x4 = coef1.x3;
|
|
|
|
|
+ coef1.x3 = coef1.x2;
|
|
|
|
|
+ coef1.x2 = coef1.x1;
|
|
|
|
|
+ coef1.x1 = channel1[i];
|
|
|
|
|
+
|
|
|
|
|
+
|
|
|
|
|
+ channel_f2 = c4*(coef2.y2)-c5*(coef2.y4)+c6*(coef2.y6)+c6*(channel1[i])-c5*(coef2.x2)+c4*(coef2.x4)-(coef2.x6);
|
|
|
|
|
+ coef2.y6 = coef2.y5;
|
|
|
|
|
+ coef2.y5 = coef2.y4;
|
|
|
|
|
+ coef2.y4 = coef2.y3;
|
|
|
|
|
+ coef2.y3 = coef2.y2;
|
|
|
|
|
+ coef2.y2 = coef2.y1;
|
|
|
|
|
+ coef2.y1 = channel_f2;
|
|
|
|
|
+ coef2.x6 = coef2.x5;
|
|
|
|
|
+ coef2.x5 = coef2.x4;
|
|
|
|
|
+ coef2.x4 = coef2.x3;
|
|
|
|
|
+ coef2.x3 = coef2.x2;
|
|
|
|
|
+ coef2.x2 = coef2.x1;
|
|
|
|
|
+ coef2.x1 = channel1[i];
|
|
|
|
|
+
|
|
|
|
|
+ channel_f3 = c1*(coef3.y2)-c2*(coef3.y4)+c3*(coef3.y6)+c3*(coef3.x1)-c2*(coef3.x3)+c1*(coef3.x5)-coef3.x7;
|
|
|
|
|
+ coef3.y6 = coef3.y5;
|
|
|
|
|
+ coef3.y5 = coef3.y4;
|
|
|
|
|
+ coef3.y4 = coef3.y3;
|
|
|
|
|
+ coef3.y3 = coef3.y2;
|
|
|
|
|
+ coef3.y2 = coef3.y1;
|
|
|
|
|
+ coef3.y1 = channel_f3;
|
|
|
|
|
+ coef3.x7 = coef3.x6;
|
|
|
|
|
+ coef3.x6 = coef3.x5;
|
|
|
|
|
+ coef3.x5 = coef3.x4;
|
|
|
|
|
+ coef3.x4 = coef3.x3;
|
|
|
|
|
+ coef3.x3 = coef3.x2;
|
|
|
|
|
+ coef3.x2 = coef3.x1;
|
|
|
|
|
+ coef3.x1 = channel2[i];
|
|
|
|
|
+
|
|
|
|
|
+
|
|
|
|
|
+ channel_f4 = c4*(coef4.y2)-c5*(coef4.y4)+c6*(coef4.y6)+c6*(channel2[i])-c5*(coef4.x2)+c4*(coef4.x4)-(coef4.x6);
|
|
|
|
|
+ coef4.y6 = coef4.y5;
|
|
|
|
|
+ coef4.y5 = coef4.y4;
|
|
|
|
|
+ coef4.y4 = coef4.y3;
|
|
|
|
|
+ coef4.y3 = coef4.y2;
|
|
|
|
|
+ coef4.y2 = coef4.y1;
|
|
|
|
|
+ coef4.y1 = channel_f4;
|
|
|
|
|
+ coef4.x6 = coef4.x5;
|
|
|
|
|
+ coef4.x5 = coef4.x4;
|
|
|
|
|
+ coef4.x4 = coef4.x3;
|
|
|
|
|
+ coef4.x3 = coef4.x2;
|
|
|
|
|
+ coef4.x2 = coef4.x1;
|
|
|
|
|
+ coef4.x1 = channel2[i];
|
|
|
|
|
+
|
|
|
|
|
+ Qval += (channel_f1*channel_f4 - channel_f2*channel_f3)/2;
|
|
|
|
|
+
|
|
|
|
|
+ }
|
|
|
|
|
+ Qval =(float) Qval/len;
|
|
|
|
|
+ // rt_printf("len = %d,qval = %f\r\n",len,Qval);
|
|
|
|
|
+ return Qval;
|
|
|
|
|
+}
|
|
|
|
|
+
|
|
|
|
|
+
|
|
|
|
|
+int filtlength(const char* name) {
|
|
|
|
|
+ int len = strlen(name);
|
|
|
|
|
+ int i = 0;
|
|
|
|
|
+ char *new_str = NULL;
|
|
|
|
|
+ int N = 0;
|
|
|
|
|
+ if (!strcmp(name,"haar") || !strcmp(name,"db1")) {
|
|
|
|
|
+ return 2;
|
|
|
|
|
+ }
|
|
|
|
|
+ else if (len > 2 && strstr(name, "db") != NULL)
|
|
|
|
|
+ {
|
|
|
|
|
+ new_str = (char*)rt_malloc(sizeof(char)*(len-2 + 1));
|
|
|
|
|
+ for (i = 2; i < len + 1; i++)
|
|
|
|
|
+ new_str[i - 2] = name[i];
|
|
|
|
|
+
|
|
|
|
|
+ N = atoi(new_str);
|
|
|
|
|
+ rt_free(new_str);
|
|
|
|
|
+ if (N>38)
|
|
|
|
|
+ {
|
|
|
|
|
+ rt_printf("\n Filter Not in Database \n");
|
|
|
|
|
+ return -1;
|
|
|
|
|
+ }
|
|
|
|
|
+
|
|
|
|
|
+ return N * 2;
|
|
|
|
|
+ }
|
|
|
|
|
+
|
|
|
|
|
+ else if (len > 4 && strstr(name, "coif") != NULL)
|
|
|
|
|
+ {
|
|
|
|
|
+ new_str = (char*)rt_malloc(sizeof(char)*(len - 4 + 1));
|
|
|
|
|
+ for (i = 4; i < len + 1; i++)
|
|
|
|
|
+ new_str[i - 4] = name[i];
|
|
|
|
|
+
|
|
|
|
|
+ N = atoi(new_str);
|
|
|
|
|
+ rt_free(new_str);
|
|
|
|
|
+ if (N>17)
|
|
|
|
|
+ {
|
|
|
|
|
+ rt_printf("\n Filter Not in Database \n");
|
|
|
|
|
+ return -1;
|
|
|
|
|
+ }
|
|
|
|
|
+
|
|
|
|
|
+ return N * 6;
|
|
|
|
|
+ }
|
|
|
|
|
+ else if (len > 3 && strstr(name, "sym") != NULL)
|
|
|
|
|
+ {
|
|
|
|
|
+ new_str = (char*)rt_malloc(sizeof(char)*(len - 3 + 1));
|
|
|
|
|
+ for (i = 3; i < len + 1; i++)
|
|
|
|
|
+ new_str[i - 3] = name[i];
|
|
|
|
|
+
|
|
|
|
|
+ N = atoi(new_str);
|
|
|
|
|
+ rt_free(new_str);
|
|
|
|
|
+ if (N>20 || N < 2)
|
|
|
|
|
+ {
|
|
|
|
|
+ rt_printf("\n Filter Not in Database \n");
|
|
|
|
|
+ return -1;
|
|
|
|
|
+ }
|
|
|
|
|
+
|
|
|
|
|
+ return N * 2;
|
|
|
|
|
+ }
|
|
|
|
|
+ else {
|
|
|
|
|
+ rt_printf("\n Filter Not in Database \n");
|
|
|
|
|
+ return -1;
|
|
|
|
|
+ }
|
|
|
|
|
+
|
|
|
|
|
+}
|
|
|
|
|
+
|
|
|
|
|
+void copy_reverse(const float *in, int N,float *out)
|
|
|
|
|
+{
|
|
|
|
|
+ int count = 0;
|
|
|
|
|
+ for (count = 0; count < N; count++)
|
|
|
|
|
+ out[count] = in[N - count - 1];
|
|
|
|
|
+}
|
|
|
|
|
+
|
|
|
|
|
+void qmf_even(const float *in, int N,float *out)
|
|
|
|
|
+{
|
|
|
|
|
+ int count = 0;
|
|
|
|
|
+ for (count = 0; count < N; count++)
|
|
|
|
|
+ {
|
|
|
|
|
+ out[count] = in[N - count - 1];
|
|
|
|
|
+ if (count % 2 != 0)
|
|
|
|
|
+ {
|
|
|
|
|
+ out[count] = -1 * out[count];
|
|
|
|
|
+ }
|
|
|
|
|
+ }
|
|
|
|
|
+}
|
|
|
|
|
+
|
|
|
|
|
+void qmf_wrev(const float *in, int N, float *out)
|
|
|
|
|
+{
|
|
|
|
|
+ float *sigOutTemp;
|
|
|
|
|
+ sigOutTemp = (float*)rt_malloc(N*sizeof(float));
|
|
|
|
|
+
|
|
|
|
|
+ qmf_even(in, N, sigOutTemp);
|
|
|
|
|
+ copy_reverse(sigOutTemp, N, out);
|
|
|
|
|
+
|
|
|
|
|
+ rt_free(sigOutTemp);
|
|
|
|
|
+ return;
|
|
|
|
|
+}
|
|
|
|
|
+
|
|
|
|
|
+void wave_copy(const float *in, int N, float *out)
|
|
|
|
|
+{
|
|
|
|
|
+ int count = 0;
|
|
|
|
|
+ for (count = 0; count < N; count++)
|
|
|
|
|
+ out[count] = in[count];
|
|
|
|
|
+}
|
|
|
|
|
+
|
|
|
|
|
+int filtcoef(const char* name, float *lp1, float *hp1, float *lp2, float *hp2) {
|
|
|
|
|
+ int i = 0;
|
|
|
|
|
+ int N = filtlength(name);
|
|
|
|
|
+ if (!strcmp(name,"haar") || !strcmp(name,"db1")) {
|
|
|
|
|
+ copy_reverse(db1, N, lp1);
|
|
|
|
|
+ qmf_wrev(db1, N, hp1);
|
|
|
|
|
+ wave_copy(db1, N, lp2);
|
|
|
|
|
+ qmf_even(db1, N, hp2);
|
|
|
|
|
+
|
|
|
|
|
+ return N;
|
|
|
|
|
+ }
|
|
|
|
|
+ else if (!strcmp(name,"db2")){
|
|
|
|
|
+ copy_reverse(db2, N, lp1);
|
|
|
|
|
+ qmf_wrev(db2, N, hp1);
|
|
|
|
|
+ wave_copy(db2, N, lp2);
|
|
|
|
|
+ qmf_even(db2, N, hp2);
|
|
|
|
|
+
|
|
|
|
|
+ return N;
|
|
|
|
|
+ }
|
|
|
|
|
+ else if (!strcmp(name,"db3")){
|
|
|
|
|
+ copy_reverse(db3, N, lp1);
|
|
|
|
|
+ qmf_wrev(db3, N, hp1);
|
|
|
|
|
+ wave_copy(db3, N, lp2);
|
|
|
|
|
+ qmf_even(db3, N, hp2);
|
|
|
|
|
+
|
|
|
|
|
+ return N;
|
|
|
|
|
+ }
|
|
|
|
|
+ else if (!strcmp(name,"db4")){
|
|
|
|
|
+ copy_reverse(db4, N, lp1);
|
|
|
|
|
+ qmf_wrev(db4, N, hp1);
|
|
|
|
|
+ wave_copy(db4, N, lp2);
|
|
|
|
|
+ qmf_even(db4, N, hp2);
|
|
|
|
|
+
|
|
|
|
|
+ return N;
|
|
|
|
|
+ }
|
|
|
|
|
+
|
|
|
|
|
+ else if (!strcmp(name,"coif1")){
|
|
|
|
|
+ float *coeffTemp;
|
|
|
|
|
+ coeffTemp = (float*)rt_malloc(N*sizeof(float));
|
|
|
|
|
+
|
|
|
|
|
+ wave_copy(coif1, N, coeffTemp);
|
|
|
|
|
+ for (i = 0; i < N; ++i) {
|
|
|
|
|
+ coeffTemp[i] *= M_SQRT2;
|
|
|
|
|
+ }
|
|
|
|
|
+
|
|
|
|
|
+ copy_reverse(coeffTemp, N, lp1);
|
|
|
|
|
+ qmf_wrev(coeffTemp, N, hp1);
|
|
|
|
|
+ wave_copy(coeffTemp, N, lp2);
|
|
|
|
|
+ qmf_even(coeffTemp, N, hp2);
|
|
|
|
|
+ rt_free(coeffTemp);
|
|
|
|
|
+
|
|
|
|
|
+ return N;
|
|
|
|
|
+ }
|
|
|
|
|
+ else if (!strcmp(name,"coif2")){
|
|
|
|
|
+ float *coeffTemp;
|
|
|
|
|
+ coeffTemp = (float*)rt_malloc(N*sizeof(float));
|
|
|
|
|
+
|
|
|
|
|
+ wave_copy(coif2, N, coeffTemp);
|
|
|
|
|
+ for (i = 0; i < N; ++i) {
|
|
|
|
|
+ coeffTemp[i] *= M_SQRT2;
|
|
|
|
|
+ }
|
|
|
|
|
+
|
|
|
|
|
+ copy_reverse(coeffTemp, N, lp1);
|
|
|
|
|
+ qmf_wrev(coeffTemp, N, hp1);
|
|
|
|
|
+ wave_copy(coeffTemp, N, lp2);
|
|
|
|
|
+ qmf_even(coeffTemp, N, hp2);
|
|
|
|
|
+ rt_free(coeffTemp);
|
|
|
|
|
+
|
|
|
|
|
+ return N;
|
|
|
|
|
+ }
|
|
|
|
|
+ else if (!strcmp(name,"coif3")){
|
|
|
|
|
+ float *coeffTemp;
|
|
|
|
|
+ coeffTemp = (float*)rt_malloc(N*sizeof(float));
|
|
|
|
|
+
|
|
|
|
|
+ wave_copy(coif3, N, coeffTemp);
|
|
|
|
|
+ for (i = 0; i < N; ++i) {
|
|
|
|
|
+ coeffTemp[i] *= M_SQRT2;
|
|
|
|
|
+ }
|
|
|
|
|
+
|
|
|
|
|
+ copy_reverse(coeffTemp, N, lp1);
|
|
|
|
|
+ qmf_wrev(coeffTemp, N, hp1);
|
|
|
|
|
+ wave_copy(coeffTemp, N, lp2);
|
|
|
|
|
+ qmf_even(coeffTemp, N, hp2);
|
|
|
|
|
+ rt_free(coeffTemp);
|
|
|
|
|
+
|
|
|
|
|
+ return N;
|
|
|
|
|
+ }
|
|
|
|
|
+ else if (!strcmp(name,"coif4")){
|
|
|
|
|
+ float *coeffTemp;
|
|
|
|
|
+ coeffTemp = (float*)rt_malloc(N*sizeof(float));
|
|
|
|
|
+
|
|
|
|
|
+ wave_copy(coif4, N, coeffTemp);
|
|
|
|
|
+ for (i = 0; i < N; ++i) {
|
|
|
|
|
+ coeffTemp[i] *= M_SQRT2;
|
|
|
|
|
+ }
|
|
|
|
|
+
|
|
|
|
|
+ copy_reverse(coeffTemp, N, lp1);
|
|
|
|
|
+ qmf_wrev(coeffTemp, N, hp1);
|
|
|
|
|
+ wave_copy(coeffTemp, N, lp2);
|
|
|
|
|
+ qmf_even(coeffTemp, N, hp2);
|
|
|
|
|
+ rt_free(coeffTemp);
|
|
|
|
|
+
|
|
|
|
|
+ return N;
|
|
|
|
|
+ }
|
|
|
|
|
+
|
|
|
|
|
+ else if (!strcmp(name,"sym2") || !strcmp(name,"sym1")){
|
|
|
|
|
+ copy_reverse(sym2, N, lp1);
|
|
|
|
|
+ qmf_wrev(sym2, N, hp1);
|
|
|
|
|
+ wave_copy(sym2, N, lp2);
|
|
|
|
|
+ qmf_even(sym2, N, hp2);
|
|
|
|
|
+ return N;
|
|
|
|
|
+ }
|
|
|
|
|
+
|
|
|
|
|
+ else if (!strcmp(name,"sym3")){
|
|
|
|
|
+ copy_reverse(sym3, N, lp1);
|
|
|
|
|
+ qmf_wrev(sym3, N, hp1);
|
|
|
|
|
+ wave_copy(sym3, N, lp2);
|
|
|
|
|
+ qmf_even(sym3, N, hp2);
|
|
|
|
|
+ return N;
|
|
|
|
|
+ }
|
|
|
|
|
+
|
|
|
|
|
+ else if (!strcmp(name,"sym4")){
|
|
|
|
|
+ copy_reverse(sym4, N, lp1);
|
|
|
|
|
+ qmf_wrev(sym4, N, hp1);
|
|
|
|
|
+ wave_copy(sym4, N, lp2);
|
|
|
|
|
+ qmf_even(sym4, N, hp2);
|
|
|
|
|
+ return N;
|
|
|
|
|
+ }
|
|
|
|
|
+
|
|
|
|
|
+ else {
|
|
|
|
|
+ rt_printf("\n Filter Not in Database \n");
|
|
|
|
|
+ return -1;
|
|
|
|
|
+ }
|
|
|
|
|
+
|
|
|
|
|
+ return 0;
|
|
|
|
|
+}
|
|
|
|
|
+
|
|
|
|
|
+
|
|
|
|
|
+wave_object wave_init(const char* wname) {
|
|
|
|
|
+ wave_object obj = NULL;
|
|
|
|
|
+ int retval;
|
|
|
|
|
+ retval = 0;
|
|
|
|
|
+
|
|
|
|
|
+ if (wname != NULL) {
|
|
|
|
|
+ retval = filtlength(wname);
|
|
|
|
|
+ }
|
|
|
|
|
+
|
|
|
|
|
+ obj = (wave_object)rt_malloc(sizeof(struct wave_set) + sizeof(float)* 4 * retval);
|
|
|
|
|
+
|
|
|
|
|
+ obj->filtlength = retval;
|
|
|
|
|
+ obj->lpd_len = obj->hpd_len = obj->lpr_len = obj->hpr_len = obj->filtlength;
|
|
|
|
|
+ strcpy(obj->wname, wname);
|
|
|
|
|
+ if (wname != NULL) {
|
|
|
|
|
+ filtcoef(wname,obj->params,obj->params+retval,obj->params+2*retval,obj->params+3*retval);
|
|
|
|
|
+ }
|
|
|
|
|
+ obj->lpd = &obj->params[0];
|
|
|
|
|
+ obj->hpd = &obj->params[retval];
|
|
|
|
|
+ obj->lpr = &obj->params[2 * retval];
|
|
|
|
|
+ obj->hpr = &obj->params[3 * retval];
|
|
|
|
|
+ return obj;
|
|
|
|
|
+}
|
|
|
|
|
+
|
|
|
|
|
+wt_object wt_init(wave_object wave,const char* method, int siglength,int J) {
|
|
|
|
|
+ int i;
|
|
|
|
|
+ wt_object obj = NULL;
|
|
|
|
|
+
|
|
|
|
|
+ if (J > 100) {
|
|
|
|
|
+ rt_printf("\n The Decomposition Iterations Cannot Exceed 100. Exiting \n");
|
|
|
|
|
+ return NULL;
|
|
|
|
|
+ }
|
|
|
|
|
+
|
|
|
|
|
+ if (!strcmp(method, "modwt") || !strcmp(method, "MODWT")) {
|
|
|
|
|
+
|
|
|
|
|
+ if (!strstr(wave->wname,"haar")) {
|
|
|
|
|
+ if (!strstr(wave->wname,"db")) {
|
|
|
|
|
+ if (!strstr(wave->wname, "sym")) {
|
|
|
|
|
+ if (!strstr(wave->wname, "coif")) {
|
|
|
|
|
+ rt_printf("\n MODWT is only implemented for orthogonal wavelet families - db, sym and coif \n");
|
|
|
|
|
+ return NULL;
|
|
|
|
|
+ }
|
|
|
|
|
+ }
|
|
|
|
|
+ }
|
|
|
|
|
+ }
|
|
|
|
|
+
|
|
|
|
|
+ // obj = (wt_object)rt_malloc(sizeof(struct wt_set) + sizeof(float)* (siglength * 2 * (J + 1)));
|
|
|
|
|
+ // obj->outlength = siglength * (J + 1); // Default
|
|
|
|
|
+ obj = (wt_object)rt_malloc(sizeof(struct wt_set) + sizeof(float)* (siglength * 2));
|
|
|
|
|
+ obj->outlength = siglength ; //
|
|
|
|
|
+ strcpy(obj->ext, "per"); // Default
|
|
|
|
|
+ }
|
|
|
|
|
+
|
|
|
|
|
+ obj->wave = wave;
|
|
|
|
|
+ obj->siglength = siglength;
|
|
|
|
|
+ obj->modwtsiglength = siglength;
|
|
|
|
|
+ obj->J = J;
|
|
|
|
|
+ // obj->MaxIter = MaxIter;
|
|
|
|
|
+ strcpy(obj->method, method);
|
|
|
|
|
+
|
|
|
|
|
+ if (siglength % 2 == 0) {
|
|
|
|
|
+ obj->even = 1;
|
|
|
|
|
+ }
|
|
|
|
|
+ else {
|
|
|
|
|
+ obj->even = 0;
|
|
|
|
|
+ }
|
|
|
|
|
+
|
|
|
|
|
+ obj->cobj = NULL;
|
|
|
|
|
+
|
|
|
|
|
+ strcpy(obj->cmethod, "direct"); // Default
|
|
|
|
|
+ obj->cfftset = 0;
|
|
|
|
|
+ obj->lenlength = J + 2;
|
|
|
|
|
+ obj->output = &obj->params[0];
|
|
|
|
|
+
|
|
|
|
|
+ if (!strcmp(method, "modwt") || !strcmp(method, "MODWT")) {
|
|
|
|
|
+ for (i = 0; i < siglength * 2; ++i) {
|
|
|
|
|
+ obj->params[i] = 0.0;
|
|
|
|
|
+ }
|
|
|
|
|
+ }
|
|
|
|
|
+
|
|
|
|
|
+ return obj;
|
|
|
|
|
+}
|
|
|
|
|
+
|
|
|
|
|
+static void modwt_per(wt_object wt, int M, float *inp, float *cA, int len_cA) {
|
|
|
|
|
+ int l, i, t, len_avg;
|
|
|
|
|
+ float s;
|
|
|
|
|
+ float *filt;
|
|
|
|
|
+ len_avg = wt->wave->lpd_len;
|
|
|
|
|
+
|
|
|
|
|
+ filt = (float*)rt_malloc(sizeof(float)* len_avg);
|
|
|
|
|
+ s = M_SQRT2;//sqrt(2.0);
|
|
|
|
|
+ for (i = 0; i < len_avg; ++i) {
|
|
|
|
|
+ filt[i] = wt->wave->lpd[i] / s;
|
|
|
|
|
+ }
|
|
|
|
|
+
|
|
|
|
|
+ for (i = 0; i < len_cA; ++i) {
|
|
|
|
|
+ t = i;
|
|
|
|
|
+ cA[i] = filt[0] * inp[t];
|
|
|
|
|
+
|
|
|
|
|
+ for (l = 1; l < len_avg; l++) {
|
|
|
|
|
+ t -= M;
|
|
|
|
|
+ while (t >= len_cA) {
|
|
|
|
|
+ t -= len_cA;
|
|
|
|
|
+ }
|
|
|
|
|
+ while (t < 0) {
|
|
|
|
|
+ // t += len_cA; //ygl 传入波形数据,有可能首尾不对称。此处首核卷积改为第一点开始
|
|
|
|
|
+ t += M;
|
|
|
|
|
+ }
|
|
|
|
|
+
|
|
|
|
|
+ cA[i] += filt[l] * inp[t];
|
|
|
|
|
+ }
|
|
|
|
|
+ }
|
|
|
|
|
+ rt_free(filt);
|
|
|
|
|
+}
|
|
|
|
|
+
|
|
|
|
|
+static void modwt_direct(wt_object wt, const float *inp) {
|
|
|
|
|
+ int i, J, temp_len, iter, M;
|
|
|
|
|
+ float *cA;
|
|
|
|
|
+
|
|
|
|
|
+
|
|
|
|
|
+ if (strcmp(wt->ext, "per")) {
|
|
|
|
|
+ rt_printf("MODWT direct method only uses periodic extension per. \n");
|
|
|
|
|
+ rt_printf(" Use MODWT fft method for symmetric extension sym \n");
|
|
|
|
|
+ return;
|
|
|
|
|
+ }
|
|
|
|
|
+
|
|
|
|
|
+ temp_len = wt->siglength;
|
|
|
|
|
+ J = wt->J;
|
|
|
|
|
+ wt->length[0] = wt->length[J] = temp_len;
|
|
|
|
|
+ wt->outlength = wt->length[J + 1] = temp_len; //只考虑近似值
|
|
|
|
|
+ M = 1;
|
|
|
|
|
+ for (iter = 1; iter < J; ++iter) {
|
|
|
|
|
+ M = 2 * M;
|
|
|
|
|
+ wt->length[iter] = temp_len;
|
|
|
|
|
+ }
|
|
|
|
|
+
|
|
|
|
|
+
|
|
|
|
|
+ cA = (float*)rt_malloc(sizeof(float)* temp_len);
|
|
|
|
|
+ M = 1;
|
|
|
|
|
+
|
|
|
|
|
+ for (i = 0; i < temp_len; ++i) {
|
|
|
|
|
+ wt->params[i] = inp[i];
|
|
|
|
|
+ }
|
|
|
|
|
+
|
|
|
|
|
+ for (iter = 0; iter < J; ++iter) {
|
|
|
|
|
+ if (iter > 0) {
|
|
|
|
|
+ M = 2 * M;
|
|
|
|
|
+ }
|
|
|
|
|
+
|
|
|
|
|
+ modwt_per(wt, M, wt->params, cA, temp_len);
|
|
|
|
|
+
|
|
|
|
|
+
|
|
|
|
|
+ for (i = 0; i < temp_len; ++i) {
|
|
|
|
|
+ wt->params[i] = cA[i];
|
|
|
|
|
+ }
|
|
|
|
|
+
|
|
|
|
|
+ }
|
|
|
|
|
+
|
|
|
|
|
+ rt_free(cA);
|
|
|
|
|
+
|
|
|
|
|
+}
|
|
|
|
|
+
|
|
|
|
|
+void modwt(wt_object wt, const float *inp) {
|
|
|
|
|
+ if (!strcmp(wt->cmethod, "direct")) {
|
|
|
|
|
+ modwt_direct(wt, inp);
|
|
|
|
|
+ }
|
|
|
|
|
+ else {
|
|
|
|
|
+ rt_printf("Error- Available Choices for this method are - direct and fft \n");
|
|
|
|
|
+ return;
|
|
|
|
|
+ }
|
|
|
|
|
+}
|
|
|
|
|
+
|
|
|
|
|
+void wave_free(wave_object object) {
|
|
|
|
|
+ rt_free(object);
|
|
|
|
|
+}
|
|
|
|
|
+
|
|
|
|
|
+void wt_free(wt_object object) {
|
|
|
|
|
+ rt_free(object);
|
|
|
|
|
+}
|
|
|
|
|
+
|
|
|
|
|
+
|
|
|
|
|
+/*------------------------------ 内部函数 -------------------------------------
|
|
|
|
|
+内部函数以下划线‘_’开头,不需要检查参数的合法性.
|
|
|
|
|
+*/
|
|
|
|
|
+
|
|
|
|
|
+
|
|
|
|
|
+/*------------------------------ 测试函数 -------------------------------------
|
|
|
|
|
+一个实体文件必须带一个本模块的测试函数来进行单元测试,如果的确不方便在本模块中
|
|
|
|
|
+进行单元测试,必须在此注明实际的测试位置(例如在哪个实体文件中使用哪个测试函数).
|
|
|
|
|
+*/
|
|
|
|
|
+
|
|
|
|
|
+// #ifndef __KERNEL__
|
|
|
|
|
+#if 0
|
|
|
|
|
+/**************************************************************************
|
|
|
|
|
+函数名称: main()
|
|
|
|
|
+函数版本:1.00
|
|
|
|
|
+作者: ygl
|
|
|
|
|
+创建日期:2022.03.10
|
|
|
|
|
+函数功能说明:小波测试函数,在vscode当前文件下测试。使用D:\\MinGW\\bin\\gcc.exe,需要配置vscode相关"launch.json"和"task.json“
|
|
|
|
|
+ 支持:haar db1~5 sym2~5 coif1~5 共四类小波函数
|
|
|
|
|
+输入参数:
|
|
|
|
|
+输出参数:
|
|
|
|
|
+返回值:
|
|
|
|
|
+***************************************************************************/
|
|
|
|
|
+int main() {
|
|
|
|
|
+ wave_object obj;
|
|
|
|
|
+ wt_object wt;
|
|
|
|
|
+ float *inp, *out, *diff;
|
|
|
|
|
+ int N, i, J;
|
|
|
|
|
+
|
|
|
|
|
+ FILE *ifp;
|
|
|
|
|
+ FILE *fpAppCoef = NULL;
|
|
|
|
|
+ float temp[10000];
|
|
|
|
|
+
|
|
|
|
|
+ char *name = "haar";
|
|
|
|
|
+ obj = wave_init(name);
|
|
|
|
|
+ // wave_summary(obj);
|
|
|
|
|
+
|
|
|
|
|
+ // ifp = fopen("signal.txt", "r");
|
|
|
|
|
+ ifp = fopen("myDWTData4.txt", "r");
|
|
|
|
|
+ fseek(ifp,0,SEEK_END);
|
|
|
|
|
+ N = ftell(ifp)/7;
|
|
|
|
|
+ fseek(ifp,0,SEEK_SET);
|
|
|
|
|
+ i = 0;
|
|
|
|
|
+ if (!ifp) {
|
|
|
|
|
+ printf("Cannot Open File");
|
|
|
|
|
+ exit(100);
|
|
|
|
|
+ }
|
|
|
|
|
+ while (!feof(ifp)) {
|
|
|
|
|
+ fscanf(ifp, "%f \n", &temp[i]);
|
|
|
|
|
+ i++;
|
|
|
|
|
+ }
|
|
|
|
|
+
|
|
|
|
|
+ //***********************************把注释取消掉,可以将变换结果存于文件****************************************
|
|
|
|
|
+ if(NULL == (fpAppCoef = fopen("DWT_AppData.txt", "w"))) //近似系数存储文件
|
|
|
|
|
+ {
|
|
|
|
|
+ printf("write file error!\n");
|
|
|
|
|
+ return -1;
|
|
|
|
|
+ }
|
|
|
|
|
+ //***********************************************************************************************************/
|
|
|
|
|
+ // N = 177;
|
|
|
|
|
+
|
|
|
|
|
+ fclose(ifp);
|
|
|
|
|
+
|
|
|
|
|
+ inp = (float*)malloc(sizeof(float)* N);
|
|
|
|
|
+ out = (float*)malloc(sizeof(float)* N);
|
|
|
|
|
+ diff = (float*)malloc(sizeof(float)* N);
|
|
|
|
|
+ //wmean = mean(temp, N);
|
|
|
|
|
+
|
|
|
|
|
+ for (i = 0; i < N; ++i) {
|
|
|
|
|
+ inp[i] = temp[i];
|
|
|
|
|
+ //printf("%g \n",inp[i]);
|
|
|
|
|
+ }
|
|
|
|
|
+ J = 4;
|
|
|
|
|
+
|
|
|
|
|
+ wt = wt_init(obj, "modwt", N, J);// Initialize the wavelet transform object
|
|
|
|
|
+
|
|
|
|
|
+ modwt(wt, inp);// Perform MODWT
|
|
|
|
|
+ //MODWT output can be accessed using wt->output vector. Use wt_summary to find out how to extract appx and detail coefficients
|
|
|
|
|
+
|
|
|
|
|
+ for (i = 0; i < wt->outlength; ++i) {
|
|
|
|
|
+ printf("%g ",wt->output[i]);
|
|
|
|
|
+ fprintf(fpAppCoef, "%f\t", wt->output[i]); //将数据流写入文件
|
|
|
|
|
+ }
|
|
|
|
|
+
|
|
|
|
|
+ // imodwt(wt, out);// Perform ISWT (if needed)
|
|
|
|
|
+ // Test Reconstruction
|
|
|
|
|
+
|
|
|
|
|
+ // for (i = 0; i < wt->siglength; ++i) {
|
|
|
|
|
+ // diff[i] = out[i] - inp[i];
|
|
|
|
|
+ // }
|
|
|
|
|
+
|
|
|
|
|
+
|
|
|
|
|
+ fclose(fpAppCoef);
|
|
|
|
|
+ wave_free(obj);
|
|
|
|
|
+ wt_free(wt);
|
|
|
|
|
+
|
|
|
|
|
+ rt_free(inp);
|
|
|
|
|
+ rt_free(out);
|
|
|
|
|
+ rt_free(diff);
|
|
|
|
|
+ return 0;
|
|
|
|
|
+}
|
|
|
|
|
+
|
|
|
|
|
+#endif
|
|
|
|
|
+/*------------------------------ 文件结束 -------------------------------------
|
|
|
|
|
+*/
|
|
|
|
|
+
|
|
|
|
|
+
|