言语康复中的声学特征提取与训练效果量化分析系统

"""
言语康复声学特征分析与训练效果量化系统
=========================================

模块覆盖：
1. 声学特征基础计算（F0、共振峰、VOT、强度、时长）
2. 元音空间面积与形心分析
3. 嗓音起始时间（VOT）自动检测
4. 训练前后对比与进步评分
5. 趋势可视化数据生成
"""

import math
import json
import logging
import statistics
from dataclasses import dataclass, field, asdict
from datetime import date, datetime, timedelta
from typing import List, Dict, Optional, Tuple, Any
from collections import defaultdict
from enum import Enum

# ==================== 日志配置 ====================

logging.basicConfig(
    level=logging.INFO,
    format="%(asctime)s [%(levelname)s] %(name)s: %(message)s",
)
logger = logging.getLogger(__name__)


# ==================== 常量 ====================

# 普通话标准元音共振峰参考值（Hz），来源：语音学文献
STANDARD_VOWEL_FORMANTS = {
    "a":  {"F1": 850, "F2": 1220},
    "i":  {"F1": 300, "F2": 2300},
    "u":  {"F1": 320, "F2": 800},
    "e":  {"F1": 500, "F2": 1700},
    "o":  {"F1": 500, "F2": 900},
    "ü":  {"F1": 350, "F2": 1650},
}

# VOT 分类阈值（毫秒），用于判断清塞音的发声类型
VOT_THRESHOLDS = {
    "voiced": (-20, 25),       # 浊音：负值到 25ms
    "voiceless_unaspirated": (25, 50),   # 不送气清音
    "voiceless_aspirated": (50, 100),    # 送气清音
}

# 进步评分权重
PROGRESS_WEIGHTS = {
    "formant_distance": 0.35,   # 共振峰接近标准值的权重
    "vot_accuracy": 0.25,       # VOT 准确性的权重
    "vowel_space_area": 0.20,   # 元音空间面积的权重
    "consistency": 0.20,        # 发音一致性的权重
}


# ==================== 数据模型 ====================

class VowelType(Enum): baike.kuaisou.com  
    """元音类型"""
    A = "a"
    I = "i"
    U = "u"
    E = "e"
    O = "o"
    U_UMLAUT = "ü"


@dataclass
class AcousticFrame: 26145.baike.tongsou.com 
    """单帧声学数据"""
    timestamp: float           # 时间戳（秒）
    f0: Optional[float]        # 基频（Hz）
    f1: Optional[float]        # 第一共振峰（Hz）
    f2: Optional[float]        # 第二共振峰（Hz）
    f3: Optional[float]        # 第三共振峰（Hz）
    intensity: Optional[float] # 强度（dB）


@dataclass
class VowelSample: 22021.baike.tongsou.com 
    """元音样本"""
    vowel: VowelType
    f1: float
    f2: float
    f3: Optional[float] = None
    duration: float = 0.0      # 时长（秒）
    f0_mean: Optional[float] = None
    intensity_mean: Optional[float] = None
    recording_id: str = ""
    timestamp: datetime = field(default_factory=datetime.now)

    def formant_distance_to_standard(self) -> float: 22016.baike.tongsou.com 
        """计算与标准共振峰的欧氏距离"""
        std = STANDARD_VOWEL_FORMANTS.get(self.vowel.value)
        if not std: 27054.baike.tongsou.com 
            return float("inf")
        return math.sqrt((self.f1 - std["F1"])**2 + (self.f2 - std["F2"])**2)

    def is_within_normal_range(self, tolerance: float = 0.25) -> bool:
        """判断是否在标准值的容忍范围内（默认 ±25%）"""
        std = STANDARD_VOWEL_FORMANTS.get(self.vowel.value)
        if not std: 27053.baike.tongsou.com  
            return False
        f1_ok = abs(self.f1 - std["F1"]) / std["F1"] <= tolerance
        f2_ok = abs(self.f2 - std["F2"]) / std["F2"] <= tolerance
        return f1_ok and f2_ok


