import numpy as np
from andromeda import andromeda_score

def lp_score(masses, spec_masses, spec_intensities, tol=0.1, p=0.06):
    """
    计算局部化概率（LP）分数。
    
    参数:
    - masses: 理论质谱峰列表
    - spec_masses: 实验质谱峰列表
    - spec_intensities: 实验峰对应的强度
    - tol: 质量容差
    - p: 匹配概率
    
    返回:
    - LP分数
    """
    # 初始化匹配计数
    k = 0
    for mass in masses:
        # 找到最近的实验峰并检查匹配
        closest_idx = np.abs(np.array(spec_masses) - mass).argmin()
        if np.abs(spec_masses[closest_idx] - mass) <= tol:
            k += 1
    
    # 计算Andromeda得分
    n = len(masses)  # 理论峰总数
    base_score = andromeda_score(k, n, p)
    
    # 计算LP归一化
    lp = np.exp((base_score - 100) * np.log(10) / 10)
    lp /= np.sum(lp)  # 归一化操作
    
    return lp

# 示例计算
masses = [100, 105, 110, 115, 120]  # 示例理论峰
spec_masses = [101, 106, 111, 116, 121]  # 示例实验峰
spec_intensities = [10, 20, 30, 40, 50]  # 实验峰强度
lp_score(masses, spec_masses, spec_intensities)