danalyzer.cs

用C#写的USB数据采集程序

using System;
using System.Collections;

namespace DAS
{
	/// <summary>
	/// 多道数据分析类
	/// </summary>
    public class DAnalyzer
    {
        #region "受保护的属性"
        /// <summary>
        /// 处理后的数据
        /// </summary>
        private double[] _data;
        #endregion

        #region "构造函数"
        /// <summary>
        /// 构造函数
        /// </summary>
        public DAnalyzer()
        {
            _data = new double[0];
        }
        #endregion

        #region "公共方法"

        #region "能谱平滑"
        /// <summary>
        /// 能谱平滑
        /// </summary>
        /// <param name="data">要平滑的能谱数据</param>
        /// <param name="type">平衡方法</param>
        /// <returns>平滑后的能谱数据</returns>
        public void Smooth(ref int[] data, SmoothType type)
        {
            switch (type)
            {
                case SmoothType.points5://5点平滑
                    if (data.Length > 5)
                    {
                        Transform(ref data, -3, 12, 17, 12, -3);
                    }
                    break;
                case SmoothType.points7://7点平滑
                    if (data.Length > 7)
                    {
                        Transform(ref data, -2, 3, 6, 7, 6, 3, -2);
                    }
                    break;
                case SmoothType.points9://9点平滑
                    if (data.Length > 9)
                    {
                        Transform(ref data, -21, 14, 39, 54, 59, 54, 39, 14, -21);
                    }
                    break;
                case SmoothType.points11://11点平滑
                    if (data.Length > 11)
                    {
                        Transform(ref data, -36, 9, 44, 69, 84, 89, 84, 69, 44, 9, -36);
                    }
                    break;
                case SmoothType.points13://13点平滑
                    if (data.Length > 13)
                    {
                        Transform(ref data, -11, 0, 9, 16, 21, 24, 25, 24, 21, 16, 9, 0, -11);
                    }
                    break;
            }
        }
        #endregion

        #region "自定义能谱平滑"
        /// <summary>
        /// 自定义能谱平滑
        /// </summary>
        /// <param name="data">要平滑的能谱数据</param>
        /// <param name="paras">平滑参数</param>
        public void Smooth(ref int[] data, params int[] paras)
        {
            Transform(ref data, paras);
        }
        #endregion

        #region "寻峰"
        /// <summary>
        /// 寻峰
        /// </summary>
        /// <param name="data">要寻峰的能谱数据</param>
        /// <param name="results">结果数据</param>
        /// <param name="type">寻峰方法</param>
        public void Peak(int[] data, out PeakResult[] results, int startIndex, int endIndex, PeakType type)
        {
            if (startIndex >= endIndex)
            {
                throw new Exception("StartIndex is over EndIndex.");
            }
            if (startIndex < 0)
            {
                throw new Exception("StartIndex is overflow.");
            }
            if (endIndex > data.Length - 1)
            {
                throw new Exception("EndIndex is overflow.");
            }
            switch (type)
            {
                case PeakType.derivative://导数法
                    double[] data1 = new double[data.Length];
                    double[] data2 = new double[data.Length];
                    Array.Copy(data, data1, data1.Length);
                    //求二阶导数
                    Transform(ref data1, 42.0, 5, 0, -3, -4, -3, 0, 5);
                    //求三阶导数
                    Array.Copy(data1, data2, data1.Length);
                    Transform(ref data2, 252.0, 22, -67, -58, 0, 58, 67, -22);

                    ArrayList tmpresults = new ArrayList();
                    PeakResult presult;

                    //寻找局部最小值，进行峰判定
                    if ((endIndex <= 4) || startIndex >= data1.Length - 4)
                    {
                        results = new PeakResult[0];
                        return;
                    }
                    if (startIndex < 4)
                    {
                        startIndex = 4;
                    }
                    if (endIndex > data1.Length - 4)
                    {
                        endIndex = data1.Length - 4;
                    }
                    for (int i = startIndex; i < endIndex; i++)
                    {
                        if ((data1[i] < data1[i - 1]) && (data1[i] < data1[i + 1]))
                        {
                            double delta = Math.Sqrt(25 * data[i - 3] + 9 * data[i - 1] + 16 * data[i] + 9 * data[i + 1] + 25 * data[i + 3]);
                            //double delta1 = Math.Sqrt(25 * data[i - 4] + 9 * data[i - 2] + 16 * data[i - 1] + 9 * data[i] + 25 * data[i + 2]);
                            //double delta2 = Math.Sqrt(25 * data[i - 2] + 9 * data[i] + 16 * data[i + 1] + 9 * data[i + 2] + 25 * data[i + 4]);

                            if (-data1[i] / delta > 0.5)
                            {
                                presult = new PeakResult();
                                //求精确峰位
                                double channel = i;
                                int j;

                                for (j = -3; j <= 3; j++)
                                {
                                    if (data2[i + j] == 0)
                                    {
                                        channel = i + j;
                                        break;
                                    }
                                    else
                                    {
                                        if ((data2[i + j] < 0) && (data2[i + j + 1] > 0))
                                        {
                                            channel = i - data2[i + j] / (data2[i + j + 1] - data2[i + j]);
                                            break;
                                        }
                                    }
                                }
                                presult.Channel = channel;
                                //求峰区域
                                for (j = 0; j < i - startIndex; j++)
                                {
                                    if (data2[i - j] == 0)
                                    {
                                        presult.RegionStartChannel = i - j;
                                        break;
                                    }
                                    else
                                    {
                                        if ((data2[i - j] < 0) && (data2[i - j - 1] > 0))
                                        {
                                            presult.RegionStartChannel = i - j - 1 + data2[i - j - 1] / (data2[i - j - 1] - data2[i - j]);
                                            break;
                                        }
                                    }
                                }
                                for (j = i; j < endIndex; j++)
                                {
                                    if (data2[j] == 0)
                                    {
                                        presult.RegionEndChannel = j;
                                        break;
                                    }
                                    else
                                    {
                                        if ((data2[j] > 0) && (data2[j + 1] < 0))
                                        {
                                            presult.RegionEndChannel = j + data2[j] / (data2[j] - data2[j + 1]);
                                            break;
                                        }
                                    }
                                }
                                double n = Math.Max(presult.RegionEndChannel - presult.Channel, presult.Channel - presult.RegionStartChannel);
                                presult.RegionStartChannel = presult.Channel - n;
                                presult.RegionEndChannel = presult.Channel + n;
                                tmpresults.Add(presult);
                            }//if
                        }
                    }
                    results = (PeakResult[])tmpresults.ToArray(typeof(PeakResult));
                    break;
                case PeakType.covariance://协方差法
                    results = new PeakResult[0];
                    break;
                case PeakType.filter://匹配滤波器法
                    results = new PeakResult[0];
                    break;
                case PeakType.linearity://线性拟合法
                    results = new PeakResult[0];
                    break;
                default:
                    results = new PeakResult[0];
                    break;
            }
        }
        #endregion

