102 lines
3.5 KiB
C#
102 lines
3.5 KiB
C#
using MathNet.Numerics;
|
||
using MathNet.Numerics.IntegralTransforms;
|
||
using System;
|
||
using System.Linq;
|
||
using System.Numerics;
|
||
|
||
namespace 测试Scottplot
|
||
{
|
||
internal class ClarkTaylorModel
|
||
{
|
||
public static void Construct()
|
||
{
|
||
// 1. 生成模拟数据
|
||
double[] time = Enumerable.Range(0, 1000).Select(i => i * 0.001).ToArray();
|
||
double[] laserPulse = GenerateGaussianPulse(time, 0.05, 0.02);
|
||
double[] idealResponse = GenerateIdealResponse(time, 0.2); // 假设真实响应
|
||
double[] measuredTemp = Convolve(idealResponse, laserPulse); // 模拟实测数据
|
||
|
||
// 2. 反卷积
|
||
double[] deconvolved = DeconvolveFFT(measuredTemp, laserPulse);
|
||
|
||
// 3. 计算热扩散率
|
||
double tHalf = FindHalfMaxTime(deconvolved, time);
|
||
double alpha = CalculateThermalDiffusivity(1e-3, tHalf);
|
||
Console.WriteLine($"Thermal Diffusivity: {alpha} m²/s");
|
||
}
|
||
//生成高斯脉冲
|
||
static double[] GenerateGaussianPulse(double[] time, double peakTime, double width)
|
||
{
|
||
return time.Select(t => Math.Exp(-Math.Pow((t - peakTime) / width, 2))).ToArray();
|
||
}
|
||
|
||
// 生成理想响应（理论模型）
|
||
static double[] GenerateIdealResponse(double[] time, double alpha)
|
||
{
|
||
double d = 1e-3; // 厚度1mm
|
||
double t0 = 0.1388 * d * d / alpha; // 理论半峰时间
|
||
return time.Select(t => t >= t0 ? Math.Exp(-0.5 * Math.Pow((t - t0) / (0.1 * t0), 2)) : 0).ToArray();
|
||
}
|
||
|
||
// 卷积运算（模拟实测信号）
|
||
static double[] Convolve(double[] signal, double[] kernel)
|
||
{
|
||
int n = signal.Length;
|
||
double[] result = new double[n];
|
||
for (int i = 0; i < n; i++)
|
||
{
|
||
for (int j = 0; j <= i; j++)
|
||
{
|
||
result[i] += signal[j] * (i - j < kernel.Length ? kernel[i - j] : 0);
|
||
}
|
||
}
|
||
return result;
|
||
}
|
||
|
||
// FFT反卷积（含维纳滤波）
|
||
static double[] DeconvolveFFT(double[] measured, double[] pulse)
|
||
{
|
||
int n = measured.Length;
|
||
Complex[] fftMeasured = new Complex[n];
|
||
Complex[] fftPulse = new Complex[n];
|
||
|
||
for (int i = 0; i < n; i++)
|
||
{
|
||
fftMeasured[i] = new Complex(measured[i], 0);
|
||
fftPulse[i] = new Complex(pulse[i], 0);
|
||
}
|
||
|
||
Fourier.Forward(fftMeasured);
|
||
Fourier.Forward(fftPulse);
|
||
|
||
double noiseLevel = 1e-6;
|
||
Complex[] fftIdeal = new Complex[n];
|
||
for (int i = 0; i < n; i++)
|
||
{
|
||
double denom = fftPulse[i].MagnitudeSquared() + noiseLevel;
|
||
fftIdeal[i] = fftMeasured[i] * Complex.Conjugate(fftPulse[i]) / denom;
|
||
}
|
||
|
||
Fourier.Inverse(fftIdeal);
|
||
return fftIdeal.Select(c => c.Real).ToArray();
|
||
}
|
||
|
||
// 找半峰时间
|
||
static double FindHalfMaxTime(double[] response, double[] time)
|
||
{
|
||
double maxTemp = response.Max();
|
||
int halfMaxIndex = Array.FindIndex(response, t => t >= maxTemp / 2);
|
||
return time[halfMaxIndex];
|
||
}
|
||
|
||
// 计算热扩散率
|
||
static double CalculateThermalDiffusivity(double thickness, double tHalf)
|
||
{
|
||
return 0.1388 * Math.Pow(thickness, 2) / tHalf;
|
||
}
|
||
}
|
||
|
||
|
||
|
||
}
|