光学仪器示波器SCPI以及相位计算 需求二

需求二：利用算法，重新计算相位差
通过SCPI指令采集光学数据
1、常用指令如下：
:WAVeform:POINts? %查询采样总点数
:WAVeform:POINts:MODE RAW %在raw模式下设置采样点数
:WAVeform:FORMat %设置数据返回的格式
:WAVeform:SOURce? %设置查询采样的通道
:TIMebase:FORMat? %设置查询当前模式x-y、y-t :TIMebase:MODE?
:WAVeform:DATA? %设置读取数据
2、计算相位
internal string StandardMethodTimer(List listX, List listY)
{
try
{
int total;

     if (listX.Count > listY.Count)
     {
         total = listY.Count;
     }
     else
     {
         total = listX.Count;
     }

     double[] point_x = new double[total];
     double[] point_y = new double[total];

     for (int i = 0; i < total-1; i++)
     {
         double a = Convert.ToSingle(listX[i]);
         double b = Convert.ToSingle(listY[i]);
         point_x[i] = a;
         point_y[i] = b;
     }

     double A1_1 = 0; double A1_2 = 0; double A1_3 = 0; double A1_4 = 0; double A1_5 = 0;
     double A2_1 = 0; double A2_2 = 0; double A2_3 = 0; double A2_4 = 0; double A2_5 = 0;
     double A3_1 = 0; double A3_2 = 0; double A3_3 = 0; double A3_4 = 0; double A3_5 = 0;
     double A4_1 = 0; double A4_2 = 0; double A4_3 = 0; double A4_4 = 0; double A4_5 = 0;
     double A5_1 = 0; double A5_2 = 0; double A5_3 = 0; double A5_4 = 0; double A5_5 = 0;

     double B1 = 0; double B2 = 0; double B3 = 0; double B4 = 0; double B5 = 0;
     for (int i = 0; i < total-1; i++)
     {
         double xx = point_x[i] * point_x[i];
         double xxx = point_x[i] * point_x[i] * point_x[i];
         double yy = point_y[i] * point_y[i];
         double yyy = point_y[i] * point_y[i] * point_y[i];
         double xy = point_x[i] * point_y[i];

         A1_1 += Math.Pow(point_x[i], 4);
         A1_2 += Math.Pow(point_x[i], 3) * point_y[i];
         A1_3 += Math.Pow(point_x[i], 2) * Math.Pow(point_y[i], 2);
         A1_4 += Math.Pow(point_x[i], 3);
         A1_5 += Math.Pow(point_x[i], 2) * point_y[i];

         A2_1 += Math.Pow(point_x[i], 3) * point_y[i];
         A2_2 += Math.Pow(point_x[i], 2) * Math.Pow(point_y[i], 2);
         A2_3 += Math.Pow(point_x[i], 1) * Math.Pow(point_y[i], 3);
         A2_4 += Math.Pow(point_x[i], 2) * Math.Pow(point_y[i], 1);
         A2_5 += Math.Pow(point_x[i], 1) * Math.Pow(point_y[i], 2);

         A3_1 += Math.Pow(point_x[i], 2) * Math.Pow(point_y[i], 2);
         A3_2 += Math.Pow(point_x[i], 1) * Math.Pow(point_y[i], 3);
         A3_3 += Math.Pow(point_y[i], 4);
         A3_4 += Math.Pow(point_x[i], 1) * Math.Pow(point_y[i], 2);
         A3_5 += Math.Pow(point_y[i], 3);

         A4_1 += Math.Pow(point_x[i], 3);
         A4_2 += Math.Pow(point_x[i], 2) * Math.Pow(point_y[i], 1);
         A4_3 += Math.Pow(point_x[i], 1) * Math.Pow(point_y[i], 2);
         A4_4 += Math.Pow(point_x[i], 2);
         A4_5 += Math.Pow(point_x[i], 1) * Math.Pow(point_y[i], 1);

         A5_1 += Math.Pow(point_x[i], 2) * Math.Pow(point_y[i], 1);
         A5_2 += Math.Pow(point_x[i], 1) * Math.Pow(point_y[i], 2);
         A5_3 += Math.Pow(point_y[i], 3);
         A5_4 += Math.Pow(point_x[i], 1) * Math.Pow(point_y[i], 1);
         A5_5 += Math.Pow(point_y[i], 2);

         B1 += Math.Pow(point_x[i], 2);
         B2 += point_x[i] * point_y[i];
         B3 += Math.Pow(point_y[i], 2);
         B4 += point_x[i];
         B5 += point_y[i];
     }

     double[,] xa = { { A1_1, A1_2,A1_3,A1_4,A1_5 },
                  { A2_1, A2_2,A2_3,A2_4,A2_5 },
                  { A3_1, A3_2,A3_3,A3_4,A3_5 },
                  { A4_1, A4_2,A4_3,A4_4,A4_5 },
                  { A5_1, A5_2,A5_3,A5_4,A5_5},
 };
     var matrix = DenseMatrix.OfArray(xa);
     var matrix_xa = matrix.Inverse();

     double[,] xa_2 = { { B1, B2, B3, B4, B5 } };
     var matrix_2 = DenseMatrix.OfArray(xa_2);

     var resValuse = matrix_2 * matrix_xa;

     double A = -resValuse.Storage[0, 0];
     double B = -resValuse.Storage[0, 1];
     double C = -resValuse.Storage[0, 2];
     double D = -resValuse.Storage[0, 3];
     double E = -resValuse.Storage[0, 4];

     double Xc = (B * E - 2 * C * D) / (4 * A * C - B * B);
     double Yc = (B * D - 2 * A * E) / (4 * A * C - B * B);

     double FA1 = 2 * (A * Math.Pow(Xc, 2) + C * Math.Pow(Yc, 2) + B * Xc * Yc - 1);
     double FA2 = Math.Sqrt(Math.Pow(A - C, 2) + B * B);

     double a0 = Math.Sqrt(FA1 / (A + C + FA2));//长轴
     double b0 = Math.Sqrt(FA1 / (A + C - FA2));//短轴

     double alpha = 0.5 * Math.Atan(B / (A - C)) * 180 / Math.PI;
     double Theta = 90 + 0.5 * (-90 - alpha);

     Double phase = 2 * Math.Atan(b0 / a0);  //弧度，需要转为角度
     double tanAngleValue2 = phase / Math.PI * 180;//求角度

     return tanAngleValue2.ToString("f1");
 }
 catch(Exception ex) 
 {
     MessageBox.Show(ex.ToString());
     return null;
 }           

}