From f58bb50384c0443c123b8acf886b36c8e7eca9ce Mon Sep 17 00:00:00 2001 From: zhaohe Date: Tue, 14 May 2024 13:45:58 +0800 Subject: [PATCH] update --- app/app.uvoptx | 24 ++++ app/app.uvprojx | 20 ++++ app/src/app_basic_service/basic/event.h | 2 +- .../app_service/ecg_service/algo/MedianFilter.c | 88 ++++++++++++++ .../app_service/ecg_service/algo/MedianFilter.h | 38 ++++++ .../ecg_service/algo/iflytop_simple_filter.c | 45 ++++++- .../ecg_service/algo/iflytop_simple_filter.h | 13 +++ .../ecg_service/algo/iflytop_simple_filter_ext.h | 5 +- app/src/app_service/ecg_service/ecg_algo.c | 114 +++++++++++------- app/src/app_service/ecg_service/ecg_algo.h | 3 +- app/src/app_service/ecg_service/ecg_service.c | 5 +- app/src/aproject_config/config.h | 5 +- app/src/ble_data_processer_utils.c | 10 +- ify_hrs_protocol | 2 +- test/MedianFilter.c | 88 ++++++++++++++ test/MedianFilter.h | 38 ++++++ test/README.md | 1 + test/test.c | 27 +++++ test/test.py | 129 +++++++++++++++++++++ 19 files changed, 595 insertions(+), 62 deletions(-) create mode 100644 app/src/app_service/ecg_service/algo/MedianFilter.c create mode 100644 app/src/app_service/ecg_service/algo/MedianFilter.h create mode 100644 test/MedianFilter.c create mode 100644 test/MedianFilter.h create mode 100644 test/README.md create mode 100644 test/test.c create mode 100644 test/test.py diff --git a/app/app.uvoptx b/app/app.uvoptx index d4029d7..af01a7a 100644 --- a/app/app.uvoptx +++ b/app/app.uvoptx @@ -1787,6 +1787,30 @@ 0 0 + + 15 + 110 + 1 + 0 + 0 + 0 + .\src\app_service\ecg_service\algo\MedianFilter.c + MedianFilter.c + 0 + 0 + + + 15 + 111 + 5 + 0 + 0 + 0 + .\src\app_service\ecg_service\algo\MedianFilter.h + MedianFilter.h + 0 + 0 + diff --git a/app/app.uvprojx b/app/app.uvprojx index c6e96d1..0991492 100644 --- a/app/app.uvprojx +++ b/app/app.uvprojx @@ -3905,6 +3905,16 @@ 1 .\src\app_service\ecg_service\algo\zdata_statistics.c + + MedianFilter.c + 1 + .\src\app_service\ecg_service\algo\MedianFilter.c + + + MedianFilter.h + 5 + .\src\app_service\ecg_service\algo\MedianFilter.h + @@ -7814,6 +7824,16 @@ 1 .\src\app_service\ecg_service\algo\zdata_statistics.c + + MedianFilter.c + 1 + .\src\app_service\ecg_service\algo\MedianFilter.c + + + MedianFilter.h + 5 + .\src\app_service\ecg_service\algo\MedianFilter.h + diff --git a/app/src/app_basic_service/basic/event.h b/app/src/app_basic_service/basic/event.h index 90d2142..5d1deaa 100644 --- a/app/src/app_basic_service/basic/event.h +++ b/app/src/app_basic_service/basic/event.h @@ -5,7 +5,7 @@ #include "aproject_config/config.h" #include "ify_hrs_protocol/heart_rate_sensor_protocol.h" -#define ECG_DATA_REPORT_FRAME_NUM 64 // ecg每次上报的帧数 +#define ECG_DATA_REPORT_FRAME_NUM 50 // ecg每次上报的帧数 typedef enum { kappevent_tmr_1s_scheduler_event, diff --git a/app/src/app_service/ecg_service/algo/MedianFilter.c b/app/src/app_service/ecg_service/algo/MedianFilter.c new file mode 100644 index 0000000..a785900 --- /dev/null +++ b/app/src/app_service/ecg_service/algo/MedianFilter.c @@ -0,0 +1,88 @@ +/* + * MedianFilter.c + * + * Created on: May 19, 2018 + * Author: alexandru.bogdan + */ + +#include "MedianFilter.h" + +int MEDIANFILTER_Init(sMedianFilter_t *medianFilter) +{ + if(medianFilter && medianFilter->medianBuffer && + (medianFilter->numNodes % 2) && (medianFilter->numNodes > 1)) + { + //initialize buffer nodes + for(unsigned int i = 0; i < medianFilter->numNodes; i++) + { + medianFilter->medianBuffer[i].value = 0; + medianFilter->medianBuffer[i].nextAge = &medianFilter->medianBuffer[(i + 1) % medianFilter->numNodes]; + medianFilter->medianBuffer[i].nextValue = &medianFilter->medianBuffer[(i + 1) % medianFilter->numNodes]; + medianFilter->medianBuffer[(i + 1) % medianFilter->numNodes].prevValue = &medianFilter->medianBuffer[i]; + } + //initialize heads + medianFilter->ageHead = medianFilter->medianBuffer; + medianFilter->valueHead = medianFilter->medianBuffer; + medianFilter->medianHead = &medianFilter->medianBuffer[medianFilter->numNodes / 2]; + + return 0; + } + + return -1; +} + +int MEDIANFILTER_Insert(sMedianFilter_t *medianFilter, int sample) +{ + unsigned int i; + sMedianNode_t *newNode, *it; + + if(medianFilter->ageHead == medianFilter->valueHead) + { //if oldest node is also the smallest node, increment value head + medianFilter->valueHead = medianFilter->valueHead->nextValue; + } + + if((medianFilter->ageHead == medianFilter->medianHead) || + (medianFilter->ageHead->value > medianFilter->medianHead->value)) + { //prepare for median correction + medianFilter->medianHead = medianFilter->medianHead->prevValue; + } + + //replace age head with new sample + newNode = medianFilter->ageHead; + newNode->value = sample; + + //remove age head from list + medianFilter->ageHead->nextValue->prevValue = medianFilter->ageHead->prevValue; + medianFilter->ageHead->prevValue->nextValue = medianFilter->ageHead->nextValue; + //increment age head + medianFilter->ageHead = medianFilter->ageHead->nextAge; + + //find new node position + it = medianFilter->valueHead; //set iterator as value head + for(i = 0; i < medianFilter->numNodes - 1; i++) + { + if(sample < it->value) + { + if(i == 0) + { //replace value head if new node is the smallest + medianFilter->valueHead = newNode; + } + break; + } + it = it->nextValue; + } + + //insert new node in list + it->prevValue->nextValue = newNode; + newNode->prevValue = it->prevValue; + it->prevValue = newNode; + newNode->nextValue = it; + + //adjust median node + if(i >= (medianFilter->numNodes / 2)) + { + medianFilter->medianHead = medianFilter->medianHead->nextValue; + } + + return medianFilter->medianHead->value; +} \ No newline at end of file diff --git a/app/src/app_service/ecg_service/algo/MedianFilter.h b/app/src/app_service/ecg_service/algo/MedianFilter.h new file mode 100644 index 0000000..e616ca5 --- /dev/null +++ b/app/src/app_service/ecg_service/algo/MedianFilter.h @@ -0,0 +1,38 @@ +/* + * MedianFilter.h + * + * Created on: May 19, 2018 + * Author: alexandru.bogdan + */ + +#ifndef MEDIANFILTER_H_ +#define MEDIANFILTER_H_ + +#ifdef __cplusplus +extern "C" { +#endif + +typedef struct sMedianNode +{ + int value; //sample value + struct sMedianNode *nextAge; //pointer to next oldest value + struct sMedianNode *nextValue; //pointer to next smallest value + struct sMedianNode *prevValue; //pointer to previous smallest value +}sMedianNode_t; + +typedef struct +{ + unsigned int numNodes; //median node buffer length + sMedianNode_t *medianBuffer; //median node buffer + sMedianNode_t *ageHead; //pointer to oldest value + sMedianNode_t *valueHead; //pointer to smallest value + sMedianNode_t *medianHead; //pointer to median value +}sMedianFilter_t; + +int MEDIANFILTER_Init(sMedianFilter_t *medianFilter); +int MEDIANFILTER_Insert(sMedianFilter_t *medianFilter, int sample); + +#ifdef __cplusplus +} +#endif +#endif \ No newline at end of file diff --git a/app/src/app_service/ecg_service/algo/iflytop_simple_filter.c b/app/src/app_service/ecg_service/algo/iflytop_simple_filter.c index fcf1fec..9f9d1fd 100644 --- a/app/src/app_service/ecg_service/algo/iflytop_simple_filter.c +++ b/app/src/app_service/ecg_service/algo/iflytop_simple_filter.c @@ -206,4 +206,47 @@ float SmoothingFilter_Update(SmoothingFilter_t *filter, float vin) { } return filter->sum / filter->datanum; -} \ No newline at end of file +} + +void BaselineDriftRemoval_Init(BaselineDriftRemoval_t *filter, int windows_size, bool enable) { + if (windows_size > MEDIAN_FILTER_SIZE) { + windows_size = MEDIAN_FILTER_SIZE; + } + + filter->windows_size = windows_size; + filter->datanum = 0; + filter->enable = enable; + for (int i = 0; i < windows_size; i++) { + filter->fifo[i] = 0; + } +} +float BaselineDriftRemoval_Update(BaselineDriftRemoval_t *filter, float vin) { + if (!filter->enable) return vin; + + // memmove(filter->fifo, filter->fifo + 1, (filter->windows_size - 1) * sizeof(float)); + filter->fifo[filter->windows_size - 1] = vin; + filter->datanum++; + if (filter->datanum > filter->windows_size) { + filter->datanum = filter->windows_size; + } + if (filter->datanum < filter->windows_size) { + return 0; + } + + // 求中值 + // memcpy(filter->processcache, filter->fifo, filter->windows_size * sizeof(float)); + for (int i = 0; i < filter->windows_size; i++) { + for (int j = i + 1; j < filter->windows_size; j++) { + if (filter->processcache[i] > filter->processcache[j]) { + float temp = filter->processcache[i]; + filter->processcache[i] = filter->processcache[j]; + filter->processcache[j] = temp; + } + } + } + + float mid = filter->processcache[filter->windows_size / 2]; + float min_data_in_fifo = filter->fifo[filter->windows_size / 2]; + + return min_data_in_fifo - mid; +} diff --git a/app/src/app_service/ecg_service/algo/iflytop_simple_filter.h b/app/src/app_service/ecg_service/algo/iflytop_simple_filter.h index ffc4f12..fdff1bf 100644 --- a/app/src/app_service/ecg_service/algo/iflytop_simple_filter.h +++ b/app/src/app_service/ecg_service/algo/iflytop_simple_filter.h @@ -94,4 +94,17 @@ void SmoothingFilter_Init(SmoothingFilter_t *filter, int windows_size, bool ena float SmoothingFilter_Update(SmoothingFilter_t *filter, float vin); +#define MEDIAN_FILTER_SIZE 300 +typedef struct { + float fifo[MEDIAN_FILTER_SIZE]; + float processcache[MEDIAN_FILTER_SIZE]; + bool enable; + + int windows_size; + int datanum; +} BaselineDriftRemoval_t; + +void BaselineDriftRemoval_Init(BaselineDriftRemoval_t *filter, int windows_size, bool enable); +float BaselineDriftRemoval_Update(BaselineDriftRemoval_t *filter, float vin); + #endif \ No newline at end of file diff --git a/app/src/app_service/ecg_service/algo/iflytop_simple_filter_ext.h b/app/src/app_service/ecg_service/algo/iflytop_simple_filter_ext.h index ad126f2..ea7e4d7 100644 --- a/app/src/app_service/ecg_service/algo/iflytop_simple_filter_ext.h +++ b/app/src/app_service/ecg_service/algo/iflytop_simple_filter_ext.h @@ -25,4 +25,7 @@ typedef struct { } NOTCHFilterExt_t; void NOTCHFilterExt_init(NOTCHFilterExt_t *filter, float centerFreqHz, float notchWidthHz, float sampleTimeS, int order, bool enable); -float NOTCHFilterExt_update(NOTCHFilterExt_t *filter, float v_in); \ No newline at end of file +float NOTCHFilterExt_update(NOTCHFilterExt_t *filter, float v_in); + + + \ No newline at end of file diff --git a/app/src/app_service/ecg_service/ecg_algo.c b/app/src/app_service/ecg_service/ecg_algo.c index 0238b8b..a60bb31 100644 --- a/app/src/app_service/ecg_service/ecg_algo.c +++ b/app/src/app_service/ecg_service/ecg_algo.c @@ -1,20 +1,30 @@ #include "ecg_algo.h" +#include "algo/MedianFilter.h" #include "algo/iflytop_simple_filter_ext.h" #include "algo/zdata_statistics.h" #include "algo/zsimple_qrs.h" #include "znordic.h" -LPFilterExt_t m_lp_filter; -HPFilterExt_t m_hp_filter; -NOTCHFilterExt_t m_notch_filter; -SmoothingFilter_t m_smoothingFilter; +#define REPORT_MEDIAN_WINDOWS_SIZE 101 // 必须是奇数 + +LPFilterExt_t m_lp_filter; +HPFilterExt_t m_hp_filter; +NOTCHFilterExt_t m_notch_filter; + +SmoothingFilter_t m_report_smoothing_filter; +median_filter_t m_display_median_filter; + +static sMedianFilter_t report_medianFilter; +static sMedianNode_t report_medianBuffer[REPORT_MEDIAN_WINDOWS_SIZE]; + +void BaselineDriftRemoval_Init(BaselineDriftRemoval_t *filter, int windows_size, bool enable); +float BaselineDriftRemoval_Update(BaselineDriftRemoval_t *filter, float vin); int32_t m_data_statistics_buf[STATISTICS_BUF_SIZE]; // 心率判断 zdata_statistics_t m_data_statistics; zsimple_qrs_t m_qrs; -int32_t data_af_high_pass_filter; int32_t reportdata; int32_t displaydata; @@ -27,42 +37,66 @@ void ecg_algo_init() { int32_t maxval; int32_t minval; int32_t nowval100; -void ecg_algo_process_data(int32_t indata) { - float data = indata; - data = HPFilterExt_update(&m_hp_filter, data); - data_af_high_pass_filter = data; - - data = LPFilterExt_update(&m_lp_filter, data); - data = NOTCHFilterExt_update(&m_notch_filter, data); - data = SmoothingFilter_Update(&m_smoothingFilter, data); + +static void display_data_processer(int32_t indata) { + // float data = indata; + + // data = BaselineDriftRemoval_Update(&m_display_drift_filter, data); // 基线漂移 + // data = median_filter_update(&m_display_median_filter, data); + + // zdata_statistics_push(&m_data_statistics, data); + // if (zdata_statistics_is_full(&m_data_statistics)) { + // zsimple_qrs_process_data(&m_qrs, data, zdata_statistics_get_min(&m_data_statistics), zdata_statistics_get_max(&m_data_statistics), + // zdata_statistics_get_avg(&m_data_statistics)); + + // maxval = zdata_statistics_get_max(&m_data_statistics); + // minval = zdata_statistics_get_min(&m_data_statistics); + + // int32_t nowvaloff = data - minval; + // nowval100 = 10 + nowvaloff * 80 / (maxval - minval); + // displaydata = nowval100; + + // } else { + // maxval = 0; + // minval = 0; + // displaydata = 50; + // } +} + +static void report_data_processer(int32_t indata) { + float data = indata; + int medianValue = MEDIANFILTER_Insert(&report_medianFilter, data); + data = indata - medianValue; + // data = BaselineDriftRemoval_Update(&m_report_drift_filter, data); // 基线漂移 + data = SmoothingFilter_Update(&m_report_smoothing_filter, data); + reportdata = data; +} - zdata_statistics_push(&m_data_statistics, data); - if (zdata_statistics_is_full(&m_data_statistics)) { - zsimple_qrs_process_data(&m_qrs, data, zdata_statistics_get_min(&m_data_statistics), zdata_statistics_get_max(&m_data_statistics), - zdata_statistics_get_avg(&m_data_statistics)); - - /** - * @brief - */ - maxval = zdata_statistics_get_max(&m_data_statistics); - minval = zdata_statistics_get_min(&m_data_statistics); - - int32_t nowvaloff = data - minval; - nowval100 = 10 + nowvaloff * 80 / (maxval - minval); - displaydata = nowval100; - - } else { - maxval = 0; - minval = 0; - displaydata = 50; - } +void ecg_algo_process_data(int32_t indata) { + float data = indata; + data = HPFilterExt_update(&m_hp_filter, data); // 高通 + data = LPFilterExt_update(&m_lp_filter, data); // 低通 + data = NOTCHFilterExt_update(&m_notch_filter, data); // 带阻 + // display_data_processer(data); + report_data_processer(data); } void ecg_algo_reset() { - LPFilterExt_init(&m_lp_filter, 40, SAMPLE_PERIOD_S, 5, true); - HPFilterExt_init(&m_hp_filter, 1, SAMPLE_PERIOD_S, 1, true); - NOTCHFilterExt_init(&m_notch_filter, 125, 2, SAMPLE_PERIOD_S, 2, true); - SmoothingFilter_Init(&m_smoothingFilter, 8, true); + LPFilterExt_init(&m_lp_filter, 40, SAMPLE_PERIOD_S, 4, true); + HPFilterExt_init(&m_hp_filter, 1, SAMPLE_PERIOD_S, 1, false); + NOTCHFilterExt_init(&m_notch_filter, 125, 2, SAMPLE_PERIOD_S, 2, false); + + SmoothingFilter_Init(&m_report_smoothing_filter, 8, true); + + { + memset(&report_medianFilter, 0, sizeof(report_medianFilter)); + memset(report_medianBuffer, 0, sizeof(report_medianBuffer)); + report_medianFilter.numNodes = REPORT_MEDIAN_WINDOWS_SIZE; + report_medianFilter.medianBuffer = report_medianBuffer; + MEDIANFILTER_Init(&report_medianFilter); + } + + median_filter_init(&m_display_median_filter, 5); zdata_statistics_clear(&m_data_statistics); zsimple_qrs_clear(); @@ -71,12 +105,6 @@ void ecg_algo_reset() { int32_t ecg_algo_get_report_data() { return reportdata; } int32_t ecg_algo_get_display_data() { return displaydata; } int32_t ecg_algo_get_heart_rate() { return zsimple_qrs_get_heartrate(&m_qrs); } -int32_t ecg_algo_get_raw_data() { return data_af_high_pass_filter; } int32_t ecg_algo_get_max_data() { return maxval; } int32_t ecg_algo_get_min_data() { return minval; } int32_t ecg_algo_get_peak2peak() { return maxval - minval; } - -int32_t ecg_algo_get_data_after_high_pass_filter_data() { return data_af_high_pass_filter; } - -// maxval -// minval \ No newline at end of file diff --git a/app/src/app_service/ecg_service/ecg_algo.h b/app/src/app_service/ecg_service/ecg_algo.h index b126661..24ecca8 100644 --- a/app/src/app_service/ecg_service/ecg_algo.h +++ b/app/src/app_service/ecg_service/ecg_algo.h @@ -5,7 +5,7 @@ #include "app_basic_service/zapp.h" #define SAMPLE_PERIOD_S ((float)(1.0 / SAMPLE_RATE)) -#define STATISTICS_BUF_SIZE ((int32_t)(1.5 / SAMPLE_PERIOD_S)) +#define STATISTICS_BUF_SIZE ((int32_t)(1 / SAMPLE_PERIOD_S)) void ecg_algo_init(); void ecg_algo_process_data(int32_t indata); @@ -15,7 +15,6 @@ int32_t ecg_algo_get_report_data(); int32_t ecg_algo_get_display_data(); int32_t ecg_algo_get_heart_rate(); -int32_t ecg_algo_get_data_after_high_pass_filter_data(); int32_t ecg_algo_get_max_data(); int32_t ecg_algo_get_min_data(); int32_t ecg_algo_get_peak2peak(); \ No newline at end of file diff --git a/app/src/app_service/ecg_service/ecg_service.c b/app/src/app_service/ecg_service/ecg_service.c index 1370b7b..4798d3d 100644 --- a/app/src/app_service/ecg_service/ecg_service.c +++ b/app/src/app_service/ecg_service/ecg_service.c @@ -244,7 +244,7 @@ static void leadoff_state_process(ads129x_capture_data_t* capture_data) { cnt = 0; } - if (cnt > 10) { + if (cnt > 3) { cnt = 0; m_leadoffstate = false; } @@ -268,7 +268,7 @@ static void ads1291_ready_pin_irq(nrfx_gpiote_pin_t pin, nrf_gpiote_polarity_t a // 上报 one_frame_t frame; if (m_report_data_in_raw_mode_flag) { - frame.data = capture_data.ch1data * 0.5; + frame.data = capture_data.ch1data; } else { frame.data = ecg_algo_get_report_data(); } @@ -434,6 +434,7 @@ void ecg_service_set_in_test_mode(bool testmode) { void ecg_service_set_report_data_in_raw_mode(bool rawmode) { m_report_data_in_raw_mode_flag = rawmode; } bool ecg_leadoff_detect() { + return false; if (m_testmode_flag) { return false; } diff --git a/app/src/aproject_config/config.h b/app/src/aproject_config/config.h index 3e7bea5..91ede93 100644 --- a/app/src/aproject_config/config.h +++ b/app/src/aproject_config/config.h @@ -7,7 +7,7 @@ #define CATEGORY "M1003" // 单导联 #define MANUFACTURER_NAME "iflytop" -#define FIRMWARE_VERSION (506) +#define FIRMWARE_VERSION (507) #define BLESTACK_VERSION 1 #define BOOTLOADER_VERSION 1 #define HARDWARE_VERSION (2) @@ -32,7 +32,7 @@ #define MAX_STORAGE_SIZE (MAX_STORAGE_TIMEOUT_S * 400) // 存储最大限制为 (256-8)kbyte #define MAX_FILE_NUM 10 // 最多存储条目数 #define SAMPLE_RATE 500 // 采样率 -#define SAMPLE_PRECISION 16 // 采样精度 +#define SAMPLE_PRECISION 32 // 采样精度 #define AUTOMATIC_SLEEP_TIME 30000 // 开机后自动休眠时间 #define SAMPLE_MIN_TIME_S (30.0) // 采样最小时间 #define LITTLE_DATA_BLOCK_FRAME_NUM 50 // 每次多少帧上报一次 @@ -60,4 +60,3 @@ ***********************************************************************************************************************/ #define ADS1291_SPI_INSTANCE 2 - diff --git a/app/src/ble_data_processer_utils.c b/app/src/ble_data_processer_utils.c index 1361a1a..5d40f89 100644 --- a/app/src/ble_data_processer_utils.c +++ b/app/src/ble_data_processer_utils.c @@ -45,14 +45,8 @@ void report_ecg_data(app_event_t* data) { reportpacket->sample_data_num = ECG_DATA_REPORT_FRAME_NUM; for (int i = 0; i < ECG_DATA_REPORT_FRAME_NUM; i++) { - int32_t frame = data->val.ecg_data_report_event.ecgData->frame[i].data; - if (frame >= INT16_MAX) { - reportpacket->frame[i] = INT16_MAX; - } else if (frame <= INT16_MIN) { - reportpacket->frame[i] = INT16_MIN; - } else { - reportpacket->frame[i] = frame; - } + int32_t frame = data->val.ecg_data_report_event.ecgData->frame[i].data; + reportpacket->frame[i] = frame; } uint8_t leadoffstate = 0; for (int i = 0; i < ECG_DATA_REPORT_FRAME_NUM; i++) { diff --git a/ify_hrs_protocol b/ify_hrs_protocol index c5a8617..286f8bb 160000 --- a/ify_hrs_protocol +++ b/ify_hrs_protocol @@ -1 +1 @@ -Subproject commit c5a86174de8f27f4c8ab6c18dd131effa7b43ea5 +Subproject commit 286f8bbe136357e72b72f7d3e2a5b0ad64e99710 diff --git a/test/MedianFilter.c b/test/MedianFilter.c new file mode 100644 index 0000000..a785900 --- /dev/null +++ b/test/MedianFilter.c @@ -0,0 +1,88 @@ +/* + * MedianFilter.c + * + * Created on: May 19, 2018 + * Author: alexandru.bogdan + */ + +#include "MedianFilter.h" + +int MEDIANFILTER_Init(sMedianFilter_t *medianFilter) +{ + if(medianFilter && medianFilter->medianBuffer && + (medianFilter->numNodes % 2) && (medianFilter->numNodes > 1)) + { + //initialize buffer nodes + for(unsigned int i = 0; i < medianFilter->numNodes; i++) + { + medianFilter->medianBuffer[i].value = 0; + medianFilter->medianBuffer[i].nextAge = &medianFilter->medianBuffer[(i + 1) % medianFilter->numNodes]; + medianFilter->medianBuffer[i].nextValue = &medianFilter->medianBuffer[(i + 1) % medianFilter->numNodes]; + medianFilter->medianBuffer[(i + 1) % medianFilter->numNodes].prevValue = &medianFilter->medianBuffer[i]; + } + //initialize heads + medianFilter->ageHead = medianFilter->medianBuffer; + medianFilter->valueHead = medianFilter->medianBuffer; + medianFilter->medianHead = &medianFilter->medianBuffer[medianFilter->numNodes / 2]; + + return 0; + } + + return -1; +} + +int MEDIANFILTER_Insert(sMedianFilter_t *medianFilter, int sample) +{ + unsigned int i; + sMedianNode_t *newNode, *it; + + if(medianFilter->ageHead == medianFilter->valueHead) + { //if oldest node is also the smallest node, increment value head + medianFilter->valueHead = medianFilter->valueHead->nextValue; + } + + if((medianFilter->ageHead == medianFilter->medianHead) || + (medianFilter->ageHead->value > medianFilter->medianHead->value)) + { //prepare for median correction + medianFilter->medianHead = medianFilter->medianHead->prevValue; + } + + //replace age head with new sample + newNode = medianFilter->ageHead; + newNode->value = sample; + + //remove age head from list + medianFilter->ageHead->nextValue->prevValue = medianFilter->ageHead->prevValue; + medianFilter->ageHead->prevValue->nextValue = medianFilter->ageHead->nextValue; + //increment age head + medianFilter->ageHead = medianFilter->ageHead->nextAge; + + //find new node position + it = medianFilter->valueHead; //set iterator as value head + for(i = 0; i < medianFilter->numNodes - 1; i++) + { + if(sample < it->value) + { + if(i == 0) + { //replace value head if new node is the smallest + medianFilter->valueHead = newNode; + } + break; + } + it = it->nextValue; + } + + //insert new node in list + it->prevValue->nextValue = newNode; + newNode->prevValue = it->prevValue; + it->prevValue = newNode; + newNode->nextValue = it; + + //adjust median node + if(i >= (medianFilter->numNodes / 2)) + { + medianFilter->medianHead = medianFilter->medianHead->nextValue; + } + + return medianFilter->medianHead->value; +} \ No newline at end of file diff --git a/test/MedianFilter.h b/test/MedianFilter.h new file mode 100644 index 0000000..e616ca5 --- /dev/null +++ b/test/MedianFilter.h @@ -0,0 +1,38 @@ +/* + * MedianFilter.h + * + * Created on: May 19, 2018 + * Author: alexandru.bogdan + */ + +#ifndef MEDIANFILTER_H_ +#define MEDIANFILTER_H_ + +#ifdef __cplusplus +extern "C" { +#endif + +typedef struct sMedianNode +{ + int value; //sample value + struct sMedianNode *nextAge; //pointer to next oldest value + struct sMedianNode *nextValue; //pointer to next smallest value + struct sMedianNode *prevValue; //pointer to previous smallest value +}sMedianNode_t; + +typedef struct +{ + unsigned int numNodes; //median node buffer length + sMedianNode_t *medianBuffer; //median node buffer + sMedianNode_t *ageHead; //pointer to oldest value + sMedianNode_t *valueHead; //pointer to smallest value + sMedianNode_t *medianHead; //pointer to median value +}sMedianFilter_t; + +int MEDIANFILTER_Init(sMedianFilter_t *medianFilter); +int MEDIANFILTER_Insert(sMedianFilter_t *medianFilter, int sample); + +#ifdef __cplusplus +} +#endif +#endif \ No newline at end of file diff --git a/test/README.md b/test/README.md new file mode 100644 index 0000000..ddbc98c --- /dev/null +++ b/test/README.md @@ -0,0 +1 @@ +pip install wfdb \ No newline at end of file diff --git a/test/test.c b/test/test.c new file mode 100644 index 0000000..b4b75e4 --- /dev/null +++ b/test/test.c @@ -0,0 +1,27 @@ +#include +#include +#include +#include + +#include "MedianFilter.h" +#define NUM_ELEMENTS 201 + +static sMedianFilter_t medianFilter; +static sMedianNode_t medianBuffer[NUM_ELEMENTS]; + +int main() { + medianFilter.numNodes = NUM_ELEMENTS; + medianFilter.medianBuffer = medianBuffer; + + MEDIANFILTER_Init(&medianFilter); + + while (1) { + int newValue = rand() % 100; + int medianValue = MEDIANFILTER_Insert(&medianFilter, newValue); + printf("New value: %d \tMedian value: %d\r\n", newValue, medianValue); + + Sleep(1); + } + + return 0; +} \ No newline at end of file diff --git a/test/test.py b/test/test.py new file mode 100644 index 0000000..66a9795 --- /dev/null +++ b/test/test.py @@ -0,0 +1,129 @@ +import csv +from datetime import datetime + +import wfdb +import matplotlib.pyplot as plt +import numpy as np +from scipy.signal import medfilt +import pywt + +fliter = int(0.8*500) +Initial_intercept_point = 0 +Final_intercept_point = 2000 +Give_up_size = int(fliter/2) + +# Get dates, high, and low temperatures from file. +filename = 'sunyan.csv' +with open(filename) as f: + reader = csv.reader(f) + header_row = next(reader) + + leaddatas= [] + for row in reader: + try: + leaddata = int(row[0]) + except ValueError: + print(leaddata, 'missing data') + else: + leaddatas.append(leaddata) +#璇诲彇鏂囦欢缁撴潫 + +#bandpass filter +from scipy import signal +fs = 500 # 閲囨牱鐜囦负500 Hz + +# 璁捐浣庨氭护娉㈠櫒 +cutoff_freq = 40 # 鎴棰戠巼 +order = 1 # 婊ゆ尝鍣ㄩ樁鏁 +nyquist = 0.5 * fs +cutoff = cutoff_freq / nyquist +b, a = signal.butter(order, cutoff, btype='low') + +# 搴旂敤浣庨氭护娉㈠櫒 +lowpass_ecg_leaddatas = signal.filtfilt(b, a, leaddatas) + +# 璁捐楂橀氭护娉㈠櫒 +cutoff_freq_high = 0.5 # 鎴棰戠巼 +order_high = 1 # 婊ゆ尝鍣ㄩ樁鏁 +#nyquist = 0.5 * fs +cutoff_high = cutoff_freq_high / nyquist +b, a = signal.butter(order_high, cutoff_high, btype='high') + +# 搴旂敤楂橀氭护娉㈠櫒 +filtered_ecg_signal_leaddatas = signal.filtfilt(b, a, lowpass_ecg_leaddatas) + +''' +# 浣跨敤灏忔尝鍙樻崲,鏁堟灉涓嶅ソ +wavelet = 'sym5' # 閫夋嫨灏忔尝绫诲瀷锛宻ym5鏄竴绉嶅绉板皬娉紝甯哥敤浜嶦CG淇″彿 +levels = 4 # 鍒嗚В绾ф暟 +coeffs = pywt.wavedec(filtered_ecg_signal_leaddatas, wavelet, level=levels) + +# 瀵瑰皬娉㈢郴鏁拌繘琛岄槇鍊煎鐞嗭紝杩欐湁鍔╀簬鍘诲櫔 +threshold = 0.2 * np.max(coeffs[-1]) +coeffs = [pywt.threshold(c, threshold, mode='soft') for c in coeffs] +reconstructed_signal = pywt.waverec(coeffs, wavelet) +''' + +''' +# 灏忔尝鍒嗚В +wavelet = 'db4' +coeffs = pywt.wavedec(filtered_ecg_signal_leaddatas, wavelet, level=4) + +# 瀵归珮棰戠郴鏁拌繘琛岄槇鍊煎鐞嗭紙杞槇鍊硷級 +sigma = np.median(np.abs(coeffs[-1])) / 0.6745 # 浼拌鍣0鏍囧噯宸 +threshold = sigma * np.sqrt(2 * np.log(len(filtered_ecg_signal_leaddatas))) + +# 瀵瑰悇灞傜郴鏁拌繘琛岄槇鍊煎鐞 +denoised_coeffs = [pywt.threshold(c, threshold, mode='soft') for c in coeffs] + +# 閲嶅缓淇″彿 +denoised_signal = pywt.waverec(denoised_coeffs, wavelet) + +''' + +# Define the window size for moving average (adjust as needed) +window_size = 2 + +# Perform moving average filtering +moving_avg_ecg1 = np.convolve(lowpass_ecg_leaddatas, np.ones(window_size)/window_size, mode='valid') + + +#record = wfdb.rdrecord('mit-bih-arrhythmia-database-1.0.0/100', sampfrom=0, sampto=25000, physical=True, channels=[0, ]) +#Original_ECG = record.p_signal[Initial_intercept_point:Final_intercept_point].flatten() +Original_ECG = moving_avg_ecg1 +ECG_baseline = medfilt(Original_ECG, fliter+1) +Totality_Bias = np.sum(ECG_baseline[Give_up_size:-Give_up_size])/(Final_intercept_point-Initial_intercept_point-fliter) +Filtered_ECG = Original_ECG - ECG_baseline +Final_Filtered_ECG = Filtered_ECG[Give_up_size:-Give_up_size]-Totality_Bias + + + +plt.figure(figsize=(100, 8)) +plt.subplot(2, 1, 1) +plt.ylabel("Original ECG signal") +plt.plot(leaddatas)#杈撳嚭鍘熷浘鍍弆eaddatas + +plt.subplot(2, 1, 2) +plt.ylabel("Filtered ECG signal") +plt.plot(Original_ECG)#鍩虹嚎婊ゆ尝缁撴灉 + +plt.show() + + +#蹇冪巼璁$畻 +import biosppy.signals.ecg as ecg + +rpeaks = ecg.hamilton_segmenter(signal=Final_Filtered_ECG, sampling_rate=500) + +def calculate_heart_rate(r_peaks, sampling_rate): + # 璁$畻鐩搁偦R娉箣闂寸殑鏃堕棿闂撮殧锛堝崟浣嶏細绉掞級 + rr_intervals = np.diff(r_peaks) / sampling_rate + # 灏嗘椂闂撮棿闅旇浆鎹负姣忓垎閽熺殑蹇冭烦鏁帮紙bpm锛 + heart_rates = 60 / rr_intervals + return heart_rates + + +# 璁$畻蹇冪巼 +heart_rates = calculate_heart_rate(rpeaks, sampling_rate=500) +print("Heart rates (bpm):", heart_rates) +