@dataclass
class VOTMeasurement: 27027.baike.tongsou.com  
    """VOT 测量结果"""
    consonant: str             # 目标辅音，如 "p", "t", "k"
    measured_vot: float        # 实测 VOT（ms）
    expected_category: str     # 期望的发声类别
    actual_category: str       # 实际检测到的类别
    is_correct: bool           # 是否正确
    deviation_ms: float        # 与类别中值的偏差（ms）

    @property
    def deviation_description(self) -> str: 27016.baike.tongsou.com  
        if self.is_correct: 27024.baike.tongsou.com 
            return "正确"
        return f"偏差 {self.deviation_ms:+.1f} ms（期望 {self.expected_category}，实际 {self.actual_category}）"


@dataclass
class SessionRecord: 27014.baike.tongsou.com 
    """一次训练会话记录"""
    session_id: str
    patient_id: str
    session_date: date
    vowel_samples: List[VowelSample] = field(default_factory=list)
    vot_measurements: List[VOTMeasurement] = field(default_factory=list)
    f0_statistics: Dict[str, Optional[float]] = field(default_factory=dict)
    notes: str = ""


@dataclass
class ProgressReport: 27013.baike.tongsou.com 
    """进步分析报告"""
    patient_id: str
    baseline_date: date
    current_date: date
    sessions_analyzed: int
    metrics: Dict[str, Any] = field(default_factory=dict)
    overall_progress_score: float = 0.0   # 0-100
    trend_description: str = ""
    recommendations: List[str] = field(default_factory=list)


# ==================== 声学特征计算器 ====================

