diff --git a/MATLAB/app_wrapper/app_io.c b/MATLAB/app_wrapper/app_io.c index 3f2a22a..0264762 100644 --- a/MATLAB/app_wrapper/app_io.c +++ b/MATLAB/app_wrapper/app_io.c @@ -175,8 +175,8 @@ void app_readInputs(const real_T* Buffer) { MB_INTERNAL.param.nominal.U = ReadInputArray(2, 2) * 10; MB_INTERNAL.param.nominal.I = ReadInputArray(2, 3) * 10; - MB_INTERNAL.param.angle.PID_Kp = ReadInputArray(2, 4) * 10000; - MB_INTERNAL.param.angle.PID_Ki = ReadInputArray(2, 5) * 10000; + MB_INTERNAL.param.angle.PID_Kp = ReadInputArray(2, 4) * 65535; + MB_INTERNAL.param.angle.PID_Ki = ReadInputArray(2, 5) * 65535; /*MB_INTERNAL.param.angle.PID_Kd*/dbg_err_limit = ReadInputArray(2, 6); MB_INTERNAL.param.angle.Angle_Max = ReadInputArray(2, 7)/180 * 65535; diff --git a/MATLAB/upp_init.m b/MATLAB/upp_init.m index 6cd838e..b7ce115 100644 --- a/MATLAB/upp_init.m +++ b/MATLAB/upp_init.m @@ -9,4 +9,6 @@ Inom = 25; Fnom = 50; Temperature1 = 2.22; % 20 градусов -Temperature2 = 2.99; % 34 градусов \ No newline at end of file +Temperature2 = 2.99; % 34 градусов + + diff --git a/MATLAB/upp_r2023.slx b/MATLAB/upp_r2023.slx index 5a4f314..dbcc52e 100644 Binary files a/MATLAB/upp_r2023.slx and b/MATLAB/upp_r2023.slx differ diff --git a/UPP/AllLibs/MyLibs b/UPP/AllLibs/MyLibs index 1f7384c..a6b27da 160000 --- a/UPP/AllLibs/MyLibs +++ b/UPP/AllLibs/MyLibs @@ -1 +1 @@ -Subproject commit 1f7384c5edc30469217ca49c6489bdab13460320 +Subproject commit a6b27da4ce1d3cf60a251225beff3e6d0af1193e diff --git a/UPP/Core/PowerMonitor/adc_tools.c b/UPP/Core/PowerMonitor/adc_tools.c index 5d16446..d96f1bb 100644 --- a/UPP/Core/PowerMonitor/adc_tools.c +++ b/UPP/Core/PowerMonitor/adc_tools.c @@ -12,6 +12,7 @@ static void ADC_EnableAllFilters(ADC_Periodic_t *adc) for(int i = 0; i < ADC_NUMB_OF_CHANNELS; i++) { Filter_Start(&adc->filter[i]); +// FilterExpInt_Process(&adc->filter[i], adc->Coefs[i].lZero); } } static void ADC_InitAllFilters(ADC_Periodic_t *adc) @@ -198,8 +199,9 @@ HAL_StatusTypeDef ADC_Handle(ADC_Periodic_t *adc) } - if(Filter_isDataReady(&adc->filter[0])) - adc->f.DataReady = 1; +// if(Filter_isDataReady(&adc->filter[0])) +// adc->f.DataReady = 1; + return HAL_OK; } diff --git a/UPP/Core/PowerMonitor/power_monitor.c b/UPP/Core/PowerMonitor/power_monitor.c index e81d8ec..8706b57 100644 --- a/UPP/Core/PowerMonitor/power_monitor.c +++ b/UPP/Core/PowerMonitor/power_monitor.c @@ -67,6 +67,16 @@ HAL_StatusTypeDef PowerMonitor_Init(PowerMonitor_t *hpm) Filter_Start(&hpm->avg[ADC_TEMP_CHANNELS_START+i]); } + + /* Инициализация фильтра для сглаживания синусоиды*/ + for(int i = 0; i < 2; i++) + { + if(FilterBandPassDerivative_Init(&hpm->ufltr[i], (50.0f*PM_FAST_PERIOD_US/1000000), 0.1)) + return HAL_ERROR; + + Filter_Start(&hpm->ufltr[i]); + } + return HAL_OK; } @@ -115,6 +125,7 @@ void PowerMonitor_SlowCalc(PowerMonitor_t *hpm) /* Расчет всякого для трех фаз отдельно */ float fmean = 0; // средняя частота по трем фазам + float imean = 0; // средний ток по трем фазам float iphase_mean = 0; // средний ток каждой фазы float uphase_mean = 0; // среднее напряжение каждой фазы // Дополнительно посчитаем значения в реальных Вольтах/Амперах @@ -140,6 +151,7 @@ void PowerMonitor_SlowCalc(PowerMonitor_t *hpm) // meas->final.I[i] = iphase_mean*PI/2/SQRT2; /*PI/2 - получить амплитудное, SQRT2 - получить действующее */ iphase_mean = Filter_Process(&hpm->rms[RMS_I+i], meas->slow.I[i]); meas->final.I[i] = Filter_Process(&hpm->rms_exp[RMS_I+i], iphase_mean); + imean += meas->final.I[i]; /* Реальные единицы измерения (Вольты/Амперы) */ meas->real.I[i] = meas->final.I[i]*i_base; @@ -160,7 +172,8 @@ void PowerMonitor_SlowCalc(PowerMonitor_t *hpm) /* Расчет амплитуд трехфазной сети */ float uamp = vector_abs_linear_calc(meas->slow.U[U_AB], meas->slow.U[U_CA])/SQRT2; /* SQRT2 - получить действующее */ - float iamp = vector_abs_phase_calc(meas->slow.I[I_A], meas->slow.I[I_C]); +// float iamp = vector_abs_phase_calc(meas->slow.I[I_A], meas->slow.I[I_C]); + float iamp = imean / 3; meas->final.Uamp = Filter_Process(&hpm->rms_exp[RMS_EXP_U], uamp); meas->final.Iamp = Filter_Process(&hpm->rms_exp[RMS_EXP_I], iamp); @@ -178,16 +191,19 @@ void PowerMonitor_SlowCalc(PowerMonitor_t *hpm) void PowerMonitor_FastCalc(PowerMonitor_t *hpm) { if(hpm == NULL) - return; - /* Считываем АЦП с пересчетами и медианой фильтрацией от выбросов */ - ADC_Handle(&hpm->adc); - - /* Заполняем Напряжения/Токи в о.е. */ + return; float u_base = u2f(PARAM_INTERNAL.nominal.U, 10); float i_base = u2f(PARAM_INTERNAL.nominal.I, 10); PowerMonitor_Measured_t *meas = &hpm->measured; - meas->fast.U[U_AB] = hpm->adc.Data[ADC_CHANNEL_UBA]/u_base; - meas->fast.U[U_CA] = hpm->adc.Data[ADC_CHANNEL_UAC]/u_base; + + /* Считываем АЦП с пересчетами и медианой фильтрацией от выбросов */ + ADC_Handle(&hpm->adc); + + /* Заполняем Напряжения/Токи в о.е. */ + float uba_fast = hpm->adc.Data[ADC_CHANNEL_UBA]/u_base; + float uac_fast = hpm->adc.Data[ADC_CHANNEL_UAC]/u_base; + meas->fast.U[U_AB] = Filter_Process(&hpm->ufltr[U_AB], uba_fast); + meas->fast.U[U_CA] = Filter_Process(&hpm->ufltr[U_CA], uac_fast); meas->fast.U[U_BC] = U_BC_calc(meas->fast.U[U_AB], meas->fast.U[U_CA]); meas->fast.I[I_C] = hpm->adc.Data[ADC_CHANNEL_IC]/i_base; diff --git a/UPP/Core/PowerMonitor/power_monitor.h b/UPP/Core/PowerMonitor/power_monitor.h index d00907c..a9fa87d 100644 --- a/UPP/Core/PowerMonitor/power_monitor.h +++ b/UPP/Core/PowerMonitor/power_monitor.h @@ -105,6 +105,7 @@ typedef struct PowerMonitor_Measured_t measured; ///< Измеренные/рассчитанные величины + FilterBandPassDerivative_t ufltr[2]; ///< Фильтры для сглаживаний напряжений в синусы FilterRMS_t rms[RMS_ALL]; ///< Фильтры для расчета действующего значения Напряжения/Токов FilterExp_t rms_exp[RMS_EXP_ALL]; ///< Фильтры для сглаживания действующего значения Напряжения/Токов +2 для результируюзих U, I FilterAverage_t avg[AVG_ALL]; ///< Фильтры для сглаживания медленных величин АЦП diff --git a/UPP/Core/UPP/angle_control.c b/UPP/Core/UPP/angle_control.c index d9b1427..78a65a2 100644 --- a/UPP/Core/UPP/angle_control.c +++ b/UPP/Core/UPP/angle_control.c @@ -70,7 +70,7 @@ HAL_StatusTypeDef Angle_PID_Init(Angle_Handle_t *hangle, float kp, float ki, flo * @param setpoint Уставка куда регулировать * @param measurement Измеренные регулируемые величины */ -void Angle_PID(Angle_Handle_t *hangle, float setpoint, float measurement) +void Angle_PID(Angle_Handle_t *hangle, float setpoint, float measurement, float Correction) { if(assert_upp(hangle)) return; @@ -81,15 +81,15 @@ void Angle_PID(Angle_Handle_t *hangle, float setpoint, float measurement) /* Ошибка регулирования = уставка - измеренное */ float err = hangle->Iref - hangle->Imeas; /* Ограничение скорости изменения */ - extern float dbg_err_limit; - if(err > dbg_err_limit) - { - err = dbg_err_limit; - } - else if (err < -dbg_err_limit) - { - err = -dbg_err_limit; - } +// extern float dbg_err_limit; +// if(err > dbg_err_limit) +// { +// err = dbg_err_limit; +// } +// else if (err < -dbg_err_limit) +// { +// err = -dbg_err_limit; +// } /* ПИД регулирование */ float open_level = arm_pid_f32(&hangle->pid, err); // 0 - открыть максимально поздно, 1 - открыть макситмально рано @@ -114,7 +114,7 @@ void Angle_PID(Angle_Handle_t *hangle, float setpoint, float measurement) float alpha = alpha_rad/PI*hangle->Config.AngleMax; // угол открытия тиристора в о.е. от максимально заданного /* Выставляем заданный уровень открытия */ - Angle_SetAlpha(hangle, alpha); + Angle_SetAlpha(hangle, alpha, Correction); } /** @@ -136,7 +136,7 @@ void Angle_PID_Reset(Angle_Handle_t *hangle) Filter_Start(&hangle->refFilter); Filter_Process(&hangle->refFilter, 0); - Angle_SetAlpha(hangle, 1); // максимально закрываем + Angle_SetAlpha(hangle, 1, 30); // максимально закрываем Angle_Reset(hangle, UPP_PHASE_A); Angle_Reset(hangle, UPP_PHASE_B); @@ -147,12 +147,13 @@ void Angle_PID_Reset(Angle_Handle_t *hangle) /** * @brief Выставление степени открытия тиристоров. * @param hangle Указатель на таймер -* @param Alpha Угол открытия тиристора в о.е. от 180 градусов: + * @param Alpha Угол открытия тиристора в о.е. от 180 градусов: - 0 - максимально закрыт, - 1 - максимально открыт + * @param Коррекция угла в градусах * @return HAL Status. */ -HAL_StatusTypeDef Angle_SetAlpha(Angle_Handle_t *hangle, float Alpha) +HAL_StatusTypeDef Angle_SetAlpha(Angle_Handle_t *hangle, float Alpha, float Correction) { if(assert_upp(hangle)) return HAL_ERROR; @@ -165,7 +166,7 @@ HAL_StatusTypeDef Angle_SetAlpha(Angle_Handle_t *hangle, float Alpha) Alpha = hangle->Config.AngleMin; // сколько надо выжидать исходя из заданного угла - hangle->alpha_real = Alpha + (30.0/180.0); // 30 градусов - сдвиг между линейными и фазными напряжениями + hangle->alpha_real = Alpha + (Correction/180.0); hangle->alpha = Alpha; diff --git a/UPP/Core/UPP/angle_control.h b/UPP/Core/UPP/angle_control.h index 89bca44..cfad508 100644 --- a/UPP/Core/UPP/angle_control.h +++ b/UPP/Core/UPP/angle_control.h @@ -53,9 +53,9 @@ HAL_StatusTypeDef Angle_PID_Init(Angle_Handle_t *hangle, float kp, float ki, flo // ====== УПРАВЛЕНИЕ ========== /* Управление углом через ПИД регулятор */ -void Angle_PID(Angle_Handle_t *hangle, float setpoint, float measurement); +void Angle_PID(Angle_Handle_t *hangle, float setpoint, float measurement, float Correction); /* Выставление текущего угла открытия тиристоров. */ -HAL_StatusTypeDef Angle_SetAlpha(Angle_Handle_t *hangle, float Angle); +HAL_StatusTypeDef Angle_SetAlpha(Angle_Handle_t *hangle, float Angle, float Correction); /* Установка угла открытия в таймер. */ HAL_StatusTypeDef Angle_Start(Angle_Handle_t *hangle, UPP_Phase_t Phase, float PeriodMs); diff --git a/UPP/Core/UPP/upp_main.c b/UPP/Core/UPP/upp_main.c index b3ab9a8..d0337b9 100644 --- a/UPP/Core/UPP/upp_main.c +++ b/UPP/Core/UPP/upp_main.c @@ -182,9 +182,14 @@ int UPP_While(void) // если пришла команда на остановку if (!upp.call->go) upp.workmode = UPP_Init; - + + // Коррекция для отсчета угла открытия + // 30 градусов - сдвиг между линейными и фазными напряжениями + // 30 градусов - фазовое смщеение эксп. фильтра АЦП для сглаживания напряжений + float Correction = 30 + 0; + // Регулирование тиристоров - Angle_PID(&upp.hangle, u2f(PARAM_PUI.Iref,100), upp.pm.measured.final.Iamp); + Angle_PID(&upp.hangle, u2f(PARAM_PUI.Iref,100), upp.pm.measured.final.Iamp, Correction); // если слишком долгий запуск if((local_time() - upp.StartTick) > (upp.PUI.params->Tdelay*1000)) diff --git a/UPP/Core/UPP/upp_params.c b/UPP/Core/UPP/upp_params.c index 677b8c2..9f36b8e 100644 --- a/UPP/Core/UPP/upp_params.c +++ b/UPP/Core/UPP/upp_params.c @@ -141,15 +141,15 @@ void UPP_Params_ControlInternal(void) { alpha_update = 1; } - if(__CheckParamF(&angle_pid_kp, PARAM_INTERNAL.angle.PID_Kp, 10000)) + if(__CheckParamF(&angle_pid_kp, PARAM_INTERNAL.angle.PID_Kp, 65535)) { alpha_update = 1; } - if(__CheckParamF(&angle_pid_ki, PARAM_INTERNAL.angle.PID_Ki, 10000)) + if(__CheckParamF(&angle_pid_ki, PARAM_INTERNAL.angle.PID_Ki, 65535)) { alpha_update = 1; } - if(__CheckParamF(&angle_pid_kd, PARAM_INTERNAL.angle.PID_Kd, 10000)) + if(__CheckParamF(&angle_pid_kd, PARAM_INTERNAL.angle.PID_Kd, 65535)) { alpha_update = 1; } @@ -279,9 +279,9 @@ HAL_StatusTypeDef UPP_Params_Init(void) /*====== ИНИЦИАЛИЗАЦИЯ МОДУЛЯ angle_control ======*/ // Инициализация ПИД if(Angle_PID_Init(&upp.hangle, - u2f(PARAM_INTERNAL.angle.PID_Kp, 10000), - u2f(PARAM_INTERNAL.angle.PID_Ki, 10000), - u2f(PARAM_INTERNAL.angle.PID_Kd, 10000), + u2f(PARAM_INTERNAL.angle.PID_Kp, 65535), + u2f(PARAM_INTERNAL.angle.PID_Ki, 65535), + u2f(PARAM_INTERNAL.angle.PID_Kd, 65535), PUI_Tnt_CalcAlpha(u2f(PARAM_PUI.Tnt, 1000), PM_SLOW_PERIOD_US)) != HAL_OK) return HAL_ERROR; @@ -456,9 +456,9 @@ void UPP_Params_SetDefault(int pui_default, int internal_default) PARAM_INTERNAL.zc.DebouneCouner = ZERO_CROSS_DEBOUNCE_CNT_DEFAULT; - PARAM_INTERNAL.angle.PID_Kp = ANGLE_PID_KP_COEF_DEFAULT*10000; - PARAM_INTERNAL.angle.PID_Ki = ANGLE_PID_KI_COEF_DEFAULT*10000; - PARAM_INTERNAL.angle.PID_Kd = ANGLE_PID_KD_COEF_DEFAULT*10000; + PARAM_INTERNAL.angle.PID_Kp = ANGLE_PID_KP_COEF_DEFAULT*65535; + PARAM_INTERNAL.angle.PID_Ki = ANGLE_PID_KI_COEF_DEFAULT*65535; + PARAM_INTERNAL.angle.PID_Kd = ANGLE_PID_KD_COEF_DEFAULT*65535; PARAM_INTERNAL.angle.Angle_Max = ANGLE_MAX_PERCENT_DEFAULT*65535; PARAM_INTERNAL.angle.Angle_Min = ANGLE_MIN_PERCENT_DEFAULT*65535; PARAM_INTERNAL.angle.PulseLengthReserve = ANGLE_PULSE_LENGTH_RESERVE_PERCENT_DEFAULT*100; diff --git a/UPP/MDK-ARM/UPP.uvoptx b/UPP/MDK-ARM/UPP.uvoptx index 0d4b4d0..c3a2b14 100644 --- a/UPP/MDK-ARM/UPP.uvoptx +++ b/UPP/MDK-ARM/UPP.uvoptx @@ -1,4 +1,4 @@ - + 1.0 @@ -45,7 +45,7 @@ 79 66 8 - + 1 @@ -104,16 +104,16 @@ 0 0 6 - - - - - - - - - - + + + + + + + + + + STLink\ST-LINKIII-KEIL_SWO.dll @@ -142,8 +142,8 @@ 0 0 ..\Core\PowerMonitor\power_monitor.c - - + + @@ -175,19 +175,19 @@ 0 0 - - + + 0 0 0 - - - - - - - - + + + + + + + + 1 1 @@ -222,7 +222,7 @@ 79 66 8 - + 1 @@ -281,16 +281,16 @@ 0 0 6 - - - - - - - - - - + + + + + + + + + + STLink\ST-LINKIII-KEIL_SWO.dll @@ -307,7 +307,7 @@ 0 ARMDBGFLAGS - + 0 @@ -325,7 +325,7 @@ -U005600373433510237363934 -O206 -SF10000 -C0 -A0 -I0 -HNlocalhost -HP7184 -P1 -N00("ARM CoreSight SW-DP (ARM Core") -D00(2BA01477) -L00(0) -TO131123 -TC168000000 -TT10000000 -TP21 -TDS800D -TDT1 -TDC1F -TIEFFFFFFFF -TIP8 -FO7 -FD20000000 -FC1000 -FN1 -FF0STM32F4xx_1024.FLM -FS08000000 -FL0100000 -FP0($$Device:STM32F417ZGTx$CMSIS\Flash\STM32F4xx_1024.FLM) -WA0 -WE0 -WVCE4 -WS2710 -WM0 -WP2 - + 0 @@ -524,19 +524,19 @@ 0 0 - - + + 0 0 0 - - - - - - - - + + + + + + + + 0 diff --git a/UPP/MDK-ARM/UPP.uvprojx b/UPP/MDK-ARM/UPP.uvprojx index 0c2212a..3def99e 100644 --- a/UPP/MDK-ARM/UPP.uvprojx +++ b/UPP/MDK-ARM/UPP.uvprojx @@ -1,7 +1,10 @@ - - + + + 2.1 +
### uVision Project, (C) Keil Software
+ UPP @@ -17,28 +20,28 @@ Keil.STM32F4xx_DFP.2.17.1 https://www.keil.com/pack/ IRAM(0x20000000-0x2002FFFF) IRAM2(0x10000000-0x1000FFFF) IROM(0x8000000-0x80FFFFF) CLOCK(25000000) FPU2 CPUTYPE("Cortex-M4") TZ - - - + + + 0 - - - - - - - - - - + + + + + + + + + + $$Device:STM32F427ZGTx$CMSIS\SVD\STM32F427x.svd 0 0 - - - - - + + + + + 0 0 @@ -53,15 +56,15 @@ 1 1 1 - + 1 0 0 0 0 - - + + 0 0 0 @@ -70,8 +73,8 @@ 0 0 - - + + 0 0 0 @@ -80,15 +83,15 @@ 0 1 - - + + 0 0 0 0 1 - + 0 @@ -102,8 +105,8 @@ 0 0 3 - - + + 0 @@ -136,11 +139,11 @@ 1 BIN\UL2V8M.DLL - - - - - + + + + + 0 @@ -173,7 +176,7 @@ 0 0 "Cortex-M4" - + 0 0 0 @@ -308,7 +311,7 @@ 0x10000 - + 1 @@ -335,9 +338,9 @@ 0 0 - + USE_HAL_DRIVER,STM32F427xx, ARM_MATH_CM4 - + ../Core/Inc;../Drivers/STM32F4xx_HAL_Driver/Inc;../Drivers/STM32F4xx_HAL_Driver/Inc/Legacy;../Drivers/CMSIS/Device/ST/STM32F4xx/Include;../Drivers/CMSIS/Include;../AllLibs/ExtMemory/Inc;../AllLibs/Modbus/Inc;../AllLibs/MyLibs/MyLibs/Inc;../AllLibs/MyLibs/RTT;../AllLibs/PeriphGeneral/Inc;../Core/Configs;../Core/PowerMonitor;../Core/Thyristors;../Core/UPP @@ -353,10 +356,10 @@ 0 1 - - - - + + + + @@ -366,15 +369,15 @@ 0 1 0 - - - - - - - - - + + + + + + + + + @@ -949,8 +952,8 @@ 2 2 11 - - + + 1 @@ -966,10 +969,10 @@ 2 0 - - - - + + + + @@ -1004,28 +1007,28 @@ Keil.STM32F4xx_DFP.2.17.1 https://www.keil.com/pack/ IRAM(0x20000000,0x00020000) IRAM2(0x10000000,0x00010000) IROM(0x08000000,0x00100000) CPUTYPE("Cortex-M4") FPU2 CLOCK(12000000) ELITTLE - - + + UL2CM3(-S0 -C0 -P0 -FD20000000 -FC1000 -FN1 -FF0STM32F4xx_1024 -FS08000000 -FL0100000 -FP0($$Device:STM32F417ZGTx$CMSIS\Flash\STM32F4xx_1024.FLM)) 0 $$Device:STM32F417ZGTx$Drivers\CMSIS\Device\ST\STM32F4xx\Include\stm32f4xx.h - - - - - - - - - + + + + + + + + + $$Device:STM32F417ZGTx$CMSIS\SVD\STM32F41x.svd 0 0 - - - - - + + + + + 0 0 @@ -1040,15 +1043,15 @@ 1 1 1 - + 1 0 0 0 0 - - + + 0 0 0 @@ -1057,8 +1060,8 @@ 0 0 - - + + 0 0 0 @@ -1067,15 +1070,15 @@ 0 0 - - + + 0 0 0 0 1 - + 0 @@ -1089,8 +1092,8 @@ 0 0 3 - - + + 0 @@ -1124,10 +1127,10 @@ 1 BIN\UL2CM3.DLL "" () - - - - + + + + 0 @@ -1160,7 +1163,7 @@ 0 0 "Cortex-M4" - + 0 0 0 @@ -1295,7 +1298,7 @@ 0x10000 - + 1 @@ -1322,9 +1325,9 @@ 0 0 - + USE_HAL_DRIVER,STM32F417xx, ARM_MATH_CM4 - + ../Core/Inc;../Drivers/STM32F4xx_HAL_Driver/Inc;../Drivers/STM32F4xx_HAL_Driver/Inc/Legacy;../Drivers/CMSIS/Device/ST/STM32F4xx/Include;../Drivers/CMSIS/Include;../AllLibs/ExtMemory/Inc;../AllLibs/Modbus/Inc;../AllLibs/MyLibs/MyLibs/Inc;../AllLibs/MyLibs/RTT;../AllLibs/PeriphGeneral/Inc;../Core/Configs;../Core/PowerMonitor;../Core/Thyristors;../Core/UPP @@ -1340,10 +1343,10 @@ 0 1 - - - - + + + + @@ -1353,15 +1356,15 @@ 0 1 0 - - - - - - - - - + + + + + + + + + @@ -1516,8 +1519,8 @@ 2 2 11 - - + + 1 @@ -1546,10 +1549,10 @@ 2 2 - - - - + + + + @@ -1577,8 +1580,8 @@ 2 2 11 - - + + 1 @@ -1607,10 +1610,10 @@ 2 2 - - - - + + + + @@ -1638,8 +1641,8 @@ 2 2 11 - - + + 1 @@ -1668,10 +1671,10 @@ 2 2 - - - - + + + + @@ -1699,8 +1702,8 @@ 2 2 11 - - + + 1 @@ -1729,10 +1732,10 @@ 2 2 - - - - + + + + @@ -1880,8 +1883,8 @@ 2 2 11 - - + + 1 @@ -1910,10 +1913,10 @@ 2 2 - - - - + + + + @@ -2196,8 +2199,8 @@ 2 2 11 - - + + 1 @@ -2213,10 +2216,10 @@ 2 0 - - - - + + + + @@ -2233,32 +2236,34 @@ + - + - + - + - + - - + + - + - - + + - + + @@ -2267,5 +2272,5 @@ -
+
diff --git a/Информация для программиста (УПП СП СЭД)/CALC/calc_biquad.m b/Информация для программиста (УПП СП СЭД)/CALC/calc_biquad.m index ff9aeb3..8d38f1d 100644 --- a/Информация для программиста (УПП СП СЭД)/CALC/calc_biquad.m +++ b/Информация для программиста (УПП СП СЭД)/CALC/calc_biquad.m @@ -2,10 +2,11 @@ clear all; close all; clc; %% Параметры моделирования -Fs = 100000; % Частота дискретизации [Гц] +Fs = 1/25e-6; % Частота дискретизации [Гц] T = 0.5; % Время моделирования [с] t = 0:1/Fs:T-1/Fs; % Временной вектор N = length(t); % Количество отсчетов +Fsrez = 50; %% Уровни шума для разных каналов noise_levels.voltage = 0.2; % 2% шума для напряжений @@ -17,28 +18,28 @@ fprintf('=== АВТОМАТИЧЕСКИЙ РАСЧЕТ КОЭФФИЦИЕНТО % 1. Полосовой фильтр 45-55 Гц для напряжений % [b_bpf, a_bpf, coeffs_bpf] = BiquadFilterDesigner.bpf(20, 10, Fs); -[b_bpf, a_bpf, coeffs_bpf] = BiquadFilterDesigner.lpf(100, Fs); +[b_bpf, a_bpf, coeffs_bpf] = BiquadFilterDesigner.lpf(Fsrez, Fs); fprintf('1. Полосовой фильтр 45-55 Гц:\n'); BiquadFilterDesigner.generate_c_code(coeffs_bpf, 'voltage_bpf'); % 2. ФНЧ 100 Гц для токов -[b_lpf_current, a_lpf_current, coeffs_lpf_current] = BiquadFilterDesigner.lpf(100, Fs); +[b_lpf_current, a_lpf_current, coeffs_lpf_current] = BiquadFilterDesigner.lpf(Fsrez, Fs); fprintf('2. ФНЧ 100 Гц (токи):\n'); BiquadFilterDesigner.generate_c_code(coeffs_lpf_current, 'current_lpf'); - -% 3. ФНЧ 10 Гц для температур -[b_lpf_temp, a_lpf_temp, coeffs_lpf_temp] = BiquadFilterDesigner.lpf(10, Fs); -fprintf('3. ФНЧ 10 Гц (температуры):\n'); -BiquadFilterDesigner.generate_c_code(coeffs_lpf_temp, 'temperature_lpf'); +% +% % 3. ФНЧ 10 Гц для температур +% [b_lpf_temp, a_lpf_temp, coeffs_lpf_temp] = BiquadFilterDesigner.lpf(10, Fs); +% fprintf('3. ФНЧ 10 Гц (температуры):\n'); +% BiquadFilterDesigner.generate_c_code(coeffs_lpf_temp, 'temperature_lpf'); % Вывод коэффициентов в консоль fprintf('\n=== РАСЧЕТНЫЕ КОЭФФИЦИЕНТЫ ===\n'); fprintf('Напряжение (BPF 45-55 Гц): b = [%.6f, %.6f, %.6f], a = [1, %.6f, %.6f]\n', ... b_bpf, a_bpf(2), a_bpf(3)); -fprintf('Ток (LPF 100 Гц): b = [%.6f, %.6f, %.6f], a = [1, %.6f, %.6f]\n', ... - b_lpf_current, a_lpf_current(2), a_lpf_current(3)); -fprintf('Температура (LPF 10 Гц): b = [%.6f, %.6f, %.6f], a = [1, %.6f, %.6f]\n\n', ... - b_lpf_temp, a_lpf_temp(2), a_lpf_temp(3)); +% fprintf('Ток (LPF 100 Гц): b = [%.6f, %.6f, %.6f], a = [1, %.6f, %.6f]\n', ... +% b_lpf_current, a_lpf_current(2), a_lpf_current(3)); +% fprintf('Температура (LPF 10 Гц): b = [%.6f, %.6f, %.6f], a = [1, %.6f, %.6f]\n\n', ... +% b_lpf_temp, a_lpf_temp(2), a_lpf_temp(3)); %% Генерация тестовых сигналов @@ -53,25 +54,25 @@ f_current = 50; current_clean = 1 * sin(2*pi*f_current*t); %.* ... %(1 + 0.2 * sin(2*pi*2*t)); % амплитудная модуляция 2 Гц -% 3. Температура (медленно меняющийся сигнал) -temperature_clean = 25 + 2 * sin(2*pi*0.1*t) + ... % медленные колебания 0.1 Гц - 0.5 * sin(2*pi*1*t); % быстрые колебания 1 Гц +% % 3. Температура (медленно меняющийся сигнал) +% temperature_clean = 25 + 2 * sin(2*pi*0.1*t) + ... % медленные колебания 0.1 Гц +% 0.5 * sin(2*pi*1*t); % быстрые колебания 1 Гц %% Добавление шума voltage_noisy = voltage_clean + noise_levels.voltage * randn(size(t)); current_noisy = current_clean + noise_levels.current * randn(size(t)); -temperature_noisy = temperature_clean + noise_levels.temperature * randn(size(t)); +% temperature_noisy = temperature_clean + noise_levels.temperature * randn(size(t)); %% Фильтрация сигналов voltage_filtered = filter(b_bpf, a_bpf, voltage_noisy); current_filtered = filter(b_lpf_current, a_lpf_current, current_noisy); -temperature_filtered = filter(b_lpf_temp, a_lpf_temp, temperature_noisy); +% temperature_filtered = filter(b_lpf_temp, a_lpf_temp, temperature_noisy); %% НОРМАЛИЗАЦИЯ УСИЛЕНИЯ (важно для правильного SNR) % Получаем АЧХ для нормализации [h_bpf, f_bpf] = freqz(b_bpf, a_bpf, 1024, Fs); [h_lpf_curr, f_lpf_curr] = freqz(b_lpf_current, a_lpf_current, 1024, Fs); -[h_lpf_temp, f_lpf_temp] = freqz(b_lpf_temp, a_lpf_temp, 1024, Fs); +% [h_lpf_temp, f_lpf_temp] = freqz(b_lpf_temp, a_lpf_temp, 1024, Fs); % Нормализация полосового фильтра на центральной частоте [~, idx_50hz] = min(abs(f_bpf - 50)); @@ -85,11 +86,11 @@ gain_lpf_current = sum(b_lpf_current) / (1 + sum(a_lpf_current(2:end))); if gain_lpf_current > 0 current_filtered = current_filtered / gain_lpf_current; end - -gain_lpf_temp = sum(b_lpf_temp) / (1 + sum(a_lpf_temp(2:end))); -if gain_lpf_temp > 0 - temperature_filtered = temperature_filtered / gain_lpf_temp; -end +% +% gain_lpf_temp = sum(b_lpf_temp) / (1 + sum(a_lpf_temp(2:end))); +% if gain_lpf_temp > 0 +% temperature_filtered = temperature_filtered / gain_lpf_temp; +% end %% ОКНО 1: НАПРЯЖЕНИЕ figure('Name', 'Анализ напряжения'); @@ -213,68 +214,68 @@ text(0.1, 0.4, sprintf('Улучшение: %.1f дБ', improvement_current), 'F text(0.1, 0.2, sprintf('Задержка 50 Гц: %.1f мс', gd_lpf_curr(idx_50hz_curr)/Fs*1000), 'FontSize', 12); title('Статистика фильтрации тока'); axis off; - -%% ОКНО 3: ТЕМПЕРАТУРА -figure('Name', 'Анализ температуры'); - -% Временные характеристики -subplot(2,3,1); -plot(t, temperature_noisy, 'b', 'LineWidth', 1); hold on; -plot(t, temperature_filtered, 'r', 'LineWidth', 2); -plot(t, temperature_clean, 'g--', 'LineWidth', 1); -title('Температура: временная область'); -legend('С шумом', 'Фильтрованный', 'Идеальный', 'Location', 'best'); -xlabel('Время [с]'); ylabel('Температура [°C]'); -grid on; - -% АЧХ фильтра -subplot(2,3,2); -plot(f_lpf_temp, 20*log10(abs(h_lpf_temp)), 'LineWidth', 2); -title('АЧХ: ФНЧ 10 Гц'); -xlabel('Частота [Гц]'); ylabel('Усиление [дБ]'); -grid on; xlim([0, 20]); - -% Спектр фильтрованного сигнала -subplot(2,3,3); -[P_temp, f_temp] = pwelch(temperature_filtered, [], [], 1024, Fs); -plot(f_temp, 10*log10(P_temp), 'LineWidth', 2); -title('Спектр фильтрованной температуры'); -xlabel('Частота [Гц]'); ylabel('Мощность [дБ]'); -grid on; xlim([0, 10]); - -% Групповая задержка -subplot(2,3,4); -[gd_lpf_temp, f_gd_temp] = grpdelay(b_lpf_temp, a_lpf_temp, 1024, Fs); -plot(f_gd_temp, gd_lpf_temp/Fs*1000, 'LineWidth', 2); -title('Групповая задержка фильтра'); -xlabel('Частота [Гц]'); ylabel('Задержка [мс]'); -grid on; xlim([0, 20]); - -% Детальный вид (последне 0.1 секунды) -idx_end_temp = max(1, length(t) - 0.1*Fs + 1):length(t); % последние 100 мс -subplot(2,3,5); -plot(t(idx_end_temp), temperature_noisy(idx_end_temp), 'b', 'LineWidth', 2); hold on; -plot(t(idx_end_temp), temperature_filtered(idx_end_temp), 'r', 'LineWidth', 2); -plot(t(idx_end_temp), temperature_clean(idx_end_temp), 'g--', 'LineWidth', 2); -title('Температура: УВЕЛИЧЕННЫЙ ВИД (900-1000 мс)'); -legend('С шумом', 'Фильтрованный', 'Идеальный', 'Location', 'best'); -xlabel('Время [с]'); ylabel('Температура [°C]'); -grid on; -xlim([t(idx_end_temp(1)) t(idx_end_temp(end))]); - -% Статистика -subplot(2,3,6); -snr_temp_in = snr(temperature_clean, temperature_noisy - temperature_clean); -snr_temp_out = snr(temperature_clean, temperature_filtered - temperature_clean); -improvement_temp = snr_temp_out - snr_temp_in; - -[~, idx_1hz] = min(abs(f_gd_temp - 1)); -text(0.1, 0.8, sprintf('SNR вход: %.1f дБ', snr_temp_in), 'FontSize', 12); -text(0.1, 0.6, sprintf('SNR выход: %.1f дБ', snr_temp_out), 'FontSize', 12); -text(0.1, 0.4, sprintf('Улучшение: %.1f дБ', improvement_temp), 'FontSize', 12); -text(0.1, 0.2, sprintf('Задержка 1 Гц: %.1f мс', gd_lpf_temp(idx_1hz)/Fs*1000), 'FontSize', 12); -title('Статистика фильтрации температуры'); -axis off; +% +% %% ОКНО 3: ТЕМПЕРАТУРА +% figure('Name', 'Анализ температуры'); +% +% % Временные характеристики +% subplot(2,3,1); +% plot(t, temperature_noisy, 'b', 'LineWidth', 1); hold on; +% plot(t, temperature_filtered, 'r', 'LineWidth', 2); +% plot(t, temperature_clean, 'g--', 'LineWidth', 1); +% title('Температура: временная область'); +% legend('С шумом', 'Фильтрованный', 'Идеальный', 'Location', 'best'); +% xlabel('Время [с]'); ylabel('Температура [°C]'); +% grid on; +% +% % АЧХ фильтра +% subplot(2,3,2); +% plot(f_lpf_temp, 20*log10(abs(h_lpf_temp)), 'LineWidth', 2); +% title('АЧХ: ФНЧ 10 Гц'); +% xlabel('Частота [Гц]'); ylabel('Усиление [дБ]'); +% grid on; xlim([0, 20]); +% +% % Спектр фильтрованного сигнала +% subplot(2,3,3); +% [P_temp, f_temp] = pwelch(temperature_filtered, [], [], 1024, Fs); +% plot(f_temp, 10*log10(P_temp), 'LineWidth', 2); +% title('Спектр фильтрованной температуры'); +% xlabel('Частота [Гц]'); ylabel('Мощность [дБ]'); +% grid on; xlim([0, 10]); +% +% % Групповая задержка +% subplot(2,3,4); +% [gd_lpf_temp, f_gd_temp] = grpdelay(b_lpf_temp, a_lpf_temp, 1024, Fs); +% plot(f_gd_temp, gd_lpf_temp/Fs*1000, 'LineWidth', 2); +% title('Групповая задержка фильтра'); +% xlabel('Частота [Гц]'); ylabel('Задержка [мс]'); +% grid on; xlim([0, 20]); +% +% % Детальный вид (последне 0.1 секунды) +% idx_end_temp = max(1, length(t) - 0.1*Fs + 1):length(t); % последние 100 мс +% subplot(2,3,5); +% plot(t(idx_end_temp), temperature_noisy(idx_end_temp), 'b', 'LineWidth', 2); hold on; +% plot(t(idx_end_temp), temperature_filtered(idx_end_temp), 'r', 'LineWidth', 2); +% plot(t(idx_end_temp), temperature_clean(idx_end_temp), 'g--', 'LineWidth', 2); +% title('Температура: УВЕЛИЧЕННЫЙ ВИД (900-1000 мс)'); +% legend('С шумом', 'Фильтрованный', 'Идеальный', 'Location', 'best'); +% xlabel('Время [с]'); ylabel('Температура [°C]'); +% grid on; +% xlim([t(idx_end_temp(1)) t(idx_end_temp(end))]); +% +% % Статистика +% subplot(2,3,6); +% snr_temp_in = snr(temperature_clean, temperature_noisy - temperature_clean); +% snr_temp_out = snr(temperature_clean, temperature_filtered - temperature_clean); +% improvement_temp = snr_temp_out - snr_temp_in; +% +% [~, idx_1hz] = min(abs(f_gd_temp - 1)); +% text(0.1, 0.8, sprintf('SNR вход: %.1f дБ', snr_temp_in), 'FontSize', 12); +% text(0.1, 0.6, sprintf('SNR выход: %.1f дБ', snr_temp_out), 'FontSize', 12); +% text(0.1, 0.4, sprintf('Улучшение: %.1f дБ', improvement_temp), 'FontSize', 12); +% text(0.1, 0.2, sprintf('Задержка 1 Гц: %.1f мс', gd_lpf_temp(idx_1hz)/Fs*1000), 'FontSize', 12); +% title('Статистика фильтрации температуры'); +% axis off; %% Вывод результатов в командное окно fprintf('\n=== ИТОГИ ФИЛЬТРАЦИИ С АВТОРАСЧЕТОМ КОЭФФИЦИЕНТОВ ===\n\n'); @@ -288,7 +289,7 @@ fprintf(' SNR вход: %.1f дБ, выход: %.1f дБ, улучшение: % snr_current_in, snr_current_out, improvement_current); fprintf(' Задержка на 50 Гц: %.1f мс\n\n', gd_lpf_curr(idx_50hz_curr)/Fs*1000); -fprintf('ТЕМПЕРАТУРА (ФНЧ 10 Гц):\n'); -fprintf(' SNR вход: %.1f дБ, выход: %.1f дБ, улучшение: %.1f дБ\n', ... - snr_temp_in, snr_temp_out, improvement_temp); -fprintf(' Задержка на 1 Гц: %.1f мс\n\n', gd_lpf_temp(idx_1hz)/Fs*1000); +% fprintf('ТЕМПЕРАТУРА (ФНЧ 10 Гц):\n'); +% fprintf(' SNR вход: %.1f дБ, выход: %.1f дБ, улучшение: %.1f дБ\n', ... +% snr_temp_in, snr_temp_out, improvement_temp); +% fprintf(' Задержка на 1 Гц: %.1f мс\n\n', gd_lpf_temp(idx_1hz)/Fs*1000); diff --git a/Информация для программиста (УПП СП СЭД)/CALC/calc_filter.m b/Информация для программиста (УПП СП СЭД)/CALC/calc_filter.m new file mode 100644 index 0000000..abedc10 --- /dev/null +++ b/Информация для программиста (УПП СП СЭД)/CALC/calc_filter.m @@ -0,0 +1,230 @@ +%% ===== РАСЧЕТ КОЭФФИЦИЕНТОВ БИКВАДРАТНОГО ФИЛЬТРА 50 Гц ===== +% Запуск: Ctrl+S (сохранить), потом F5 (запустить) + +clear all; close all; clc; + +%% 1. ПАРАМЕТРЫ (МЕНЯТЬ ЗДЕСЬ) +fs = 1/25e-6; % Частота дискретизации (Гц) +fc = 100; % Частота среза 50 Гц +Q = 0.707; % Добротность (0.707 = Баттерворт) +filter_type = 'lpf'; % 'lpf', 'hpf', 'bpf', 'notch' + +test_freq = 55; + +%% 2. РАСЧЕТ КОЭФФИЦИЕНТОВ +w0 = 2 * pi * fc / fs; +alpha = sin(w0) / (2 * Q); +cos_w0 = cos(w0); + +switch filter_type + case 'lpf' % ФНЧ + b0 = (1 - cos_w0) / 2; + b1 = 1 - cos_w0; + b2 = b0; + a0 = 1 + alpha; + a1 = -2 * cos_w0; + a2 = 1 - alpha; + + case 'hpf' % ФВЧ + b0 = (1 + cos_w0) / 2; + b1 = -(1 + cos_w0); + b2 = b0; + a0 = 1 + alpha; + a1 = -2 * cos_w0; + a2 = 1 - alpha; + + case 'bpf' % Полосовой + b0 = alpha; + b1 = 0; + b2 = -alpha; + a0 = 1 + alpha; + a1 = -2 * cos_w0; + a2 = 1 - alpha; + + case 'notch' % Режекторный (50 Гц) + b0 = 1; + b1 = -2 * cos_w0; + b2 = 1; + a0 = 1 + alpha; + a1 = -2 * cos_w0; + a2 = 1 - alpha; +end + +% Нормализация +b = [b0, b1, b2] / a0; +a = [1, a1/a0, a2/a0]; + +% Проверка полюсов +poles = roots([1, a(2), a(3)]); +if any(abs(poles) >= 1) + error('Фильтр НЕУСТОЙЧИВ! Полюса: %s', mat2str(poles)); +end +%% 3. ВЫВОД РЕЗУЛЬТАТОВ +fprintf('\n══════════════════════════════════════════════════\n'); +fprintf('ПАРАМЕТРЫ ФИЛЬТРА:\n'); +fprintf('Тип: %s\n', upper(filter_type)); +fprintf('Частота среза: %.1f Гц\n', fc); +fprintf('Частота дискретизации: %.0f Гц\n', fs); +fprintf('Добротность Q: %.3f\n', Q); +fprintf('══════════════════════════════════════════════════\n'); + +fprintf('\nКОЭФФИЦИЕНТЫ ДЛЯ CMSIS-DSP:\n'); +fprintf('float32_t coeffs[5] = {\n'); +fprintf(' %ff, // b0\n', b(1)); +fprintf(' %ff, // b1\n', b(2)); +fprintf(' %ff, // b2\n', b(3)); +fprintf(' %ff, // a1\n', a(2)); +fprintf(' %ff // a2\n', a(3)); +fprintf('};\n'); + +fprintf('\nФОРМАТ ДЛЯ Biquad_InitDirect:\n'); +fprintf('Biquad_InitDirect(&filter, %.6ff, %.6ff, %.6ff, %.6ff, %.6ff);\n', ... + b(1), b(2), b(3), a(2), a(3)); + +%% 4. ПРОВЕРКА ЧАСТОТНОЙ ХАРАКТЕРИСТИКИ +[h, f] = freqz(b, a, 2048, fs); + +figure(); + +% АЧХ +subplot(1, 2, 1); +plot(f, 20*log10(abs(h)), 'LineWidth', 2); +grid on; hold on; +xline(fc, '--r', 'Частота среза', 'LineWidth', 1.5); +xlim([0, fs/2]); +ylim([-60, 5]); +xlabel('Частота (Гц)'); +ylabel('Усиление (дБ)'); +title(['АЧХ: ' upper(filter_type) ' фильтр ' num2str(fc) ' Гц']); + +% ФЧХ +subplot(1, 2, 2); +plot(f, angle(h)*180/pi, 'LineWidth', 2); +grid on; hold on; +xline(fc, '--r', 'Частота среза', 'LineWidth', 1.5); +xlim([0, fs/2]); +xlabel('Частота (Гц)'); +ylabel('Фазовый сдвиг (градусы)'); +title(['ФЧХ: ' upper(filter_type) ' фильтр ' num2str(fc) ' Гц']); + +%% 4b. ГРАФИК ЗАДЕРЖКИ (40-60 Гц) +figure(); + +% Рассчитываем на нужном диапазоне частот +f_range = 40:0.1:60; % Частоты с шагом 0.1 Гц +h_range = freqz(b, a, f_range, fs); + +% Групповая задержка (численная производная) +phase_unwrapped = unwrap(angle(h_range)); +w_range = 2*pi*f_range/fs; +group_delay = -gradient(phase_unwrapped) ./ gradient(w_range); + +% Фазовая задержка +phase_delay = -angle(h_range) ./ w_range; + +% График 1: Групповая задержка +subplot(1, 2, 1); +plot(f_range, group_delay * 1000, 'LineWidth', 3); +grid on; hold on; +xline(test_freq, '--r', '50 Гц', 'LineWidth', 2, 'FontSize', 12); +xlim([40, 60]); +ylim([0, max(group_delay*1000)*1.1]); +xlabel('Частота (Гц)', 'FontSize', 12); +ylabel('Групповая задержка (мс)', 'FontSize', 12); +title(['Групповая задержка: ' upper(filter_type) ' фильтр'], 'FontSize', 14); + +% Значение на 50 Гц +idx_50 = find(f_range >= test_freq, 1); +if ~isempty(idx_50) + delay_50hz = group_delay(idx_50) * 1000; + plot(test_freq, delay_50hz, 'ro', 'MarkerSize', 15, 'LineWidth', 3); + text(test_freq + 0.5, delay_50hz*1.05, sprintf('%.2f мс', delay_50hz), ... + 'FontSize', 14, 'FontWeight', 'bold', 'BackgroundColor', 'white'); +end + +% График 2: Фазовая задержка +subplot(1, 2, 2); +plot(f_range, phase_delay * 1000, 'LineWidth', 3, 'Color', [0, 0.5, 0]); +grid on; hold on; +xline(test_freq, '--r', '50 Гц', 'LineWidth', 2, 'FontSize', 12); +xlim([40, 60]); +xlabel('Частота (Гц)', 'FontSize', 12); +ylabel('Фазовая задержка (мс)', 'FontSize', 12); +title(['Фазовая задержка: ' upper(filter_type) ' фильтр'], 'FontSize', 14); + +% Значение на 50 Гц +if ~isempty(idx_50) + phase_delay_50hz = phase_delay(idx_50) * 1000; + plot(test_freq, phase_delay_50hz, 'ro', 'MarkerSize', 15, 'LineWidth', 3); + text(test_freq + 0.5, phase_delay_50hz*1.05, sprintf('%.2f мс', phase_delay_50hz), ... + 'FontSize', 14, 'FontWeight', 'bold', 'BackgroundColor', 'white'); +end + +%% 4c. ГРАФИК УСИЛЕНИЯ (40-60 Гц) +figure(); + +% Усиление в дБ +gain_db = 20*log10(abs(h_range)); + +plot(f_range, gain_db, 'LineWidth', 3, 'Color', [0.8, 0, 0]); +grid on; hold on; +xline(test_freq, '--r', '50 Гц', 'LineWidth', 2, 'FontSize', 12); +xlim([40, 60]); +ylim([min(gain_db)-1, max(gain_db)+1]); +xlabel('Частота (Гц)', 'FontSize', 12); +ylabel('Усиление (дБ)', 'FontSize', 12); +title(['Усиление: ' upper(filter_type) ' фильтр (40-60 Гц)'], 'FontSize', 14); + +% Значение на 50 Гц +if ~isempty(idx_50) + gain_50hz = gain_db(idx_50); + plot(test_freq, gain_50hz, 'ro', 'MarkerSize', 15, 'LineWidth', 3); + text(test_freq +0.5, gain_50hz, sprintf('%.2f дБ', gain_50hz), ... + 'FontSize', 14, 'FontWeight', 'bold', 'BackgroundColor', 'white'); +end + +%% 5. ТЕСТОВЫЙ СИГНАЛ +t = 0:1/fs:0.2; % 200 мс +f_test = test_freq; % Тест на 50 Гц +signal = sin(2*pi*f_test*t); + +% Фильтрация +filtered = filter(b, a, signal); + +figure(); +subplot(2, 1, 1); +plot(t, signal, 'b', 'LineWidth', 1.5); +hold on; grid on; +plot(t, filtered, 'r', 'LineWidth', 2); +xlabel('Время (с)'); +ylabel('Амплитуда'); +legend('Исходный', 'После фильтра', 'Location', 'best'); +title(['Тест фильтра: ' num2str(f_test) ' Гц']); + +% Спектр +subplot(2, 1, 2); +[P1, f1] = pwelch(signal, [], [], [], fs); +[P2, f2] = pwelch(filtered, [], [], [], fs); +plot(f1, 10*log10(P1), 'b', 'LineWidth', 1.5); +hold on; grid on; +plot(f2, 10*log10(P2), 'r', 'LineWidth', 2); +xlim([0, 200]); +xlabel('Частота (Гц)'); +ylabel('Мощность (дБ)'); +legend('Исходный', 'После фильтра', 'Location', 'best'); +title('Спектр сигнала'); + + +fprintf('\n══════════════════════════════════════════════════\n'); +fprintf('ЗНАЧЕНИЯ НА %d Гц:\n', f_test); + +if ~isempty(idx_50) + fprintf('Групповая задержка: %.2f мс\n', delay_50hz); + fprintf('Фазовая задержка: %.2f мс\n', phase_delay_50hz); + fprintf('Усиление: %.2f дБ (%.3f в линейном масштабе)\n', ... + gain_50hz, abs(h_range(idx_50))); +else + fprintf('Не удалось рассчитать значения на 50 Гц\n'); +end + +fprintf('══════════════════════════════════════════════════\n'); \ No newline at end of file diff --git a/Информация для программиста (УПП СП СЭД)/CALC/test_filt.m b/Информация для программиста (УПП СП СЭД)/CALC/test_filt.m new file mode 100644 index 0000000..67a2b51 --- /dev/null +++ b/Информация для программиста (УПП СП СЭД)/CALC/test_filt.m @@ -0,0 +1,199 @@ +%% МОДЕЛИРОВАНИЕ ПОЛОСОВОГО ФИЛЬТРА С ДИФФЕРЕНЦИАТОРОМ +clear all; close all; clc; + +%% 1. ПАРАМЕТРЫ ФИЛЬТРА (подставь свои значения) +b0 = 0.000392540911; +b1 = 0.0; +b2 = -0.000392540911; +a1 = -1.99915338; +a2 = 0.999214947; + +% Или рассчитай новые: +center_freq = 50; % Центральная частота (Гц) +sample_freq = 40000; % Частота дискретизации (Гц) +bandwidth = 5; % Ширина полосы (Гц) + +% %% 2. РАСЧЕТ КОЭФФИЦИЕНТОВ (если нужно) +% if 1 % Поставь 0 если используешь готовые коэффициенты выше +% % Отношение частот +% fc_ratio = center_freq / sample_freq; +% bandwidth_ratio = bandwidth / center_freq; +% +% % Расчет коэффициентов полосового фильтра +% w0 = 2 * pi * fc_ratio; +% Q = 1 / bandwidth_ratio; +% alpha = sin(w0) / (2 * Q); +% cos_w0 = cos(w0); +% +% % Коэффициенты биквадратного полосового фильтра +% b0_bp = alpha; +% b1_bp = 0; +% b2_bp = -alpha; +% a0_bp = 1 + alpha; +% a1_bp = -2 * cos_w0; +% a2_bp = 1 - alpha; +% +% % Нормализация (a0 = 1) +% b0 = b0_bp / a0_bp; +% b1 = b1_bp / a0_bp; +% b2 = b2_bp / a0_bp; +% a1 = a1_bp / a0_bp; +% a2 = a2_bp / a0_bp; +% end + +fprintf('Коэффициенты фильтра:\n'); +fprintf('b0 = %.6f\n', b0); +fprintf('b1 = %.6f\n', b1); +fprintf('b2 = %.6f\n', b2); +fprintf('a1 = %.6f\n', a1); +fprintf('a2 = %.6f\n', a2); + +%% 3. ПРОВЕРКА УСТОЙЧИВОСТИ +poles = roots([1, a1, a2]); +fprintf('\nПолюса фильтра:\n'); +for i = 1:length(poles) + fprintf(' pole%d = %.6f %+.6fj (|pole| = %.6f)\n', ... + i, real(poles(i)), imag(poles(i)), abs(poles(i))); +end + +if any(abs(poles) >= 1) + fprintf('⚠️ ВНИМАНИЕ: Фильтр НЕУСТОЙЧИВ!\n'); +else + fprintf('✅ Фильтр устойчив\n'); +end + +%% 5. ЧАСТОТНАЯ ХАРАКТЕРИСТИКА ПОЛНОГО ФИЛЬТРА (с дифференциатором) +% Дифференциатор: H_diff(z) = 1 - z^-1 +b_diff = [1, -1]; +a_diff = 1; + +% % Каскадное соединение: дифференциатор + полосовой фильтр +% b_total = conv(b_diff, b_bp); % Перемножение полиномов +% a_total = conv(a_diff, a_bp); +% полосовой фильтр +b_total = [b0, b1, b2]; +a_total = [1, a1, a2]; + +[h_total, f_total] = freqz(b_total, a_total, 2048, sample_freq); + +figure(); +subplot(1,2,1); +plot(f_total, 20*log10(abs(h_total)), 'r', 'LineWidth', 2); +grid on; hold on; +xline(center_freq, '--r', sprintf('%d Гц', center_freq), 'LineWidth', 1.5); +xlim([0, sample_freq/2]); +ylim([-60, 20]); +xlabel('Частота (Гц)'); +ylabel('Усиление (дБ)'); +title('АЧХ полного фильтра (дифференциатор + полосовой)'); + +subplot(1,2,2); +plot(f_total, angle(h_total)*180/pi, 'r', 'LineWidth', 2); +grid on; hold on; +xline(center_freq, '--r', sprintf('%d Гц', center_freq), 'LineWidth', 1.5); +xlim([0, sample_freq/2]); +xlabel('Частота (Гц)'); +ylabel('Фаза (градусы)'); +title('ФЧХ полного фильтра (дифференциатор + полосовой)'); + +%% 6. МОДЕЛИРОВАНИЕ ВО ВРЕМЕННОЙ ОБЛАСТИ +% Создаем тестовый сигнал +duration = 0.2; % 200 мс +t = 0:1/sample_freq:duration; + +% Сигнал 50 Гц + шум +signal_freq = 50; +signal = sin(2*pi*signal_freq*t) + 0.1*randn(size(t)); + +% Имитация работы твоего кода на C +% Прямая форма II с дифференциатором +x1 = 0; x2 = 0; % Состояния входа фильтра +y1 = 0; y2 = 0; % Состояния выхода фильтра +prev_input = 0; % Для дифференциатора + +filtered_signal = zeros(size(signal)); + +for i = 1:length(signal) + % 1. Дифференциатор + diff = signal(i); + prev_input = signal(i); + + % 2. Полосовой фильтр (прямая форма II) + y = b0*diff + b1*x1 + b2*x2 - a1*y1 - a2*y2; + + % 3. Обновление состояний + x2 = x1; + x1 = diff; + y2 = y1; + y1 = y; + + filtered_signal(i) = y; +end + +% Встроенная функция MATLAB для сравнения +filtered_matlab = filter(b_total, a_total, signal); + +figure(); + +% Сигнал и выход +subplot(2,1,1); +plot(t, signal, 'b', 'LineWidth', 1); +hold on; grid on; +plot(t, filtered_signal, 'r', 'LineWidth', 2); +plot(t, filtered_matlab, 'g--', 'LineWidth', 1); +xlabel('Время (с)'); +ylabel('Амплитуда'); +legend('Исходный сигнал', 'Наш фильтр (имитация C)', 'MATLAB filter()', ... + 'Location', 'best'); +title(sprintf('Тестовый сигнал: %d Гц + помеха 150 Гц + шум', signal_freq)); + +% Ошибка между нашей реализацией и MATLAB +subplot(2,1,2); +error = filtered_signal - filtered_matlab; +plot(t, error, 'k', 'LineWidth', 1); +grid on; +xlabel('Время (с)'); +ylabel('Ошибка'); +title('Разница между нашей реализацией и MATLAB filter()'); +fprintf('\nМаксимальная ошибка: %.2e\n', max(abs(error))); + +%% 7. АНАЛИЗ НА ЧАСТОТЕ 50 Гц +freq_test = 50; +w_test = 2*pi*freq_test/sample_freq; +z = exp(1i*w_test); + +% Передаточная функция полного фильтра +H_z = (b0 + b1*z^-1 + b2*z^-2) / (1 + a1*z^-1 + a2*z^-2) * (1 - z^-1); + +fprintf('\n══════════════════════════════════════════════════\n'); +fprintf('АНАЛИЗ НА %d Гц:\n', freq_test); +fprintf('Усиление: %.3f (%.2f дБ)\n', abs(H_z), 20*log10(abs(H_z))); +fprintf('Фазовый сдвиг: %.1f градусов\n', angle(H_z)*180/pi); +fprintf('Задержка: %.2f мс (фазовая)\n', -angle(H_z)*1000/(2*pi*freq_test)); +fprintf('══════════════════════════════════════════════════\n'); + +%% 8. ГРАФИК ПОЛЮСОВ И НУЛЕЙ +figure(); + +% Единичная окружность +theta = linspace(0, 2*pi, 100); +plot(cos(theta), sin(theta), 'k--', 'LineWidth', 1); +hold on; grid on; axis equal; + +% Полюса (красные кресты) +plot(real(poles), imag(poles), 'rx', 'MarkerSize', 15, 'LineWidth', 2); + +% Нули полосового фильтра +zeros_bp = roots([b0, b1, b2]); +plot(real(zeros_bp), imag(zeros_bp), 'bo', 'MarkerSize', 10, 'LineWidth', 2); + +% Нули дифференциатора (z=1) +plot(1, 0, 'go', 'MarkerSize', 10, 'LineWidth', 2); + +xlim([-1.2, 1.2]); +ylim([-1.2, 1.2]); +xlabel('Re(z)'); +ylabel('Im(z)'); +title('Диаграмма полюсов и нулей фильтра'); +legend('Единичная окружность', 'Полюса', 'Нули полосового фильтра', ... + 'Нуль дифференциатора (z=1)', 'Location', 'best'); \ No newline at end of file