        #region "最小二乘法曲线拟合"
        /// <summary>
        /// 最小二乘法
        /// </summary>
        /// <param name="data"></param>
        /// <returns></returns>
        public LeastSquareResult LeastSquare(MultinomialType type, params LeastSquareData[] data)
        {
            //只有1个点则产生异常
            if (data.Length < 2)
            {
                throw new Exception("Parameters is too less.");
            }
            double a1 = 0;
            double a2 = 0;
            double a3 = 0;
            double a4 = 0;
            double a5 = 0;
            double a6 = 0;
            double a7 = 0;
            double a8 = 0;
            double b1 = 0;
            double b2 = 0;
            double b3 = 0;
            double b4 = 0;
            double b5 = 0;
            double x1 = 0;
            double x2 = 0;
            double x3 = 0;
            LeastSquareResult result = new LeastSquareResult();

            switch (type)
            {
                case MultinomialType.Linear:
                    //*****要计算的矩阵*****
                    //|a1 a2||x1|   |b1|
                    //|a4 0 ||x2|   |b3|
                    //**********************
                    //计算矩阵系数
                    a1 = data.Length;
                    foreach (LeastSquareData d in data)
                    {
                        a2 += d.X;
                        a3 += d.X * d.X;
                        b1 += d.Fx;
                        b2 += d.Fx * d.X;
                    }
                    a4 = a2 - a1 * a3 / a2;
                    b3 = b2 - b1 * a3 / a2;
                    //计算结果
                    x1 = b3 / a4;
                    x2 = (b1 - a1 * x1) / a2;
                    result.A = x1;
                    result.B = x2;
                    result.C = 0;
                    break;
                case MultinomialType.Square:
                    //如果只有两个点则作直线拟合
                    if (data.Length < 3)
                    {
                        double b = 0;
                        double a = 0;
                        b = (data[0].Fx - data[1].Fx) / (data[0].X - data[1].X);
                        a = data[0].Fx - b * data[0].X;
                        result.A = a;
                        result.B = b;
                        result.C = 0;
                        return result;
                    }
                    else
                    {
                        //*****要计算的矩阵*****
                        //|a1 a2 a3||x1|   |b1|
                        //|a6 a7 0 ||x2| = |b4|
                        //|a8 0  0 ||x3|   |b5|
                        //**********************
                        //计算矩阵系数
                        a1 = data.Length;
                        foreach (LeastSquareData d in data)
                        {
                            a2 += d.X;