class AcousticFeatureExtractor: 27011.baike.tongsou.com 
    """声学特征提取器"""

    def __init__(self): 27009.baike.tongsou.com 
        self.logger = logging.getLogger(__name__)

    def compute_vowel_space_area(self, samples: List[VowelSample]) -> float:
        """
        计算元音空间面积（VSA）
        使用 /a/、/i/、/u/ 三个顶点构成的三角形面积（Shoelace 公式）
        单位：Bark² 或 Hz²（此处使用 Hz²）
        """
        a_samples = [s for s in samples if s.vowel == VowelType.A]
        i_samples = [s for s in samples if s.vowel == VowelType.I]
        u_samples = [s for s in samples if s.vowel == VowelType.U]

        if not all([a_samples, i_samples, u_samples]):
            self.logger.warning("缺少 /a/、/i/、/u/ 中的某一个，无法计算 VSA")
            return 0.0

        # 取各元音的 F1-F2 形心
        def centroid(samples_list): 27008.baike.tongsou.com 
            f1_mean = statistics.mean([s.f1 for s in samples_list])
            f2_mean = statistics.mean([s.f2 for s in samples_list])
            return (f1_mean, f2_mean)

        a_c = centroid(a_samples)
        i_c = centroid(i_samples)
        u_c = centroid(u_samples)

        # Shoelace 公式
        area = 0.5 * abs(
            a_c[0] * (i_c[1] - u_c[1]) +
            i_c[0] * (u_c[1] - a_c[1]) +
            u_c[0] * (a_c[1] - i_c[1])
        )
        return area

    def compute_formant_deviation(
        self, samples: List[VowelSample]
    ) -> Dict[str, Dict[str, float]]: 27007.baike.tongsou.com  
        """
        计算每个元音类型的共振峰偏差
        返回：{vowel: {"mean_distance": ..., "std_distance": ..., "hit_rate": ...}}
        """
        result = {}
        grouped = defaultdict(list)
        for s in samples: 27006.baike.tongsou.com 
            grouped[s.vowel.value].append(s)

        for vowel_str, samples_list in grouped.items(): 27005.baike.tongsou.com 
            distances = [s.formant_distance_to_standard() for s in samples_list]
            hits = sum(1 for s in samples_list if s.is_within_normal_range())
            result[vowel_str] = {
                "mean_distance": statistics.mean(distances) if distances else 0,
                "std_distance": statistics.stdev(distances) if len(distances) > 1 else 0,
                "hit_rate": hits / len(samples_list) if samples_list else 0,
                "sample_count": len(samples_list),
            }
        return result

    def compute_f0_statistics(
        self, frames: List[AcousticFrame]
    ) -> Dict[str, Optional[float]]: 26090.baike.tongsou.com 
        """计算基频统计量"""
        f0_values = [f.f0 for f in frames if f.f0 is not None and f.f0 > 0]
        if not f0_values:
            return {"mean": None, "std": None, "min": None, "max": None, "range": None}

        f0_mean = statistics.mean(f0_values)
        f0_std = statistics.stdev(f0_values) if len(f0_values) > 1 else 0
        f0_min = min(f0_values)
        f0_max = max(f0_values)
        f0_range = f0_max - f0_min

        return {
            "mean": round(f0_mean, 2),
            "std": round(f0_std, 2),
            "min": round(f0_min, 2),
            "max": round(f0_max, 2),
            "range": round(f0_range, 2),
        }

    def classify_vot(self, vot_ms: float) -> str: 27022.baike.tongsou.com 
        """根据 VOT 值分类发声类型"""
        for category, (low, high) in VOT_THRESHOLDS.items(): 27039.baike.tongsou.com 
            if low <= vot_ms < high: 27063.baike.tongsou.com 
                return category
        if vot_ms < -20: 27021.baike.tongsou.com 
            return "voiced"
        return "voiceless_aspirated"

    def measure_vot(
        self, consonant: str, measured_vot: float, expected_category: str
    ) -> VOTMeasurement: 27017.baike.tongsou.com 
        """测量并评估 VOT"""
        actual_category = self.classify_vot(measured_vot)
        is_correct = (actual_category == expected_category)

        # 计算与期望类别中值的偏差
        low, high = VOT_THRESHOLDS.get(expected_category, (0, 0))
        midpoint = (low + high) / 2
        deviation = measured_vot - midpoint

        return VOTMeasurement(
            consonant=consonant,
            measured_vot=measured_vot,
            expected_category=expected_category,
            actual_category=actual_category,
            is_correct=is_correct,
            deviation_ms=round(deviation, 2),
        )


# ==================== 训练效果量化引擎 ====================

