一般的讲数字信号处理的书中都会提到窗函数。大多数只会提及其中的几种。这里我把这些窗都用C语言实现了一下,都不复杂,但如果要自己去弄也挺费时间。所有函数都用Matlab验证了。包括以下窗:
1 /*窗类型*/ 2 typedef enum 3 { 4 Bartlett = 0, 5 BartLettHann, 6 BlackMan, 7 BlackManHarris, 8 Bohman, 9 Chebyshev, 10 FlatTop, 11 Gaussian, 12 Hamming, 13 Hann, 14 Kaiser, 15 Nuttal, 16 Parzen, 17 Rectangular, 18 Taylor, 19 Triangular, 20 Tukey 21 }winType;
别的不多说了,直接上干货。
1 /* 2 *file WindowFunction.h 3 *author Vincent Cui 4 *e-mail [email protected] 5 *version 0.3 6 *data 31-Oct-2014 7 *brief 各种窗函数的C语言实现 8 */ 9 10 11 #ifndef _WINDOWFUNCTION_H_ 12 #define _WINDOWFUNCTION_H_ 13 14 #include "GeneralConfig.h" 15 16 #define BESSELI_K_LENGTH 10 17 18 #define FLATTOPWIN_A0 0.215578995 19 #define FLATTOPWIN_A1 0.41663158 20 #define FLATTOPWIN_A2 0.277263158 21 #define FLATTOPWIN_A3 0.083578947 22 #define FLATTOPWIN_A4 0.006947368 23 24 #define NUTTALL_A0 0.3635819 25 #define NUTTALL_A1 0.4891775 26 #define NUTTALL_A2 0.1365995 27 #define NUTTALL_A3 0.0106411 28 29 #define BLACKMANHARRIS_A0 0.35875 30 #define BLACKMANHARRIS_A1 0.48829 31 #define BLACKMANHARRIS_A2 0.14128 32 #define BLACKMANHARRIS_A3 0.01168 33 34 35 36 dspErrorStatus taylorWin(dspUint_16 N, dspUint_16 nbar, dspDouble sll, dspDouble **w); 37 dspErrorStatus triangularWin(dspUint_16 N, dspDouble **w); 38 dspErrorStatus tukeyWin(dspUint_16 N, dspDouble r, dspDouble **w); 39 dspErrorStatus bartlettWin(dspUint_16 N, dspDouble **w); 40 dspErrorStatus bartLettHannWin(dspUint_16 N, dspDouble **w); 41 dspErrorStatus blackManWin(dspUint_16 N, dspDouble **w); 42 dspErrorStatus blackManHarrisWin(dspUint_16 N, dspDouble **w); 43 dspErrorStatus bohmanWin(dspUint_16 N, dspDouble **w); 44 dspErrorStatus chebyshevWin(dspUint_16 N, dspDouble r, dspDouble **w); 45 dspErrorStatus flatTopWin(dspUint_16 N, dspDouble **w); 46 dspErrorStatus gaussianWin(dspUint_16 N, dspDouble alpha, dspDouble **w); 47 dspErrorStatus hammingWin(dspUint_16 N, dspDouble **w); 48 dspErrorStatus hannWin(dspUint_16 N, dspDouble **w); 49 dspErrorStatus kaiserWin(dspUint_16 N, dspDouble beta, dspDouble **w); 50 dspErrorStatus nuttalWin(dspUint_16 N, dspDouble **w); 51 dspErrorStatus parzenWin(dspUint_16 N, dspDouble **w); 52 dspErrorStatus rectangularWin(dspUint_16 N, dspDouble **w); 53 54 55 56 #endif
WindowFunction.h
1 /* 2 *file WindowFunction.c 3 *author Vincent Cui 4 *e-mail [email protected] 5 *version 0.3 6 *data 31-Oct-2014 7 *brief 各种窗函数的C语言实现 8 */ 9 10 11 #include "WindowFunction.h" 12 #include "GeneralConfig.h" 13 #include "MathReplenish.h" 14 #include "math.h" 15 #include <stdlib.h> 16 #include <stdio.h> 17 #include <string.h> 18 19 /*函数名:taylorWin 20 *说明:计算泰勒窗。泰勒加权函数 21 *输入: 22 *输出: 23 *返回: 24 *调用:prod()连乘函数 25 *其它:用过以后,需要手动释放掉*w的内存空间 26 * 调用示例:ret = taylorWin(99, 4, 40, &w); 注意此处的40是正数 表示-40dB 27 */ 28 dspErrorStatus taylorWin(dspUint_16 N, dspUint_16 nbar, dspDouble sll, dspDouble **w) 29 { 30 dspDouble A; 31 dspDouble *retDspDouble; 32 dspDouble *sf; 33 dspDouble *result; 34 dspDouble alpha,beta,theta; 35 dspUint_16 i,j; 36 37 /*A = R cosh(PI, A) = R*/ 38 A = (dspDouble)acosh(pow((dspDouble)10.0,(dspDouble)sll/20.0)) / PI; 39 A = A * A; 40 41 /*开出存放系数的空间*/ 42 retDspDouble = (dspDouble *)malloc(sizeof(dspDouble) * (nbar - 1)); 43 if(retDspDouble == NULL) 44 return DSP_ERROR; 45 sf = retDspDouble; 46 47 /*开出存放系数的空间*/ 48 retDspDouble = (dspDouble *)malloc(sizeof(dspDouble) * N); 49 if(retDspDouble == NULL) 50 return DSP_ERROR; 51 result = retDspDouble; 52 53 alpha = prod(1, 1, (nbar - 1)); 54 alpha *= alpha; 55 beta = (dspDouble)nbar / sqrt( A + pow((nbar - 0.5), 2) ); 56 for(i = 1; i <= (nbar - 1); i++) 57 { 58 *(sf + i - 1) = prod(1,1,(nbar -1 + i)) * prod(1,1,(nbar -1 - i)); 59 theta = 1; 60 for(j = 1; j <= (nbar - 1); j++) 61 { 62 theta *= 1 - (dspDouble)(i * i) / ( beta * beta * ( A + (j - 0.5) * (j - 0.5)) ); 63 } 64 *(sf + i - 1) = alpha * (dspDouble)theta / (*(sf + i - 1)); 65 } 66 67 /*奇数阶*/ 68 if((N % 2) == 1) 69 { 70 for(i = 0; i < N; i++) 71 { 72 alpha = 0; 73 for(j = 1; j <= (nbar - 1); j++) 74 { 75 alpha += (*(sf + j - 1)) * cos( 2 * PI * j * (dspDouble)(i - ((N-1)/2))/N ); 76 } 77 *(result + i) = 1 + 2 * alpha; 78 } 79 } 80 /*偶数阶*/ 81 else 82 { 83 for(i = 0; i < N; i++) 84 { 85 alpha = 0; 86 for(j = 1; j <= (nbar - 1); j++) 87 { 88 alpha += (*(sf + j - 1)) * cos( PI * j * (dspDouble)(2 * (i - (N/2)) + 1) / N ); 89 } 90 *(result + i) = 1 + 2 * alpha; 91 92 } 93 } 94 *w = result; 95 free(sf); 96 97 return DSP_SUCESS; 98 99 } 100 101 102 /* 103 *函数名:triangularWin 104 *说明:计算三角窗函数 105 *输入: 106 *输出: 107 *返回: 108 *调用: 109 *其它:用过以后,需要手动释放掉*w的内存空间 110 * 调用示例:ret = triangularWin(99, &w); 111 */ 112 dspErrorStatus triangularWin(dspUint_16 N, dspDouble **w) 113 { 114 dspDouble *ptr; 115 dspUint_16 i; 116 117 ptr = (dspDouble *)malloc(N * sizeof(dspDouble)); 118 if(ptr == NULL) 119 return DSP_ERROR; 120 121 122 /*阶数为奇*/ 123 if((N % 2) == 1) 124 { 125 for(i = 0; i < ((N - 1)/2); i++) 126 { 127 *(ptr + i) = 2 * (dspDouble)(i + 1) / (N + 1); 128 } 129 for(i = ((N - 1)/2); i < N; i++) 130 { 131 *(ptr + i) = 2 * (dspDouble)(N - i) / (N + 1); 132 } 133 } 134 /*阶数为偶*/ 135 else 136 { 137 for(i = 0; i < (N/2); i++) 138 { 139 *(ptr + i) = (i + i + 1) * (dspDouble)1 / N; 140 } 141 for(i = (N/2); i < N; i++) 142 { 143 *(ptr + i) = *(ptr + N - 1 - i); 144 } 145 } 146 *w = ptr; 147 148 return DSP_SUCESS; 149 } 150 151 /* 152 *函数名:tukeyWin 153 *说明:计算tukey窗函数 154 *输入: 155 *输出: 156 *返回:linSpace() 157 *调用: 158 *其它:用过以后,需要手动释放掉*w的内存空间 159 * 调用示例:ret = tukeyWin(99, 0.5, &w); 160 */ 161 dspErrorStatus tukeyWin(dspUint_16 N, dspDouble r, dspDouble **w) 162 { 163 dspErrorStatus retErrorStatus; 164 dspUint_16 index; 165 dspDouble *x,*result,*retPtr; 166 dspDouble alpha; 167 168 retErrorStatus = linSpace(0, 1, N, &x); 169 if(retErrorStatus == DSP_ERROR) 170 return DSP_ERROR; 171 172 result = (dspDouble *)malloc(N * sizeof(dspDouble)); 173 if(result == NULL) 174 return DSP_ERROR; 175 176 /*r <= 0 就是矩形窗*/ 177 if(r <= 0) 178 { 179 retErrorStatus = rectangularWin(N, &retPtr); 180 if(retErrorStatus == DSP_ERROR) 181 return DSP_ERROR; 182 /*将数据拷出来以后,释放调用的窗函数的空间*/ 183 memcpy(result, retPtr, ( N * sizeof(dspDouble))); 184 free(retPtr); 185 } 186 /*r >= 1 就是汉宁窗*/ 187 else if(r >= 1) 188 { 189 retErrorStatus = hannWin(N, &retPtr); 190 if(retErrorStatus == DSP_ERROR) 191 return DSP_ERROR; 192 /*将数据拷出来以后,释放调用的窗函数的空间*/ 193 memcpy(result, retPtr, ( N * sizeof(dspDouble))); 194 free(retPtr); 195 } 196 else 197 { 198 for(index = 0; index < N; index++) 199 { 200 alpha = *(x + index); 201 if(alpha < (r/2)) 202 { 203 *(result + index) = (dspDouble)(1 + cos( 2 * PI * (dspDouble)(alpha - (dspDouble)r/2)/r))/2; 204 } 205 else if((alpha >= (r/2)) && (alpha <(1 - r/2))) 206 { 207 *(result + index) = 1; 208 } 209 else 210 { 211 *(result + index) = (dspDouble)(1 + cos( 2 * PI * (dspDouble)(alpha - 1 + (dspDouble)r/2)/r))/2; 212 } 213 214 } 215 } 216 217 free(x); 218 219 *w = result; 220 221 return DSP_SUCESS; 222 223 } 224 225 /* 226 *函数名:bartlettWin 227 *说明:计算bartlettWin窗函数 228 *输入: 229 *输出: 230 *返回: 231 *调用: 232 *其它:用过以后,需要手动释放掉*w的内存空间 233 * 调用示例:ret = bartlettWin(99, &w); 234 */ 235 dspErrorStatus bartlettWin(dspUint_16 N, dspDouble **w) 236 { 237 dspDouble *ret; 238 dspUint_16 n; 239 240 ret = (dspDouble *)malloc(N * sizeof(dspDouble)); 241 if(ret == NULL) 242 return DSP_ERROR; 243 244 for(n = 0; n < ( N - 1 ) / 2; n++) 245 { 246 *(ret + n) = 2 * (dspDouble)n / (N - 1); 247 } 248 249 for(n = ( N - 1 ) / 2; n < N; n++) 250 { 251 *(ret + n) = 2 - 2 * (dspDouble)n / (( N - 1 )); 252 } 253 254 *w = ret; 255 256 return DSP_SUCESS; 257 } 258 259 /* 260 *函数名:bartLettHannWin 261 *说明:计算bartLettHannWin窗函数 262 *输入: 263 *输出: 264 *返回: 265 *调用: 266 *其它:用过以后,需要手动释放掉*w的内存空间 267 * 调用示例:ret = bartLettHannWin(99, &w); 268 */ 269 dspErrorStatus bartLettHannWin(dspUint_16 N, dspDouble **w) 270 { 271 dspUint_16 n; 272 dspDouble *ret; 273 274 ret = (dspDouble *)malloc(N * sizeof(dspDouble)); 275 if(ret == NULL) 276 return DSP_ERROR; 277 /*奇*/ 278 if(( N % 2 ) == 1) 279 { 280 for(n = 0; n < N; n++) 281 { 282 *(ret + n) = 0.62 - 0.48 * myAbs( ( (dspDouble)n / ( N - 1 ) ) - 0.5 ) + 0.38 * cos( 2 * PI * ( ((dspDouble)n / ( N - 1 ) ) - 0.5 ) ); 283 } 284 for(n = 0; n < (N-1)/2; n++) 285 { 286 *(ret + n) = *(ret + N - 1 - n); 287 } 288 } 289 /*偶*/ 290 else 291 { 292 for(n = 0; n < N; n++) 293 { 294 *(ret + n) = 0.62 - 0.48 * myAbs( ( (dspDouble)n / ( N - 1 ) ) - 0.5 ) + 0.38 * cos( 2 * PI * ( ((dspDouble)n / ( N - 1 ) ) - 0.5 ) ); 295 } 296 for(n = 0; n < N/2; n++) 297 { 298 *(ret + n) = *(ret + N -1 - n); 299 } 300 } 301 302 *w = ret; 303 304 return DSP_SUCESS; 305 306 } 307 308 /* 309 *函数名:blackManWin 310 *说明:计算blackManWin窗函数 311 *输入: 312 *输出: 313 *返回: 314 *调用: 315 *其它:用过以后,需要手动释放掉*w的内存空间 316 * 调用示例:ret = blackManWin(99, &w); 317 */ 318 dspErrorStatus blackManWin(dspUint_16 N, dspDouble **w) 319 { 320 dspUint_16 n; 321 dspDouble *ret; 322 ret = (dspDouble *)malloc(N * sizeof(dspDouble)); 323 if(ret == NULL) 324 return DSP_ERROR; 325 326 for(n = 0; n < N; n++) 327 { 328 *(ret + n) = 0.42 - 0.5 * cos(2 * PI * (dspDouble)n / ( N - 1 )) + 0.08 * cos( 4 * PI * ( dspDouble )n / ( N - 1 ) ); 329 } 330 331 *w = ret; 332 333 return DSP_SUCESS; 334 } 335 336 /* 337 *函数名:blackManHarrisWin 338 *说明:计算blackManHarrisWin窗函数 339 *输入: 340 *输出: 341 *返回: 342 *调用: 343 *其它:用过以后,需要手动释放掉*w的内存空间 344 * 调用示例:ret = blackManHarrisWin(99, &w); 345 * minimum 4-term Blackman-harris window -- From Matlab 346 */ 347 dspErrorStatus blackManHarrisWin(dspUint_16 N, dspDouble **w) 348 { 349 dspUint_16 n; 350 dspDouble *ret; 351 352 ret = (dspDouble *)malloc(N * sizeof(dspDouble)); 353 if(ret == NULL) 354 return DSP_ERROR; 355 356 for(n = 0; n < N; n++) 357 { 358 *(ret + n) = BLACKMANHARRIS_A0 - BLACKMANHARRIS_A1 * cos( 2 * PI * (dspDouble)n / (N) ) + 359 BLACKMANHARRIS_A2 * cos( 4 * PI * (dspDouble)n/ (N) ) - 360 BLACKMANHARRIS_A3 * cos( 6 * PI * (dspDouble)n/ (N) ); 361 } 362 363 *w = ret; 364 365 return DSP_SUCESS; 366 } 367 368 /* 369 *函数名:bohmanWin 370 *说明:计算bohmanWin窗函数 371 *输入: 372 *输出: 373 *返回: 374 *调用: 375 *其它:用过以后,需要手动释放掉*w的内存空间 376 * 调用示例:ret = bohmanWin(99, &w); 377 */ 378 dspErrorStatus bohmanWin(dspUint_16 N, dspDouble **w) 379 { 380 dspUint_16 n; 381 dspDouble *ret; 382 dspDouble x; 383 ret = (dspDouble *)malloc(N * sizeof(dspDouble)); 384 if(ret == NULL) 385 return DSP_ERROR; 386 387 for(n = 0; n < N; n++) 388 { 389 x = -1 + n * (dspDouble)2 / ( N - 1 ) ; 390 /*取绝对值*/ 391 x = x >= 0 ? x : ( x * ( -1 ) ); 392 *(ret + n) = ( 1 - x ) * cos( PI * x) + (dspDouble)(1 / PI) * sin( PI * x); 393 } 394 395 *w = ret; 396 397 return DSP_SUCESS; 398 } 399 400 /* 401 *函数名:chebyshevWin 402 *说明:计算chebyshevWin窗函数 403 *输入: 404 *输出: 405 *返回: 406 *调用: 407 *其它:用过以后,需要手动释放掉*w的内存空间 408 * 调用示例:ret = chebyshevWin(99,100, &w); 409 */ 410 dspErrorStatus chebyshevWin(dspUint_16 N, dspDouble r, dspDouble **w) 411 { 412 dspUint_16 n,index; 413 dspDouble *ret; 414 dspDouble x, alpha, beta, theta, gama; 415 416 ret = (dspDouble *)malloc(N * sizeof(dspDouble)); 417 if(ret == NULL) 418 return DSP_ERROR; 419 420 421 /*10^(r/20)*/ 422 theta = pow((dspDouble)10, (dspDouble)(myAbs(r)/20)); 423 beta = pow(cosh(acosh(theta)/(N - 1)),2); 424 alpha = 1 - (dspDouble)1 / beta; 425 426 if((N % 2) == 1) 427 { 428 /*计算一半的区间*/ 429 for( n = 1; n < ( N + 1 ) / 2; n++ ) 430 { 431 gama = 1; 432 for(index = 1; index < n; index++) 433 { 434 x = index * (dspDouble)( N - 1 - 2 * n + index) /(( n - index ) * (n + 1 -index)); 435 gama = gama * alpha * x + 1; 436 } 437 *(ret + n) = (N - 1) * alpha * gama; 438 } 439 440 theta = *( ret + (N - 1)/2 ); 441 *ret = 1; 442 443 for(n = 0; n < ( N + 1 ) / 2; n++ ) 444 { 445 *(ret + n) = (dspDouble)(*(ret + n)) / theta; 446 } 447 448 /*填充另一半*/ 449 for(; n < N; n++) 450 { 451 *(ret + n) = ret[N - n - 1]; 452 } 453 } 454 else 455 { 456 /*计算一半的区间*/ 457 for( n = 1; n < ( N + 1 ) / 2; n++ ) 458 { 459 gama = 1; 460 for(index = 1; index < n; index++) 461 { 462 x = index * (dspDouble)( N - 1 - 2 * n + index) /(( n - index ) * (n + 1 -index)); 463 gama = gama * alpha * x + 1; 464 } 465 *(ret + n) = (N - 1) * alpha * gama; 466 } 467 468 theta = *( ret + (N/2) - 1); 469 *ret = 1; 470 471 for(n = 0; n < ( N + 1 ) / 2; n++ ) 472 { 473 *(ret + n) = (dspDouble)(*(ret + n)) / theta; 474 } 475 476 /*填充另一半*/ 477 for(; n < N; n++) 478 { 479 *(ret + n) = ret[N - n - 1]; 480 } 481 } 482 483 484 *w = ret; 485 486 return DSP_SUCESS; 487 } 488 489 /* 490 *函数名:flatTopWin 491 *说明:计算flatTopWin窗函数 492 *输入: 493 *输出: 494 *返回: 495 *调用: 496 *其它:用过以后,需要手动释放掉*w的内存空间 497 * 调用示例:ret = flatTopWin(99, &w); 498 */ 499 dspErrorStatus flatTopWin(dspUint_16 N, dspDouble **w) 500 { 501 dspUint_16 n; 502 dspDouble *ret; 503 ret = (dspDouble *)malloc(N * sizeof(dspDouble)); 504 if(ret == NULL) 505 return DSP_ERROR; 506 507 for(n = 0; n < N; n++) 508 { 509 *(ret + n) = FLATTOPWIN_A0 - FLATTOPWIN_A1 * cos(2 * PI * (dspDouble)n / (N - 1)) +510 FLATTOPWIN_A2 * cos(4 * PI * (dspDouble)n / (N - 1)) -511 FLATTOPWIN_A3 * cos(6 * PI * (dspDouble)n / (N - 1)) +512 FLATTOPWIN_A4 * cos(8 * PI * (dspDouble)n / (N - 1)); 513 } 514 515 *w = ret; 516 517 return DSP_SUCESS; 518 } 519 520 521 /* 522 *函数名:gaussianWin 523 *说明:计算gaussianWin窗函数 524 *输入: 525 *输出: 526 *返回: 527 *调用: 528 *其它:用过以后,需要手动释放掉*w的内存空间 529 * 调用示例:ret = gaussianWin(99,2.5, &w); 530 */ 531 dspErrorStatus gaussianWin(dspUint_16 N, dspDouble alpha, dspDouble **w) 532 { 533 dspUint_16 n; 534 dspDouble k, beta, theta; 535 dspDouble *ret; 536 537 ret = (dspDouble *)malloc(N * sizeof(dspDouble)); 538 if(ret == NULL) 539 return DSP_ERROR; 540 541 for(n =0; n < N; n++) 542 { 543 if((N % 2) == 1) 544 { 545 k = n - (N - 1)/2; 546 beta = 2 * alpha * (dspDouble)k / (N - 1); 547 } 548 else 549 { 550 k = n - (N)/2; 551 beta = 2 * alpha * (dspDouble)k / (N - 1); 552 } 553 554 theta = pow(beta, 2); 555 *(ret + n) = exp((-1) * (dspDouble)theta / 2); 556 } 557 558 *w = ret; 559 560 return DSP_SUCESS; 561 } 562 563 /* 564 *函数名:hammingWin 565 *说明:计算hammingWin窗函数 566 *输入: 567 *输出: 568 *返回: 569 *调用: 570 *其它:用过以后,需要手动释放掉*w的内存空间 571 * 调用示例:ret = hammingWin(99, &w); 572 */ 573 dspErrorStatus hammingWin(dspUint_16 N, dspDouble **w) 574 { 575 dspUint_16 n; 576 dspDouble *ret; 577 ret = (dspDouble *)malloc(N * sizeof(dspDouble)); 578 if(ret == NULL) 579 return DSP_ERROR; 580 581 for(n = 0; n < N; n++) 582 { 583 *(ret + n) = 0.54 - 0.46 * cos (2 * PI * ( dspDouble )n / ( N - 1 ) ); 584 } 585 586 *w = ret; 587 588 return DSP_SUCESS; 589 } 590 591 /* 592 *函数名:hannWin 593 *说明:计算hannWin窗函数 594 *输入: 595 *输出: 596 *返回: 597 *调用: 598 *其它:用过以后,需要手动释放掉*w的内存空间 599 * 调用示例:ret = hannWin(99, &w); 600 */ 601 dspErrorStatus hannWin(dspUint_16 N, dspDouble **w) 602 { 603 dspUint_16 n; 604 dspDouble *ret; 605 ret = (dspDouble *)malloc(N * sizeof(dspDouble)); 606 if(ret == NULL) 607 return DSP_ERROR; 608 609 for(n = 0; n < N; n++) 610 { 611 *(ret + n) = 0.5 * ( 1 - cos( 2 * PI * (dspDouble)n / (N - 1))); 612 } 613 614 *w = ret; 615 616 return DSP_SUCESS; 617 } 618 619 /* 620 *函数名:kaiserWin 621 *说明:计算kaiserWin窗函数 622 *输入: 623 *输出: 624 *返回: 625 *调用:besseli()第一类修正贝塞尔函数 626 *其它:用过以后,需要手动释放掉*w的内存空间 627 * 调用示例:ret = kaiserWin(99, 5, &w); 628 */ 629 dspErrorStatus kaiserWin(dspUint_16 N, dspDouble beta, dspDouble **w) 630 { 631 dspUint_16 n; 632 dspDouble *ret; 633 dspDouble theta; 634 635 ret = (dspDouble *)malloc(N * sizeof(dspDouble)); 636 if(ret == NULL) 637 return DSP_ERROR; 638 639 for(n = 0; n < N; n++) 640 { 641 theta = beta * sqrt( 1 - pow( ( (2 * (dspDouble)n/(N -1)) - 1),2 ) ); 642 *(ret + n) = (dspDouble)besseli(0, theta, BESSELI_K_LENGTH) / besseli(0, beta, BESSELI_K_LENGTH); 643 } 644 645 *w = ret; 646 647 return DSP_SUCESS; 648 } 649 650 /* 651 *函数名:nuttalWin 652 *说明:计算nuttalWin窗函数 653 *输入: 654 *输出: 655 *返回: 656 *调用: 657 *其它:用过以后,需要手动释放掉*w的内存空间 658 * 调用示例:ret = nuttalWin(99, &w); 659 */ 660 dspErrorStatus nuttalWin(dspUint_16 N, dspDouble **w) 661 { 662 dspUint_16 n; 663 dspDouble *ret; 664 665 ret = (dspDouble *)malloc(N * sizeof(dspDouble)); 666 if(ret == NULL) 667 return DSP_ERROR; 668 669 for(n = 0; n < N; n++) 670 { 671 *(ret + n) =NUTTALL_A0 - NUTTALL_A1 * cos(2 * PI * (dspDouble)n / (N - 1)) +672 NUTTALL_A2 * cos(4 * PI * (dspDouble)n / (N - 1)) -673 NUTTALL_A3 * cos(6 * PI * (dspDouble)n / (N - 1)); 674 675 } 676 677 *w = ret; 678 679 return DSP_SUCESS; 680 } 681 682 /* 683 *函数名:parzenWin 684 *说明:计算parzenWin窗函数 685 *输入: 686 *输出: 687 *返回: 688 *调用: 689 *其它:用过以后,需要手动释放掉*w的内存空间 690 * 调用示例:ret = parzenWin(99, &w); 691 */ 692 dspErrorStatus parzenWin(dspUint_16 N, dspDouble **w) 693 { 694 dspUint_16 n; 695 dspDouble *ret; 696 dspDouble alpha,k; 697 698 ret = (dspDouble *)malloc(N * sizeof(dspDouble)); 699 if(ret == NULL) 700 return DSP_ERROR; 701 702 if(( N % 2) == 1) 703 { 704 for(n = 0; n < N; n++) 705 { 706 k = n - (N - 1) / 2; 707 alpha = 2 * (dspDouble)myAbs(k) / N; 708 if(myAbs(k) <= (N - 1) / 4) 709 { 710 *(ret + n) = 1 - 6 * pow(alpha,2) + 6 * pow(alpha, 3); 711 } 712 else 713 { 714 *(ret + n) = 2 * pow( (1 - alpha), 3 ); 715 } 716 717 } 718 } 719 else 720 { 721 for(n = 0; n < N; n++) 722 { 723 k = n - (N - 1) / 2; 724 alpha = 2 * (dspDouble)myAbs(k) / N; 725 if(myAbs(k) <= (dspDouble)(N -1) / 4) 726 { 727 *(ret + n) = 1 - 6 * pow(alpha,2) + 6 * pow(alpha, 3); 728 } 729 else 730 { 731 *(ret + n) = 2 * pow( (1 - alpha), 3 ); 732 } 733 734 } 735 } 736 737 738 739 *w = ret; 740 741 return DSP_SUCESS; 742 } 743 744 /* 745 *函数名:rectangularWin 746 *说明:计算rectangularWin窗函数 747 *输入: 748 *输出: 749 *返回: 750 *调用: 751 *其它:用过以后,需要手动释放掉*w的内存空间 752 * 调用示例:ret = rectangularWin(99, &w); 753 */ 754 dspErrorStatus rectangularWin(dspUint_16 N, dspDouble **w) 755 { 756 dspUint_16 n; 757 dspDouble *ret; 758 759 ret = (dspDouble *)malloc(N * sizeof(dspDouble)); 760 if(ret == NULL) 761 return DSP_ERROR; 762 763 for(n = 0; n < N; n++) 764 { 765 *(ret + n) = 1; 766 } 767 768 *w = ret; 769 770 return DSP_SUCESS; 771 }
WindowFunction.c
欢迎多交流!
时间: 2024-10-09 08:23:51