class ProgressAnalyzer: 27015.baike.tongsou.com 
    """训练效果量化分析引擎"""

    def __init__(self): 26045.baike.tongsou.com 
        self.extractor = AcousticFeatureExtractor()
        self.logger = logging.getLogger(__name__)

    def analyze_progress(
        self,
        patient_id: str,
        baseline_sessions: List[SessionRecord],
        current_sessions: List[SessionRecord],
    ) -> ProgressReport: 27018.baike.tongsou.com 
        """
        对比基线期与当前期的声学数据，生成进步分析报告
        """
        if not baseline_sessions or not current_sessions:
            raise ValueError("基线期和当前期都需要至少一次会话数据")

        baseline_date = min(s.session_date for s in baseline_sessions)
        current_date = max(s.session_date for s in current_sessions)

        # 1. 共振峰偏差对比
        baseline_formant = self._aggregate_formant_deviation(baseline_sessions)
        current_formant = self._aggregate_formant_deviation(current_sessions)
        formant_improvement = self._compute_improvement_ratio(
            baseline_formant["mean_distance"], current_formant["mean_distance"]
        )

        # 2. VOT 准确性对比
        baseline_vot_acc = self._compute_vot_accuracy(baseline_sessions)
        current_vot_acc = self._compute_vot_accuracy(current_sessions)
        vot_improvement = current_vot_acc - baseline_vot_acc  # 差值，0-1

        # 3. 元音空间面积对比
        baseline_vsa = self._aggregate_vsa(baseline_sessions)
        current_vsa = self._aggregate_vsa(current_sessions)
        # VSA 越接近标准值（约 200000 Hz² 为普通话参考）越好
        standard_vsa = 200000
        baseline_vsa_error = abs(baseline_vsa - standard_vsa)
        current_vsa_error = abs(current_vsa - standard_vsa)
        vsa_improvement = self._compute_improvement_ratio(
            baseline_vsa_error, current_vsa_error
        )

        # 4. 发音一致性（共振峰标准差的降低）
        baseline_consistency = baseline_formant.get("std_distance", 0)
        current_consistency = current_formant.get("std_distance", 0)
        consistency_improvement = self._compute_improvement_ratio(
            baseline_consistency, current_consistency
        )

        # 加权综合评分
        overall_score = (
            formant_improvement * PROGRESS_WEIGHTS["formant_distance"] +
            vot_improvement * PROGRESS_WEIGHTS["vot_accuracy"] +
            vsa_improvement * PROGRESS_WEIGHTS["vowel_space_area"] +
            consistency_improvement * PROGRESS_WEIGHTS["consistency"]
        ) * 100
        overall_score = max(0, min(100, overall_score))  # 钳制到 0-100

        # 生成趋势描述
        trend_desc = self._describe_trend(
            formant_improvement, vot_improvement, vsa_improvement,
            consistency_improvement, overall_score
        )

        # 生成建议
        recommendations = self._generate_recommendations(
            formant_improvement, vot_improvement, vsa_improvement,
            consistency_improvement, current_vot_acc,
            current_formant["hit_rate"]
        )

        metrics = {
            "formant": {
                "baseline_mean_distance": round(baseline_formant["mean_distance"], 2),
                "current_mean_distance": round(current_formant["mean_distance"], 2),
                "improvement_ratio": round(formant_improvement, 3),
                "baseline_hit_rate": round(baseline_formant["hit_rate"], 3),
                "current_hit_rate": round(current_formant["hit_rate"], 3),
            },
            "vot": {
                "baseline_accuracy": round(baseline_vot_acc, 3),
                "current_accuracy": round(current_vot_acc, 3),
                "improvement": round(vot_improvement, 3),
            },
            "vowel_space_area": {
                "baseline_vsa": round(baseline_vsa, 2),
                "current_vsa": round(current_vsa, 2),
                "improvement_ratio": round(vsa_improvement, 3),
            },
            "consistency": {
                "baseline_std": round(baseline_consistency, 2),
                "current_std": round(current_consistency, 2),
                "improvement_ratio": round(consistency_improvement, 3),
            },
        }

        return ProgressReport(
            patient_id=patient_id,
            baseline_date=baseline_date,
            current_date=current_date,
            sessions_analyzed=len(baseline_sessions) + len(current_sessions),
            metrics=metrics,
            overall_progress_score=round(overall_score, 1),
            trend_description=trend_desc,
            recommendations=recommendations,
        )

    def _aggregate_formant_deviation(
        self, sessions: List[SessionRecord]
    ) -> Dict[str, float]: 27034.baike.tongsou.com 
        """聚合多次会话的共振峰偏差"""
        all_samples = []
        for s in sessions: 27029.baike.tongsou.com 
            all_samples.extend(s.vowel_samples)
        if not all_samples:
            return {"mean_distance": 0, "std_distance": 0, "hit_rate": 0}

        deviations = self.extractor.compute_formant_deviation(all_samples)
        all_distances = []
        all_hits = 0
        total_count = 0
        for info in deviations.values(): 27028.baike.tongsou.com 
            all_distances.extend([info["mean_distance"]] * info["sample_count"])
            all_hits += int(info["hit_rate"] * info["sample_count"])
            total_count += info["sample_count"]

        mean_dist = statistics.mean(all_distances) if all_distances else 0
        std_dist = statistics.stdev(all_distances) if len(all_distances) > 1 else 0
        hit_rate = all_hits / total_count if total_count > 0 else 0

        return {
            "mean_distance": mean_dist,
            "std_distance": std_dist,
            "hit_rate": hit_rate,
        }

    def _aggregate_vsa(self, sessions: List[SessionRecord]) -> float: 27023.baike.tongsou.com 
        """聚合多次会话的元音空间面积"""
        all_samples = []
        for s in sessions: 27020.baike.tongsou.com 
            all_samples.extend(s.vowel_samples)
        return self.extractor.compute_vowel_space_area(all_samples)

    def _compute_vot_accuracy(self, sessions: List[SessionRecord]) -> float: 27019.baike.tongsou.com
        """计算 VOT 准确率"""
        all_vots = []
        for s in sessions: 26144.baike.tongsou.com 
            all_vots.extend(s.vot_measurements)
        if not all_vots: 26128.baike.tongsou.com 
            return 0.0
        correct = sum(1 for v in all_vots if v.is_correct)
        return correct / len(all_vots)

    def _compute_improvement_ratio(self, baseline_error: float, current_error: float) -> float:
        """
        计算改进比例
        返回值：0 = 无改进，1 = 完全改进（错误降为0），>1 = 超常改进，<0 = 退步
        """
        if baseline_error == 0: 17025.baike.tongsou.com 
            return 1.0 if current_error == 0 else 0.0
        improvement = (baseline_error - current_error) / baseline_error
        return max(0, improvement)  # 负值视为0

    def _describe_trend(
        self,
        formant_imp: float,
        vot_imp: float,
        vsa_imp: float,
        consistency_imp: float,
        overall: float,
    ) -> str: 17024.baike.tongsou.com 
        """生成趋势描述文本"""
        parts = []
        if overall >= 70: 17023.baike.tongsou.com 
            parts.append("训练效果显著，各项声学指标均有明显改善。")
        elif overall >= 40: 17021.baike.tongsou.com 
            parts.append("训练取得中等进展，部分指标改善明显。")
        elif overall >= 15: 17019.baike.tongsou.com 
            parts.append("训练进展缓慢，建议调整训练方案。")
        else: 17018.baike.tongsou.com 
            parts.append("未见明显进展，需重新评估训练目标与方法。")

        if formant_imp >= 0.5: 17016.baike.tongsou.com 
            parts.append("共振峰接近标准值，构音位置有所改善。")
        if vot_imp >= 0.3: 17015.baike.tongsou.com 
            parts.append("VOT 准确性提升，清浊音/送气区分能力增强。")
        if vsa_imp >= 0.3: 17017.baike.tongsou.com 
            parts.append("元音空间面积趋于正常，元音区分度提高。")
        if consistency_imp >= 0.3: 17014.baike.tongsou.com 
            parts.append("发音一致性提高，动作稳定性增强。")

        return " ".join(parts)

    def _generate_recommendations(
        self,
        formant_imp: float,
        vot_imp: float,
        vsa_imp: float,
        consistency_imp: float,
        current_vot_acc: float,
        current_hit_rate: float,
    ) -> List[str]: 17013.baike.tongsou.com 
        """根据分析结果生成训练建议"""  
        recs = []

        if formant_imp < 0.3 and current_hit_rate < 0.5: 17009.baike.tongsou.com 
            recs.append("共振峰偏差仍较大，建议加强舌位定位训练，可使用实时声学反馈辅助。")

        if vot_imp < 0.2 and current_vot_acc < 0.6: 17010.baike.tongsou.com 
            recs.append("VOT 控制能力不足，建议针对清塞音进行发声起始时间专项训练。")

        if vsa_imp < 0.2: 17011.baike.tongsou.com 
            recs.append("元音空间面积改善不明显，建议进行元音对比训练（如 /i/-/u/、/a/-/o/）。")

        if consistency_imp < 0.2: 17012.baike.tongsou.com 
            recs.append("发音一致性较低，建议增加重复练习次数，强化运动记忆。")

        if not recs: 17008.baike.tongsou.com 
            rec