diff --git a/Crystallography.Controls/CheckLocation.cs b/Crystallography.Controls/CheckLocation.cs
new file mode 100644
index 0000000..23cb1aa
--- /dev/null
+++ b/Crystallography.Controls/CheckLocation.cs
@@ -0,0 +1,29 @@
+using System.Drawing;
+using System.Linq;
+using System.Windows.Forms;
+
+namespace Crystallography.Controls
+{
+ public static class CheckLocation
+ {
+ public static void Adjust(Form form)
+ {
+ double fX = form.Location.X, fY = form.Location.Y, fW = form.Size.Width, fH = form.Size.Height;
+
+ var needsAdjust = true;
+ foreach (var scr in Screen.AllScreens)
+ {
+ double bX = scr.Bounds.X, bY = scr.Bounds.Y, bW = scr.Bounds.Width, bH = scr.Bounds.Height;
+ //どれかのウィンドウの中に20%以上含まれていたら、needsAdjustをfalse (=更新の必要なし) にする。
+ if ((fX - bX) / fW > -0.8 && (bX + bW - fX - fW) / fW > -0.8 && (fY - bY) / fH > -0.8 && (bY + bH - fY - fH) / fH > -0.8)
+ needsAdjust = false;
+ }
+
+ if (needsAdjust)
+ {
+ var scr = Screen.AllScreens.First(e => e.DeviceName == Screen.FromControl(form).DeviceName);
+ form.Location = new Point(scr.Bounds.X + 100, scr.Bounds.Y + 100);
+ }
+ }
+ }
+}
diff --git a/Crystallography.Controls/Crystallography.Controls.csproj b/Crystallography.Controls/Crystallography.Controls.csproj
index df1b76b..289d76a 100644
--- a/Crystallography.Controls/Crystallography.Controls.csproj
+++ b/Crystallography.Controls/Crystallography.Controls.csproj
@@ -4,8 +4,8 @@
Library
net7.0-windows
true
- 2023.3.8.0528
- 2023.3.8.0528
+ 2023.3.18.0710
+ 2023.3.18.0710
PerMonitorV2
true
true
diff --git a/Crystallography.Controls/Numeric/NumericBox.cs b/Crystallography.Controls/Numeric/NumericBox.cs
index db4a6e8..19ee6d5 100644
--- a/Crystallography.Controls/Numeric/NumericBox.cs
+++ b/Crystallography.Controls/Numeric/NumericBox.cs
@@ -162,15 +162,14 @@ public string ToolTip
public string HeaderText { set => labelHeader.Text = value; get => labelHeader.Text; }
[Category("Font && Color")]
+ [Localizable(true)]
[DefaultValue(typeof(Padding), "0,0,0,0")]
-
public Padding HeaderMargin { set => labelHeader.Margin = value; get => labelHeader.Margin; }
[Localizable(true)]
[DefaultValue(typeof(Font), "Segoe UI Symbol, 9.75pt")]
[Category("Font && Color")]
-
public Font HeaderFont { set => labelHeader.Font = value; get => labelHeader.Font; }
[DefaultValue(typeof(Color), "ControlText")]
@@ -218,6 +217,7 @@ public string ToolTip
public Color TextBoxBackColor { set => textBox.BackColor = value; get => textBox.BackColor; }
[DefaultValue(typeof(Font), "Segoe UI Symbol, 9.75pt")]
+ [Localizable(true)]
[Category("Font && Color")]
///
/// font
diff --git a/Crystallography/Atom/AtomStatic.cs b/Crystallography/Atom/AtomStatic.cs
index 76d3fc8..a7e5b0e 100644
--- a/Crystallography/Atom/AtomStatic.cs
+++ b/Crystallography/Atom/AtomStatic.cs
@@ -2667,8 +2667,8 @@ public double FactorImaginaryAnnular(double kV, Vector3DBase g, double m, double
f_kPlusG += A * Math.Exp(-kPlusG * B / 100);
}
return f_kMinusG * f_kPlusG * (1 - Math.Exp(m * (gLen2 - kMinusG - kPlusG))); ;// * sinThetaを外に出して、少しでも早く
- }, 0, 2 * Math.PI, 30) * sinθ;
- }, inner, outer, 80);
+ }, 0, 2 * Math.PI, 20) * sinθ;
+ }, inner, outer, 60);
return gamma * k0 / 2 * result * 0.01;
}
diff --git a/Crystallography/Atom/Atoms.cs b/Crystallography/Atom/Atoms.cs
index 734eea5..edfbff0 100644
--- a/Crystallography/Atom/Atoms.cs
+++ b/Crystallography/Atom/Atoms.cs
@@ -2,7 +2,6 @@
using System.Drawing;
using System.Numerics;
using System.Xml.Serialization;
-using static System.Net.Mime.MediaTypeNames;
namespace Crystallography;
@@ -554,6 +553,11 @@ public enum Type { U, B }
///
public double Biso000 => (B11 * a2 + B22 * b2 + B33 * c2 + 2 * B12 * ab + 2 * B23 * bc + 2 * B31 * ca) * 4.0 / 3.0;
+ ///
+ /// xq[̏ꍇtrue
+ ///
+ public bool IsZero => UseIso ? Biso == 0 : B11 == 0 && B22 == 0 && B33 == 0 && B12 == 0 && B23 == 0 && B31 == 0;
+
///
/// unit: nm^2
///
@@ -607,28 +611,68 @@ public enum Type { U, B }
///
public double B31_err => OriginalType == Type.B ? Aniso31_err : Aniso31_err * coeff31;
-
-
#endregion
#region U type. Get̂
+ ///
+ /// Pʂ nm^2
+ ///
public double Uiso => OriginalType == Type.U ? Iso : Iso / PI2 / 8;
+ ///
+ /// Pʂ nm^2
+ ///
public double Uiso_err => OriginalType == Type.U ? Iso_err : Iso_err / PI2 / 8;
+ ///
+ /// Pʂ nm^2
+ ///
public double U11 => OriginalType == Type.U ? Aniso11 : Aniso11 / coeff11;
+ ///
+ /// Pʂ nm^2
+ ///
public double U22 => OriginalType == Type.U ? Aniso22 : Aniso22 / coeff22;
+ ///
+ /// Pʂ nm^2
+ ///
public double U33 => OriginalType == Type.U ? Aniso33 : Aniso33 / coeff33;
+ ///
+ /// Pʂ nm^2
+ ///
public double U12 => OriginalType == Type.U ? Aniso12 : Aniso12 / coeff12;
+ ///
+ /// Pʂ nm^2
+ ///
public double U23 => OriginalType == Type.U ? Aniso23 : Aniso23 / coeff23;
+ ///
+ /// Pʂ nm^2
+ ///
public double U31 => OriginalType == Type.U ? Aniso31 : Aniso31 / coeff31;
+ ///
+ /// Pʂ nm^2
+ ///
public double U11_err => OriginalType == Type.U ? Aniso11_err : Aniso11_err / coeff11;
+ ///
+ /// Pʂ nm^2
+ ///
public double U22_err => OriginalType == Type.U ? Aniso22_err : Aniso22_err / coeff22;
+ ///
+ /// Pʂ nm^2
+ ///
public double U33_err => OriginalType == Type.U ? Aniso33_err : Aniso33_err / coeff33;
+ ///
+ /// Pʂ nm^2
+ ///
public double U12_err => OriginalType == Type.U ? Aniso12_err : Aniso12_err / coeff12;
+ ///
+ /// Pʂ nm^2
+ ///
public double U23_err => OriginalType == Type.U ? Aniso23_err : Aniso23_err / coeff23;
+ ///
+ /// Pʂ nm^2
+ ///
public double U31_err => OriginalType == Type.U ? Aniso31_err : Aniso31_err / coeff31;
#endregion
- #region IWi̒l
+ #region IWi̒l
///
/// Pʂ nm^2
///
diff --git a/Crystallography/BetheMethod.cs b/Crystallography/BetheMethod.cs
index 0ca2d27..3b363be 100644
--- a/Crystallography/BetheMethod.cs
+++ b/Crystallography/BetheMethod.cs
@@ -1,30 +1,25 @@
#region using
-
using MathNet.Numerics;
using MathNet.Numerics.LinearAlgebra;
using MathNet.Numerics.LinearAlgebra.Complex;
-using OpenTK.Graphics.OpenGL;
+using OpenTK.Platform.Windows;
using System;
using System.Buffers;
-using System.CodeDom.Compiler;
using System.Collections.Generic;
using System.ComponentModel;
-using System.Diagnostics;
using System.Drawing;
-using System.Globalization;
using System.Linq;
using System.Numerics;
-using System.Reflection;
using System.Runtime.InteropServices;
-using System.Text;
using System.Threading;
using System.Threading.Tasks;
-using System.Windows.Forms;
using System.Xml.Serialization;
+
using static System.Buffers.ArrayPool;
using static System.Numerics.Complex;
using DMat = MathNet.Numerics.LinearAlgebra.Complex.DenseMatrix;
using DVec = MathNet.Numerics.LinearAlgebra.Complex.DenseVector;
+
#endregion
namespace Crystallography;
@@ -119,9 +114,9 @@ public enum Solver { Eigen_MKL, Eigen_Eigen, MtxExp_MKL, MtxExp_Eigen, Auto }
static BetheMethod()
{
EigenEnabled = NativeWrapper.Enabled;
- BlasEnabled = MathNet.Numerics.Control.TryUseNativeOpenBLAS();
- MklEnabled = MathNet.Numerics.Control.TryUseNativeMKL();
- CudaEnabled = MathNet.Numerics.Control.TryUseNativeCUDA();
+ BlasEnabled = Control.TryUseNativeOpenBLAS();
+ MklEnabled = Control.TryUseNativeMKL();
+ CudaEnabled = Control.TryUseNativeCUDA();
}
public BetheMethod(Crystal crystal)
{
@@ -166,8 +161,8 @@ public void CancelCBED()
/// 基準となる方位
/// 厚みの配列
/// 基準となる方位に乗算する方位配列
- public void RunCBED(int maxNumOfBloch, double voltage, Matrix3D rotation,
- double[] thickness, Vector3DBase[] beamDirections,bool LACBED, Solver solver = Solver.Auto, int thread = 1)
+ public void RunCBED(int maxNumOfBloch, double voltage, Matrix3D rotation,
+ double[] thickness, Vector3DBase[] beamDirections, bool LACBED, Solver solver = Solver.Auto, int thread = 1)
{
MaxNumOfBloch = maxNumOfBloch;
AccVoltage = voltage;
@@ -348,21 +343,21 @@ private unsafe void cbed_DoWork(object sender, DoWorkEventArgs e)
//まず、各ディスクを構成するピクセルの座標を計算
var diskTemp = new (RectangleD Rect, PointD[] Pos)[Beams.Length];
Parallel.For(0, Beams.Length, g =>
+ {
+ if (!bwCBED.CancellationPending)
{
- if (!bwCBED.CancellationPending)
+ var pos = new PointD[BeamDirections.Length];
+ for (int r = 0; r < pos.Length; r++)
{
- var pos = new PointD[BeamDirections.Length];
- for (int r = 0; r < pos.Length; r++)
- {
- //Ewald球中心(試料)から見た、逆格子ベクトルの方向
- var vec = kvac * BeamDirections[r] + Disks[0][g].G;
- //var vec = BeamDirections[r] * (new Vector3DBase(0, 0, kvac) - Disks[0][g].G);
- //var vec = BeamDirections[r] - Disks[0][g].G;
- pos[r] = new PointD(vec.X / vec.Z, vec.Y / vec.Z); //カメラ長 1 を想定した検出器上のピクセルの座標値を格納
- }
- diskTemp[g] = (new RectangleD(new PointD(pos.Min(p => p.X), pos.Min(p => p.Y)), new PointD(pos.Max(p => p.X), pos.Max(p => p.Y))), pos);
+ //Ewald球中心(試料)から見た、逆格子ベクトルの方向
+ var vec = kvac * BeamDirections[r] + Disks[0][g].G;
+ //var vec = BeamDirections[r] * (new Vector3DBase(0, 0, kvac) - Disks[0][g].G);
+ //var vec = BeamDirections[r] - Disks[0][g].G;
+ pos[r] = new PointD(vec.X / vec.Z, vec.Y / vec.Z); //カメラ長 1 を想定した検出器上のピクセルの座標値を格納
}
- });
+ diskTemp[g] = (new RectangleD(new PointD(pos.Min(p => p.X), pos.Min(p => p.Y)), new PointD(pos.Max(p => p.X), pos.Max(p => p.Y))), pos);
+ }
+ });
//g1のディスク中のピクセル(r1)に対して、他のディスク(g2)の重なるピクセル(r2)を足し合わせていく。
Parallel.For(0, Beams.Length, g1 =>
@@ -713,31 +708,31 @@ public Beam[] GetPrecessionElectronDiffraction(int maxNumOfBloch, double voltage
gDic.Clear();
stepP.ForAll(k =>
- {
- var rotAngle = 2.0 * Math.PI * k / step;
- var beamRotation = Matrix3D.Rot(new Vector3DBase(Math.Cos(rotAngle), Math.Sin(rotAngle), 0), SemianglePED);
- //計算対象のg-Vectorsを決める。
- var potentialMatrix = Array.Empty();
- var vecK0 = getVecK0(kvac, u0, beamRotation * new Vector3D(0, 0, -1));
- BeamsPED[k] = Find_gVectors(BaseRotation, vecK0, MaxNumOfBloch, true);
- var len = BeamsPED[k].Length;
- potentialMatrix = getEigenMatrix(BeamsPED[k]);
- var dim = BeamsPED[k].Length;
- //A行列に関する固有値、固有ベクトルを取得
- if (useEigen)
- {//Eigenを使う場合
- var (val, vec) = NativeWrapper.EigenSolver(dim, potentialMatrix);
- (EigenValuesPED[k], EigenVectorsPED[k]) = (new DVec(val), new DMat(len, len, vec));
- EigenVectorsInversePED[k] = new DMat(len, len, NativeWrapper.Inverse(len, EigenVectorsPED[k].Values));
- }
- else
- {//MKLを使う場合
- var evd = new DMat(dim, dim, potentialMatrix).Evd(Symmetricity.Asymmetric);
- EigenValuesPED[k] = (DVec)evd.EigenValues;
- EigenVectorsPED[k] = (DMat)evd.EigenVectors;
- EigenVectorsInversePED[k] = (DMat)EigenVectorsPED[k].Inverse();
- }
- });
+ {
+ var rotAngle = 2.0 * Math.PI * k / step;
+ var beamRotation = Matrix3D.Rot(new Vector3DBase(Math.Cos(rotAngle), Math.Sin(rotAngle), 0), SemianglePED);
+ //計算対象のg-Vectorsを決める。
+ var potentialMatrix = Array.Empty();
+ var vecK0 = getVecK0(kvac, u0, beamRotation * new Vector3D(0, 0, -1));
+ BeamsPED[k] = Find_gVectors(BaseRotation, vecK0, MaxNumOfBloch, true);
+ var len = BeamsPED[k].Length;
+ potentialMatrix = getEigenMatrix(BeamsPED[k]);
+ var dim = BeamsPED[k].Length;
+ //A行列に関する固有値、固有ベクトルを取得
+ if (useEigen)
+ {//Eigenを使う場合
+ var (val, vec) = NativeWrapper.EigenSolver(dim, potentialMatrix);
+ (EigenValuesPED[k], EigenVectorsPED[k]) = (new DVec(val), new DMat(len, len, vec));
+ EigenVectorsInversePED[k] = new DMat(len, len, NativeWrapper.Inverse(len, EigenVectorsPED[k].Values));
+ }
+ else
+ {//MKLを使う場合
+ var evd = new DMat(dim, dim, potentialMatrix).Evd(Symmetricity.Asymmetric);
+ EigenValuesPED[k] = (DVec)evd.EigenValues;
+ EigenVectorsPED[k] = (DMat)evd.EigenVectors;
+ EigenVectorsInversePED[k] = (DMat)EigenVectorsPED[k].Inverse();
+ }
+ });
}
//各方向でのbeamの振幅を求める
@@ -804,7 +799,7 @@ public void CancelSTEM()
if (bwSTEM.IsBusy)
bwSTEM.CancelAsync();
}
- public void RunSTEM(int maxNumOfBloch, double voltage, double cs, double delta, double sliceThickness, Size imageSize, double resolution,
+ public void RunSTEM(int maxNumOfBloch, double voltage, double cs, double delta, double sliceThickness, Size imageSize, double resolution, double sourceSize,
Matrix3D baseRotation, double[] thicknesses, double[] defocusses,
Vector3DBase[] beamDirections, double convergenceAngle, double detAngleInner, double detAngleOuter,
bool calcElas, bool calcInel, Solver solver = Solver.Auto, int thread = 1)
@@ -817,14 +812,15 @@ public void RunSTEM(int maxNumOfBloch, double voltage, double cs, double delta,
BaseRotation = new Matrix3D(baseRotation);
BeamDirections = beamDirections;
Thicknesses = thicknesses;
- if(!bwSTEM.IsBusy)
- bwSTEM.RunWorkerAsync((solver, thread, cs, delta, sliceThickness, convergenceAngle, detAngleInner, detAngleOuter, thicknesses, defocusses, imageSize, resolution, calcElas, calcInel));
+ if (!bwSTEM.IsBusy)
+ bwSTEM.RunWorkerAsync((solver, thread, cs, delta, sliceThickness, convergenceAngle, detAngleInner, detAngleOuter, thicknesses, defocusses, imageSize, resolution, sourceSize, calcElas, calcInel));
}
- public void stem_DoWork(object sender, DoWorkEventArgs e)
+ public unsafe void stem_DoWork(object sender, DoWorkEventArgs e)
{
//MathNetの行列の内部は、1列目の要素、2列目の要素、という順番で格納されている
- var (solver, thread, cs, delta, sliceThickness, convergenceAngle, detAngleInner, detAngleOuter, thicknesses, defocusses, imageSize, resolution, calcElas, calcInel)
- = ((Solver, int, double, double, double, double, double, double, double[], double[], Size, double, bool, bool))e.Argument;
+
+ var (solver, thread, cs, delta, sliceThickness, convergenceAngle, detAngleInner, detAngleOuter, thicknesses, defocusses, imageSize, resolution, sourceSize, calcElas, calcInel)
+ = ((Solver, int, double, double, double, double, double, double, double[], double[], Size, double, double, bool, bool))e.Argument;
var diameterPix = (int)Math.Sqrt(BeamDirections.Length);
var radiusPix = diameterPix / 2.0;
@@ -877,7 +873,7 @@ public void stem_DoWork(object sender, DoWorkEventArgs e)
int dLen = defocusses.Length, tLen = thicknesses.Length, bLen = Beams.Length;
//入射面での波動関数を定義
- var psi0 = Enumerable.Range(0, Beams.Length).ToList().Select(g => g == 0 ? One : 0).ToArray();
+ var psi0 = Enumerable.Range(0, Beams.Length).ToList().Select(g => g == 0 ? -One : 0).ToArray();
//ポテンシャルマトリックスを初期化
uDictionary.Clear();
var potentialMatrix = getPotentialMatrix(Beams);
@@ -903,16 +899,17 @@ public void stem_DoWork(object sender, DoWorkEventArgs e)
var kg_z = new double[BeamDirections.Length][];
//進捗状況報告用の各種定数を初期化
int count = 0;
- //Transmission coefficient tc[k][t][g]
var tc = BeamDirections.AsParallel().WithDegreeOfParallelism(thread).Select((beamDirection, i) =>
{
var vecK0 = getVecK0(kvac, u0, beamDirection);
k_vec[i] = vecK0;
- if (!inside(i) || bwSTEM.CancellationPending) return null;
- if (Interlocked.Increment(ref count) % 10 == 0) bwSTEM.ReportProgress(count, reportString);//進捗状況を報告
+ if (Interlocked.Increment(ref count) % 10 == 0) bwSTEM.ReportProgress((int)(1E6 * count / BeamDirections.Count()), reportString);//進捗状況を報告
+ if (bwSTEM.CancellationPending) { e.Cancel = true; return null; }
- var eigenMatrix = new Complex[bLen * bLen];// Shared.Rent(bLen * bLen);
+ if (!inside(i)) return null;
+
+ var eigenMatrix = GC.AllocateUninitializedArray(bLen * bLen);// Shared.Rent(bLen * bLen);
var beams = ArrayPool.Shared.Rent(bLen);
try
{
@@ -934,7 +931,7 @@ public void stem_DoWork(object sender, DoWorkEventArgs e)
var tg = new DMat(bLen, tLen);
for (int t = 0; t < tLen; t++)
{
- var gammmaAlpha = new DVec(eVal[i].Select((ev, g) => Exp(TwoPiI * ev * thicknesses[t] ) * α[i][g]).ToArray());//ガンマの対称行列×アルファを作成
+ var gammmaAlpha = new DVec(eVal[i].Select((ev, g) => Exp(TwoPiI * ev * thicknesses[t]) * α[i][g]).ToArray());//ガンマの対称行列×アルファを作成
tg.SetColumn(t, evd.EigenVectors * gammmaAlpha);//深さtにおけるψを求める
}
result = tg.Values;
@@ -950,13 +947,11 @@ public void stem_DoWork(object sender, DoWorkEventArgs e)
_tc[t][g] = result[t * bLen + g] * Exp(TwoPiI * kg_z[i][g] * thicknesses[t]);
return _tc;
}
- finally { /*Shared.Return(eigenMatrix);*/ ArrayPool.Shared.Return(beams); }
+ finally { ArrayPool.Shared.Return(beams); }
}).ToArray();
#endregion
- if (bwSTEM.CancellationPending) { e.Cancel = true; return; }
-
var k_xy = k_vec.Select(e => e.ToPointD).ToArray();
var k_z = k_vec.Select(e => e.Z).ToArray();
@@ -991,7 +986,7 @@ Complex Lenz(in PointD k, in PointD kq, in double defocus)
var mat = BaseRotation * Crystal.MatrixInverse.Transpose();
var qList = Beams.AsParallel().SelectMany(e1 => Beams.Select(e2 => (e1 - e2).Index)).Distinct()
.Select(e => new Beam(e, mat * e)).Where(e => k_xy.Any(e2 => A(e2) && A(e2 + e.Vec.ToPointD))).OrderBy(e => e.Vec.Length2).ToList();
-
+
//if(qList.Count > Beams.Length)
// qList.RemoveRange(Beams.Length, qList.Count - Beams.Length);
@@ -1030,24 +1025,6 @@ Complex Lenz(in PointD k, in PointD kq, in double defocus)
#endregion
- #region U行列の計算
- count = 0;
- var U = new Complex[qList.Count][];
- uDictionary.Clear();
- if (calcInel) //U行列を作成
- Parallel.For(0, qList.Count, m =>
- {
- U[m] = new Complex[bLen * bLen];
- int k = 0;
- for (int i = 0; i < bLen; i++)
- for (int j = 0; j < bLen; j++)
- U[m][k++] = getU(AccVoltage, qList[m] + Beams[i] - Beams[j], null, detAngleInner, detAngleOuter).Imag;//局所形式の場合 i, jの順番が正解
- //U[m][k++] = getU(AccVoltage, qList[m], -Beams[i] + Beams[j], detAngleInner, detAngleOuter).Imag;//非局所形式の場合
- bwSTEM.ReportProgress((int)(1000.0 * Interlocked.Increment(ref count) / qList.Count), "Calculating U matrix");//状況を報告
- if (bwSTEM.CancellationPending) { e.Cancel = true; return; }
- });
- #endregion
-
//必要な情報だけを追加してParallelにしたtcP
var tcP = tc.AsParallel().Select((e, i) => (index: i, result: e, xy: k_xy[i])).Where(e => e.result != null && A(e.xy)).Select(e => e.index);//.WithDegreeOfParallelism(1);
@@ -1077,7 +1054,9 @@ Complex Lenz(in PointD k, in PointD kq, in double defocus)
});
#endregion
+
#region 弾性散乱 の計算
+ bwSTEM.ReportProgress(0, "Calculating I_elastic(Q)");//状況を報告
Complex[,,] I_Elas = new Complex[qList.Count, tLen, dLen];
count = 0;
if (calcElas)
@@ -1096,20 +1075,58 @@ Complex Lenz(in PointD k, in PointD kq, in double defocus)
I_Elas[qIndex, t, d] += i_Elas[t] * lenz[d];
}
if (bwSTEM.CancellationPending) return;
- if (Interlocked.Increment(ref count) % 10 == 0) bwSTEM.ReportProgress(count, "Calculating I_elastic(Q)");//状況を報告
+ if (Interlocked.Increment(ref count) % 10 == 0) bwSTEM.ReportProgress((int)(1E6 * count / tcP.Count()), "Calculating I_elastic(Q)");//状況を報告
});
#endregion
#region 非弾性散乱を計算する場合
var I_Inel = new Complex[qList.Count, tLen, dLen];
- var PiecewiseQuadrature = true;
-
if (calcInel)
{
+ bwSTEM.ReportProgress(0, "Calculating U matrix");//状況を報告
+ var bLen2 = bLen * bLen;
+ #region U行列の計算
+ count = 0;
+ //var U = new Complex[qList.Count][];
+ var U = new Complex[qList.Count * bLen2];
+ uDictionary.Clear();
+
+ //マルチスレッドの効率を上げるため、まずqList[qIndex] + Beams[i] - Beams[j]の重複を除く
+ var tmpDic = new Dictionary<(int h, int k, int l), (Beam b, int q, int i, int j)>();
+ for (int q = 0; q < qList.Count; q++)
+ for (int j = 0; j < bLen; j++)
+ for (int i = 0; i < bLen; i++)
+ {
+ var b = qList[q] + Beams[i] - Beams[j];
+ tmpDic.TryAdd(b.Index, (b, q, i, j));
+ }
+ tmpDic.AsParallel().ForAll(d =>
+ {
+ getU(AccVoltage, d.Value.b, null, detAngleInner, detAngleOuter);//共役とると、なぜかいい感じ。
+ if (Interlocked.Increment(ref count) % 10 == 0) bwSTEM.ReportProgress((int)(1E6 * count / tmpDic.Count), "Calculating U matrix");//状況を報告
+ if (bwSTEM.CancellationPending) { e.Cancel = true; return; }
+ });
+
+ Parallel.For(0, qList.Count, qIndex =>
+ {
+ //U[qIndex] = GC.AllocateUninitializedArray(bLen * bLen);
+ for (int j = 0; j < bLen; j++)
+ for (int i = 0; i < bLen; i++)
+ {
+ //U[qIndex][j * bLen + i] = getU(AccVoltage, qList[qIndex] + Beams[i] - Beams[j], null, detAngleInner, detAngleOuter).Imag.Conjugate();//共役とると、なぜかいい感じ。
+ U[qIndex* bLen2 + j * bLen + i] = getU(AccVoltage, qList[qIndex] + Beams[i] - Beams[j], null, detAngleInner, detAngleOuter).Imag.Conjugate();//共役とると、なぜかいい感じ。
+ //U[m][k++] = getU(AccVoltage, qList[m], -Beams[i] + Beams[j], detAngleInner, detAngleOuter).Imag;//非局所形式の場合
+ }
+ });
+ #endregion
+
+ var PiecewiseQuadrature = true;
if (PiecewiseQuadrature)
#region 区分求積法アルゴリズム
{
+ bwSTEM.ReportProgress(0, "Calculating I_inelastic(Q)");
+
#region 各厚みを、指定された厚み程度で切り分ける
var _thick = new double[Thicknesses.Length][];
var tStep = new double[Thicknesses.Length];
@@ -1133,18 +1150,22 @@ Complex Lenz(in PointD k, in PointD kq, in double defocus)
#endregion
#region 各種変数の設定
- var sum = new Complex[qList.Count, dLen];
+
+
+ var tc_k = GC.AllocateUninitializedArray(tc.Length * bLen);
- var tc_k = tc.Select(e => new Complex[bLen]).ToArray();
var validTc = list.Where(e1 => e1 != null).SelectMany(e2 => e2.SelectMany(e3 => e3.N)).Distinct().ToList().AsParallel();
var total = _thick.Sum(e => e.Length) * tcP.Count();
+
+
count = 0;
#endregion
#region メインのループ
for (int t = 0; t < Thicknesses.Length; t++)
{
+ var sum = new Complex[qList.Count * dLen];//ゼロ初期化が必要
foreach (var thickness in _thick[t])
{
if (bwSTEM.CancellationPending) return;
@@ -1165,146 +1186,152 @@ Complex Lenz(in PointD k, in PointD kq, in double defocus)
// for (int j = 0; j < bLen; j++)
// tc_k[kIndex][g] += eVec[kIndex][j * bLen + g] * exp_kgz[g] * exp_λ[j];
#endregion
- NativeWrapper.GenerateTC(bLen, thickness, kg_z[kIndex], eVal[kIndex], eVec[kIndex], ref tc_k[kIndex]);
+ fixed (Complex* _tc_k = tc_k, _eVal = eVal[kIndex], _eVec = eVec[kIndex])
+ fixed (double* _kg_z = kg_z[kIndex])
+ NativeWrapper.GenerateTC1(bLen, thickness, _kg_z, _eVal, _eVec, _tc_k + kIndex * bLen);
});
#endregion
tcP.ForAll(kIndex =>
{
- Complex[] TgThU2 = Shared.Rent(list[kIndex].Count * dLen), tc_kq = Shared.Rent(bLen);//, u_tck = Shared.Rent(list[kIndex].Count * bLen);
+ Complex[] sumTmp = Shared.Rent(list[kIndex].Count * dLen), tc_kq = Shared.Rent(bLen);
try
{
- for (int i = 0, j = 0; i < list[kIndex].Count; i++)
- {
- var (qIndex, n, r, lenz) = list[kIndex][i];
- //厚み_thick[t][_t]における透過係数_tc_kqを計算
- NativeWrapper.BlendAndConjugate(bLen, tc_k[n[0]], tc_k[n[1]], tc_k[n[2]], tc_k[n[3]], r[0], r[1], r[2], r[3], ref tc_kq);
- #region この内容をNativeコードで実行 4倍速い
- //var sum1 = new Complex(0, 0);
- //for (int h = 0; h < bLen; h++)
- // for (int g = 0; g < bLen; g++)
- // sum1 += tc_k[kIndex][g] * tc_kq[h] * U[qIndex][g * bLen + h];
- #endregion
- var TgThU1 = NativeWrapper.RowVec_SqMat_ColVec(bLen, tc_kq, U[qIndex], tc_k[kIndex]);
- for (int d = 0; d < dLen; d++)
- TgThU2[j++] = TgThU1 * lenz[d];
- }
+ fixed (Complex* _tc_k = tc_k, _U = U, _tc_kq = tc_kq)
+ for (int i = 0; i < list[kIndex].Count; i++)
+ {
+ var (qIndex, n, r, lenz) = list[kIndex][i];
+ //厚み_thick[t][_t]における透過係数_tc_kqを計算
+ NativeWrapper.BlendAndConjugate(bLen, _tc_k + n[0] * bLen, _tc_k + n[1] * bLen, _tc_k + n[2] * bLen, _tc_k + n[3] * bLen, r[0], r[1], r[2], r[3], _tc_kq);
+ var tmp = NativeWrapper.RowVec_SqMat_ColVec(bLen, _tc_kq, _U + qIndex * bLen2, _tc_k + kIndex * bLen);
+
+ for (int dIndex = 0; dIndex < dLen; dIndex++)
+ sumTmp[i * dLen + dIndex] = tmp * lenz[dIndex];
+ }
lock (lockObj1)
- for (int i = 0, j = 0; i < list[kIndex].Count; i++)
+ for (int i = 0; i < list[kIndex].Count; i++)
for (int d = 0; d < dLen; d++)
- sum[list[kIndex][i].qIndex, d] += TgThU2[j++] * tStep[t];
+ sum[list[kIndex][i].qIndex * dLen + d] += sumTmp[i * dLen + d] ;
}
- finally { Shared.Return(tc_kq); Shared.Return(TgThU2); }
+ finally { Shared.Return(sumTmp); Shared.Return(tc_kq); }
- if (Interlocked.Increment(ref count) % 1000 == 0) bwSTEM.ReportProgress((int)(1000000.0 / total * count), "Calculating I_inelastic(Q)");//状況を報告
+ if (Interlocked.Increment(ref count) % 1000 == 0) bwSTEM.ReportProgress((int)(1E6 / total * count), "Calculating I_inelastic(Q)");//状況を報告
});
+
}
- for (int qIndex = 0; qIndex < qList.Count; qIndex++)
- for (int d = 0; d < dLen; d++)
- I_Inel[qIndex, t, d] = sum[qIndex, d] / kvac;
+ var coeff = 2 * Math.PI / kvac * tStep[t];
+ Parallel.For(0, qList.Count, qIndex =>
+ {
+ for (int dIndex = 0; dIndex < dLen; dIndex++)
+ {
+ I_Inel[qIndex, t, dIndex] = sum[qIndex * dLen + dIndex] * coeff;
+ if (t > 0)
+ I_Inel[qIndex, t, dIndex] += I_Inel[qIndex, t - 1, dIndex];
+ }
+ });
}
+
#endregion
}
#endregion
- else
#region 解析的に非弾性を計算する場合
- {
- //ZOLZのみだったらこれでいいが、HOLZや晶帯軸から傾いている場合に対応できていない。
- //対応しようとするとものすごい計算コスト。あきらめるか。
-
- //固有値・ベクトルをブレンドするのではなく、最後にブレンドする。
-
- //最初にeVecにαを掛けておく
- Parallel.For(0, tc.Length, kIndex =>
- {
- if (eVal[kIndex] != null)
- for (int col = 0; col < bLen; col++)
- for (int row = 0; row < bLen; row++)
- eVec[kIndex][col * bLen + row] *= α[kIndex][col];
- });
-
- //複素共役なC, λ, αを用意
- Complex[][] C = eVec, _C = new Complex[tc.Length][], λ = eVal, _λ = new Complex[tc.Length][];
- Complex[][] exp = new Complex[tc.Length][], _exp= new Complex[tc.Length][];
- list.AsParallel().Where(e1 => e1 != null).SelectMany(e2 => e2.SelectMany(e3 => e3.N)).Distinct().ForAll(kIndex =>
- {
- _C[kIndex] = (new DMat(bLen, bLen, C[kIndex]).ConjugateTranspose() as DMat).Values;
- _λ[kIndex] = (new DVec(λ[kIndex]).Conjugate() as DVec).Values;
-
- exp[kIndex] = new Complex[tLen * bLen];
- _exp[kIndex] = new Complex[tLen * bLen];
- for (int t = 0; t < tLen; t++)
- for (int j = 0; j < bLen; j++)
- {
- exp[kIndex][t * bLen + j] = Exp(TwoPiI * (λ[kIndex][j] + kg_z[kIndex][j]) * Thicknesses[t]);
- _exp[kIndex][t * bLen + j] = exp[kIndex][t * bLen + j].Conjugate();
- }
- });
-
-
- var total = tcP.Count();
- count = 0;
- tcP.ForAll(kIndex =>
- {
- if (bwSTEM.CancellationPending) return;
-
- Complex[] λ_k = λ[kIndex], exp_k = exp[kIndex];
- double[] kz_k = kg_z[kIndex];
-
- Complex[] TDS = Shared.Rent(bLen * bLen), tmpMat = Shared.Rent(bLen * bLen);
- try
- {
- foreach (var (qIndex, n, r, lenz) in list[kIndex])
- {
- NativeWrapper.MultiplyMxM(bLen, U[qIndex], C[kIndex], ref tmpMat);
- var tmpSum = new Complex[tLen];
- for (int m = 0; m < n.Length; m++)
- {
- NativeWrapper.MultiplyMxM(bLen, _C[n[m]], tmpMat, ref TDS);
-
- Complex[] λ_kq = _λ[n[m]], exp_kq = _exp[n[m]];
- double[] kz_kq = kg_z[n[m]];
-
- //B行列の中身を計算し、アダマール積を取る
- for (int t = 0; t < tLen; t++)
- {
- int l = 0;
- for (int j = 0; j < bLen; j++)
- for (int i = 0; i < bLen; i++)
- tmpSum[t] += r[m] * ((exp_k[t * bLen + j] * exp_kq[t * bLen + i] - 1) / TwoPiI / (kz_k[j] - kz_kq[i] + λ_k[j] - λ_kq[i])) * TDS[l++];//B行列は作らず、直接アダマール積を取る
- }
- }
- lock (lockObj2)
- for (int t = 0; t < tLen; t++)
- for (int d = 0; d < dLen; d++)
- I_Inel[qIndex, t, d] += tmpSum[t] / kvac * lenz[d];
- }
- }
- finally { Shared.Return(TDS); Shared.Return(tmpMat); }
- if (Interlocked.Increment(ref count) % 10 == 0) bwSTEM.ReportProgress((int)(1000000.0 / total * count), "Calculating I_inelastic(Q)");//状況を報告
- });
-
- }
+ //else
+ //{
+ // //ZOLZのみだったらこれでいいが、HOLZや晶帯軸から傾いている場合に対応できていない。
+ // //対応しようとするとものすごい計算コスト。あきらめるか。
+
+ // //固有値・ベクトルをブレンドするのではなく、最後にブレンドする。
+
+ // //最初にeVecにαを掛けておく
+ // Parallel.For(0, tc.Length, kIndex =>
+ // {
+ // if (eVal[kIndex] != null)
+ // for (int col = 0; col < bLen; col++)
+ // for (int row = 0; row < bLen; row++)
+ // eVec[kIndex][col * bLen + row] *= α[kIndex][col];
+ // });
+
+ // //複素共役なC, λ, αを用意
+ // Complex[][] C = eVec, _C = new Complex[tc.Length][], λ = eVal, _λ = new Complex[tc.Length][];
+ // Complex[][] exp = new Complex[tc.Length][], _exp = new Complex[tc.Length][];
+ // list.AsParallel().Where(e1 => e1 != null).SelectMany(e2 => e2.SelectMany(e3 => e3.N)).Distinct().ForAll(kIndex =>
+ // {
+ // _C[kIndex] = (new DMat(bLen, bLen, C[kIndex]).ConjugateTranspose() as DMat).Values;
+ // _λ[kIndex] = (new DVec(λ[kIndex]).Conjugate() as DVec).Values;
+
+ // exp[kIndex] = new Complex[tLen * bLen];
+ // _exp[kIndex] = new Complex[tLen * bLen];
+ // for (int t = 0; t < tLen; t++)
+ // for (int j = 0; j < bLen; j++)
+ // {
+ // exp[kIndex][t * bLen + j] = Exp(TwoPiI * (λ[kIndex][j] + kg_z[kIndex][j]) * Thicknesses[t]);
+ // _exp[kIndex][t * bLen + j] = exp[kIndex][t * bLen + j].Conjugate();
+ // }
+ // });
+
+
+ // var total = tcP.Count();
+ // count = 0;
+ // tcP.ForAll(kIndex =>
+ // {
+ // if (bwSTEM.CancellationPending) return;
+
+ // Complex[] λ_k = λ[kIndex], exp_k = exp[kIndex];
+ // double[] kz_k = kg_z[kIndex];
+
+ // Complex[] TDS = Shared.Rent(bLen * bLen), tmpMat = Shared.Rent(bLen * bLen);
+ // try
+ // {
+ // foreach (var (qIndex, n, r, lenz) in list[kIndex])
+ // {
+ // NativeWrapper.MultiplyMxM(bLen, U[qIndex], C[kIndex], ref tmpMat);
+ // var tmpSum = new Complex[tLen];
+ // for (int m = 0; m < n.Length; m++)
+ // {
+ // NativeWrapper.MultiplyMxM(bLen, _C[n[m]], tmpMat, ref TDS);
+
+ // Complex[] λ_kq = _λ[n[m]], exp_kq = _exp[n[m]];
+ // double[] kz_kq = kg_z[n[m]];
+
+ // //B行列の中身を計算し、アダマール積を取る
+ // for (int t = 0; t < tLen; t++)
+ // {
+ // int l = 0;
+ // for (int j = 0; j < bLen; j++)
+ // for (int i = 0; i < bLen; i++)
+ // tmpSum[t] += r[m] * ((exp_k[t * bLen + j] * exp_kq[t * bLen + i] - 1) / TwoPiI / (kz_k[j] - kz_kq[i] + λ_k[j] - λ_kq[i])) * TDS[l++];//B行列は作らず、直接アダマール積を取る
+ // }
+ // }
+ // lock (lockObj2)
+ // for (int t = 0; t < tLen; t++)
+ // for (int d = 0; d < dLen; d++)
+ // I_Inel[qIndex, t, d] += tmpSum[t] / kvac * lenz[d];
+ // }
+ // }
+ // finally { Shared.Return(TDS); Shared.Return(tmpMat); }
+ // if (Interlocked.Increment(ref count) % 10 == 0) bwSTEM.ReportProgress((int)(1000000.0 / total * count), "Calculating I_inelastic(Q)");//状況を報告
+ // });
+
+ //}
#endregion
}
- #endregion
-
if (bwSTEM.CancellationPending) { e.Cancel = true; return; }
+ #endregion
+ #region 各ピクセルの計算
//imagesを初期化
int width = imageSize.Width, height = imageSize.Height;
- STEM_Image = Thicknesses.Select(e => defocusses.Select(e2 => new double[width * height]).ToArray()).ToArray();
+ STEM_Image = Thicknesses.Select(e => defocusses.Select(e2 => GC.AllocateUninitializedArray(width * height)).ToArray()).ToArray();
- #region 各ピクセルの計算
double cX = width / 2.0, cY = height / 2.0;
var shift = (Crystal.RotationMatrix * (Crystal.A_Axis + Crystal.B_Axis + Crystal.C_Axis) / 2).ToPointD;
Parallel.For(0, height, y =>
{
for (int x = 0; x < width; x++)
{
- var rVec = new PointD(-resolution * (x - cX), -resolution * (height - y - 1 - cY)) + shift;
+ var rVec = new PointD(-resolution * (x - cX), resolution * (y - cY)) + shift;//X座標はマイナス。
for (int t = 0; t < Thicknesses.Length; t++)
for (int d = 0; d < dLen; d++)
{
@@ -1319,6 +1346,13 @@ Complex Lenz(in PointD k, in PointD kq, in double defocus)
}
}
});
+
+ //ガウスブラーを適用
+ if (sourceSize > 0)
+ for (int t = 0; t < Thicknesses.Length; t++)
+ for (int d = 0; d < dLen; d++)
+ ImageProcess.GaussianBlurFast(ref STEM_Image[t][d], width, sourceSize / resolution);
+
#endregion
return;
@@ -1445,7 +1479,7 @@ public double[][] GetHRTEMImage(Beam[] beams, Size size, double res, double cs,
Parallel.For(0, divTotal, div =>
{
int start = step * div, count = div == divTotal - 1 ? width * height - start : step;
- var rVec = Enumerable.Range(start, count).SelectMany(n => new[] { -res * (n % width - cX) + shift.X, -res * (height - n / width - 1 - cY) + shift.Y }).ToArray();
+ var rVec = Enumerable.Range(start, count).SelectMany(n => new[] { -res * (n % width - cX) + shift.X, res * (n / width - cY) + shift.Y }).ToArray();//X座標はマイナス。
var results = NativeWrapper.HRTEM_Solver(gPsi, gVec, gLenz, rVec, quasiMode);
for (var i = 0; i < defLen; i++)
Array.Copy(results, i * count, images[i], start, count);
@@ -1455,7 +1489,7 @@ public double[][] GetHRTEMImage(Beam[] beams, Size size, double res, double cs,
{
Parallel.For(0, width * height, n =>
{
- PointD r = new(-(n % width - cX) * res + shift.X, -(height - n / width - 1 - cY) * res + shift.Y), _vec = new(double.NaN, double.NaN);
+ PointD r = new(-(n % width - cX) * res + shift.X, (n / width - cY) * res + shift.Y), _vec = new(double.NaN, double.NaN);//X座標はマイナス。
var sums = new Complex[defLen];
var exp = new Complex(0, 0);
foreach (var (Psi, Vec, Lenz) in gList)
@@ -1472,11 +1506,6 @@ public double[][] GetHRTEMImage(Beam[] beams, Size size, double res, double cs,
images[i][n] = quasiMode ? sums[i].MagnitudeSquared() : Math.Abs(sums[i].Real);
});
}
-
- //20220519 上下左右が反転しているみたいなので、その対処
- for (int i = 0; i < images.Length; i++)
- Array.Reverse(images[i]);
-
return images;
}
#endregion Image Simulation
@@ -1507,49 +1536,36 @@ public double[][] GetHRTEMImage(Beam[] beams, Size size, double res, double cs,
Complex fReal = 0, fImag = 0;
foreach (var atoms in Crystal.Atoms)
{
- //var es = AtomStatic.ElectronScatteringPeng[atoms.AtomicNumber][atoms.SubNumberElectron];//5 gaussian
- var es = AtomStatic.ElectronScatteringEightGaussian[atoms.AtomicNumber];//8 gaussian
+ var es = AtomStatic.ElectronScatteringPeng[atoms.AtomicNumber][atoms.SubNumberElectron];//5 gaussian
+ //var es = AtomStatic.ElectronScatteringEightGaussian[atoms.AtomicNumber];//8 gaussian
var real = es.Factor(s2);//弾性散乱因子
- #region お蔵
- //var m = atoms.Dsf.UseIso || index == (0, 0, 0) ? atoms.Dsf.Biso : 0;
- //if (!atoms.Dsf.UseIso && double.IsNaN(m) && index == (0, 0, 0))// 非等方でg = 0、かつmがNaNの時 Acta Cryst. (1959). 12, 609 , Hamilton の式に従って、Bisoを計算
- // m = (atoms.Dsf.B11 * a * a + atoms.Dsf.B22 * b * b + atoms.Dsf.B33 * c * c + 2 * atoms.Dsf.B12 * a * b + 2 * atoms.Dsf.B23 * b * c + 2 * atoms.Dsf.B31 * c * a) * 4.0 / 3.0;
-
- //if (atoms.Dsf.UseIso || index == (0, 0, 0))
- // imag = (double.IsNaN(inner * outer)) ? es.FactorImaginary(kV, s2, m) :
- // h == null ? es.FactorImaginaryAnnular(kV, g.Vec, m, inner, outer) : es.FactorImaginaryAnnular(kV, g.Vec, h.Vec, m, inner, outer);//非弾性散乱因子 答えは無次元
- #endregion
- double imag = double.NaN, m = double.NaN;
+ var dsf = atoms.Dsf;
+ var zero = dsf.IsZero;
+ double imag = zero ? 0 : double.NaN, m = zero ? 0 : double.NaN;
foreach (var atom in atoms.Atom)
{
- if((!atoms.Dsf.UseIso && index != (0, 0, 0)) || double.IsNaN(imag))//非等方でg≠0の時、あるいは初めての時
+ if (!zero)
{
- if (atoms.Dsf.UseIso)
- m = atoms.Dsf.Biso;
- else if (index == (0, 0, 0))
- m = double.IsNaN(atoms.Dsf.Biso) ? atoms.Dsf.Biso000 : atoms.Dsf.Biso;
- else
+ if ((!dsf.UseIso && index != (0, 0, 0)) || double.IsNaN(imag))//非等方でg≠0の時、あるいは初めての時
{
- var (H, K, L) = atom.Operation.ConvertPlaneIndex(index);
- m = (atoms.Dsf.B11 * H * H + atoms.Dsf.B22 * K * K + atoms.Dsf.B33 * L * L + 2 * atoms.Dsf.B12 * H * K + 2 * atoms.Dsf.B23 * K * L + 2 * atoms.Dsf.B31 * L * H) / s2;
+ if (dsf.UseIso)
+ m = dsf.Biso;
+ else if (index == (0, 0, 0))
+ m = double.IsNaN(dsf.Biso) ? dsf.Biso000 : dsf.Biso;
+ else
+ {
+ var (H, K, L) = atom.Operation.ConvertPlaneIndex(index);
+ m = (dsf.B11 * H * H + dsf.B22 * K * K + dsf.B33 * L * L + 2 * dsf.B12 * H * K + 2 * dsf.B23 * K * L + 2 * dsf.B31 * L * H) / s2;
+ }
+ if (double.IsNaN(m))
+ m = 0;
+
+ imag = m == 0 ? 0 : (double.IsNaN(inner * outer)) ? es.FactorImaginary(kV, s2, m) :
+ h == null ? es.FactorImaginaryAnnular(kV, g.Vec, m, inner, outer) : es.FactorImaginaryAnnular(kV, g.Vec, h.Vec, m, inner, outer);//非弾性散乱因子 答えは無次元
}
- if (double.IsNaN(m))
- m = 0;
- imag = (double.IsNaN(inner * outer)) ? es.FactorImaginary(kV, s2, m) :
- h == null ? es.FactorImaginaryAnnular(kV, g.Vec, m, inner, outer) : es.FactorImaginaryAnnular(kV, g.Vec, h.Vec, m, inner, outer);//非弾性散乱因子 答えは無次元
}
- #region お蔵
- //if (!atoms.Dsf.UseIso && index != (0, 0, 0))
- //{
- // var (H, K, L) = atom.Operation.ConvertPlaneIndex(index);
- // m = (atoms.Dsf.B11 * H * H + atoms.Dsf.B22 * K * K + atoms.Dsf.B33 * L * L + 2 * atoms.Dsf.B12 * H * K + 2 * atoms.Dsf.B23 * K * L + 2 * atoms.Dsf.B31 * L * H) / s2;
- // imag = (double.IsNaN(inner * outer)) ? es.FactorImaginary(kV, s2, m) :
- // h == null ? es.FactorImaginaryAnnular(kV, g.Vec, m, inner, outer) : es.FactorImaginaryAnnular(kV, g.Vec, h.Vec, m, inner, outer);//非弾性散乱因子 答えは無次元
- //}
- //if (double.IsNaN(m)) m = 0;
- #endregion
- var d = Exp(-m * s2 - TwoPiI * (atom * index)) * atoms.Occ;
+ var d = Exp(-m * s2 - TwoPiI * (atom * index)) * atoms.Occ;//expの中がマイナスなのが、U'マトリックスを転置させなければいけない理由かも。違うか。。
fReal += real * d;
fImag += imag * d;
}
@@ -1695,8 +1711,8 @@ public Beam[] Find_gVectors(Matrix3D baseRotation, Vector3DBase vecK0, int maxNu
Vector3DBase g;
while (beams.Count < maxNumOfBloch * 20 && whole.Count < 1000000 && outer.Count > 0)
{
- var min = outer[0].gLen + shift;
- var end = outer.FindLastIndex(o => o.gLen - min < shift * 2);
+ var min = outer[0].gLen + shift;
+ var end = outer.FindLastIndex(o => o.gLen - min < shift * 2);
foreach (var (key, gLen) in CollectionsMarshal.AsSpan(outer)[..(end + 1)])
{
diff --git a/Crystallography/Crystal/ConvertCrystalData.cs b/Crystallography/Crystal/ConvertCrystalData.cs
index 93ce922..04e0181 100644
--- a/Crystallography/Crystal/ConvertCrystalData.cs
+++ b/Crystallography/Crystal/ConvertCrystalData.cs
@@ -3,7 +3,6 @@
using System.ComponentModel;
using System.Data;
using System.Diagnostics;
-using System.Diagnostics.Contracts;
using System.Drawing;
using System.IO;
using System.Linq;
@@ -13,13 +12,11 @@
using System.Threading;
using V3 = OpenTK.Vector3d;
-
-
namespace Crystallography;
public class ConvertCrystalData
{
- static System.StringComparison Ord = System.StringComparison.Ordinal;
+ static readonly System.StringComparison Ord = System.StringComparison.Ordinal;
#region CrystalList(xml`)̓ǂݍ/
public static bool SaveCrystalListXml(Crystal[] crystals, string filename)
@@ -1619,8 +1616,8 @@ private static int SearchSGseriesNumberForCIF(string SgNameHM, string SgNameHall
public static string ConvertToCIF(Crystal crystal)
{
var sb = new StringBuilder();
- sb.AppendLine("# This file is exported from \"" + System.Diagnostics.Process.GetCurrentProcess().ProcessName + "\"");
- sb.AppendLine("# http://pmsl.planet.sci.kobe-u.ac.jp/~seto");
+ sb.AppendLine("# This file is exported from \"" + Process.GetCurrentProcess().ProcessName + "\"");
+ sb.AppendLine("# https://github.com/seto77/");
sb.AppendLine("data_global");
sb.AppendLine("_chemical_name '" + crystal.Name + "'");
@@ -1650,6 +1647,7 @@ public static string ConvertToCIF(Crystal crystal)
sb.AppendLine(";");
#endregion
+ #region iq萔AΏ̐
sb.AppendLine("_chemical_formula_sum '" + crystal.ChemicalFormulaSum + "'");
sb.AppendLine("_cell_length_a " + (crystal.A * 10).ToString("f6"));
sb.AppendLine("_cell_length_b " + (crystal.B * 10).ToString("f6"));
@@ -1668,6 +1666,7 @@ public static string ConvertToCIF(Crystal crystal)
hm = hm.Replace("Rho", "");
sb.AppendLine("_symmetry_space_group_name_H-M '" + hm + "'");
sb.AppendLine("_symmetry_space_group_name_Hall '" + sym.SpaceGroupHallStr + "'");
+ #endregion
#region q̓ʒu
sb.AppendLine("loop_");
@@ -1733,31 +1732,32 @@ public static string ConvertToCIF(Crystal crystal)
foreach (var a in crystal.Atoms)
{
- var u = double.IsNaN(a.Dsf.Uiso) ? 0 : a.Dsf.Uiso;
- sb.AppendLine($"{a.Label} {AtomStatic.AtomicName(a.AtomicNumber)} {a.X:f5} {a.Y:f5} {a.Z:f5} {a.Occ:f5} {u:f5}");
+ var u = double.IsNaN(a.Dsf.Uiso) ? 0 : a.Dsf.Uiso * 100;
+ sb.AppendLine($"{a.Label} {AtomStatic.AtomicName(a.AtomicNumber)} {a.X:f6} {a.Y:f6} {a.Z:f6} {a.Occ:f6} {u:f6}");
}
- sb.AppendLine("loop_");
- sb.AppendLine("_atom_site_aniso_label");
- sb.AppendLine("_atom_site_aniso_U_11");
- sb.AppendLine("_atom_site_aniso_U_22");
- sb.AppendLine("_atom_site_aniso_U_33");
- sb.AppendLine("_atom_site_aniso_U_23");
- sb.AppendLine("_atom_site_aniso_U_13");
- sb.AppendLine("_atom_site_aniso_U_12");
- foreach (var a in crystal.Atoms)
+
+ if (crystal.Atoms.Any(a => !a.Dsf.UseIso))
{
- if (!a.Dsf.UseIso)
+ sb.AppendLine("loop_");
+ sb.AppendLine("_atom_site_aniso_label");
+ sb.AppendLine("_atom_site_aniso_U_11");
+ sb.AppendLine("_atom_site_aniso_U_22");
+ sb.AppendLine("_atom_site_aniso_U_33");
+ sb.AppendLine("_atom_site_aniso_U_23");
+ sb.AppendLine("_atom_site_aniso_U_13");
+ sb.AppendLine("_atom_site_aniso_U_12");
+ foreach (var a in crystal.Atoms.Where(e => !e.Dsf.UseIso))
{
- var u11 = double.IsNaN(a.Dsf.U11) ? 0 : a.Dsf.U11;
- var u22 = double.IsNaN(a.Dsf.U22) ? 0 : a.Dsf.U22;
- var u33 = double.IsNaN(a.Dsf.U33) ? 0 : a.Dsf.U33;
- var u23 = double.IsNaN(a.Dsf.U23) ? 0 : a.Dsf.U23;
- var u31 = double.IsNaN(a.Dsf.U31) ? 0 : a.Dsf.U31;
- var u12 = double.IsNaN(a.Dsf.U12) ? 0 : a.Dsf.U12;
- sb.AppendLine($"{a.Label} {u11:f5} {u22:f5} {u33:f5} {u23:f5} {u31:f5} {u12:f5}");
+ var u11 = double.IsNaN(a.Dsf.U11) ? 0 : a.Dsf.U11 * 100;
+ var u22 = double.IsNaN(a.Dsf.U22) ? 0 : a.Dsf.U22 * 100;
+ var u33 = double.IsNaN(a.Dsf.U33) ? 0 : a.Dsf.U33 * 100;
+ var u23 = double.IsNaN(a.Dsf.U23) ? 0 : a.Dsf.U23 * 100;
+ var u31 = double.IsNaN(a.Dsf.U31) ? 0 : a.Dsf.U31 * 100;
+ var u12 = double.IsNaN(a.Dsf.U12) ? 0 : a.Dsf.U12 * 100;
+ sb.AppendLine($"{a.Label} {u11:f6} {u22:f6} {u33:f6} {u23:f6} {u31:f6} {u12:f6}");
}
- }
+ }
#endregion
return sb.ToString();
diff --git a/Crystallography/Crystal/CrystalMinimum.cs b/Crystallography/Crystal/CrystalMinimum.cs
index 31542e2..7654c00 100644
--- a/Crystallography/Crystal/CrystalMinimum.cs
+++ b/Crystallography/Crystal/CrystalMinimum.cs
@@ -22,22 +22,16 @@ public partial class Crystal2
public float density;
- //[BrotliStringFormatter]
public string name;
- //[BrotliStringFormatter]
public string note;
- //[BrotliStringFormatter]
public string jour;
- //[BrotliStringFormatter]
public string auth;
- //[BrotliStringFormatter]
public string sect;
- //[BrotliStringFormatter]
public string formula;//vZ\ȏꍇ́B
public short sym;
@@ -46,7 +40,6 @@ public partial class Crystal2
public float[] d;//x8ʂ܂łdl
- //[BrotliStringFormatter]
public string fileName;
#endregion
@@ -166,12 +159,8 @@ public static Crystal GetCrystal(Crystal2 c)
var bonds = Bonds.GetVestaBonds(atom.Select(a => a.AtomicNumber));
- var crystal = new Crystal(cell.Values, cell.Errors,
- c.sym, c.name, System.Drawing.Color.FromArgb(c.argb), new Matrix3D(),
- atom.ToArray(), (c.note, c.auth, c.jour, c.sect), bonds);
-
- return crystal;
- }
+ return new Crystal(cell.Values, cell.Errors, c.sym, c.name, System.Drawing.Color.FromArgb(c.argb), new Matrix3D(), atom.ToArray(), (c.note, c.auth, c.jour, c.sect), bonds);
+ }
public static Crystal2 FromCrystal(Crystal c)
{
@@ -228,7 +217,11 @@ public static Crystal2 FromCrystal(Crystal c)
}
[MemoryPackIgnore]
- static readonly char[] toCharDic = new[] { '0', '1', '2', '3', '4', '5', '6', '7', '8', '9', '.', '/', '-', '|', 'E' };
+ private static readonly CultureInfo culture = CultureInfo.InvariantCulture;
+ [MemoryPackIgnore]
+ private static readonly NumberStyles style = NumberStyles.Number;
+ [MemoryPackIgnore]
+ private static readonly StringComparison Ord = StringComparison.Ordinal;
[MemoryPackIgnore]
static readonly string[] toStringDic = new string[]
{
@@ -539,14 +532,14 @@ public static byte[] ToBytes(string s)
return new[] { (byte)(240 + 0) };
else
{
- if (s.StartsWith("0.", StringComparison.Ordinal))
+ if (s.StartsWith("0.", Ord))
{
if (s == "0.")
s = "0";
else
s = s[1..];
}
- else if (s.StartsWith("-0.", StringComparison.Ordinal))
+ else if (s.StartsWith("-0.", Ord))
s = s.Replace("-0.", "-.");
try
{
@@ -583,18 +576,9 @@ public static string ToString(byte[] bytes)
return sb.ToString();
}
-
-
-
-
private static (double Value, double Error) Decompose(string str) => Decompose(str, false);
public static (double Value, double Error) Decompose(string str, int sgnum) => Decompose(str, sgnum >= 430 && sgnum <= 488);
- [MemoryPackIgnore]
- private static readonly CultureInfo culture = CultureInfo.InvariantCulture;
- [MemoryPackIgnore]
- private static readonly NumberStyles style = NumberStyles.Number;
- [MemoryPackIgnore]
- private static readonly StringComparison Ord = StringComparison.Ordinal;
+
///
/// 9.726|5|, 1.234|12|E-6 ̂悤ȕAValueErrorɕă^vŕԂ.
diff --git a/Crystallography/Crystallography.csproj b/Crystallography/Crystallography.csproj
index 77da986..0609116 100644
--- a/Crystallography/Crystallography.csproj
+++ b/Crystallography/Crystallography.csproj
@@ -4,8 +4,8 @@
Library
net7.0-windows
true
- 2023.3.8.0528
- 2023.3.8.0528
+ 2023.3.18.0710
+ 2023.3.18.0710
diff --git a/Crystallography/Images/ImageProcessor.cs b/Crystallography/Images/ImageProcessor.cs
index 83a6777..016b462 100644
--- a/Crystallography/Images/ImageProcessor.cs
+++ b/Crystallography/Images/ImageProcessor.cs
@@ -8,6 +8,7 @@
namespace Crystallography;
public static class ImageProcess
{
+ #region Gaussian Blur
///
/// GaussianBlurを施す。
///
@@ -72,6 +73,31 @@ unsafe static public double[] GaussianBlur(double[] pixels, int width, double ra
return blurTemp[0];
}
+ ///
+ /// GaussianBlurを施す。横方向と縦方向を分離するため、高速。
+ ///
+ /// 画像データの一次元配列
+ /// 画像の幅
+ /// ピクセル単位でのフィルムにじみ半値半幅
+ ///
+ static public double[] GaussianBlurFast(double[] pixels, int width, double hwhm)
+ {
+ var results = GC.AllocateUninitializedArray(pixels.Length);
+ GaussianBlurFast(pixels, width, hwhm, results);
+ return results;
+ }
+
+ ///
+ /// GaussianBlurを施す。横方向と縦方向を分離するため、高速。pixelsが直接書き換えられる。
+ ///
+ /// 画像データの一次元配列
+ /// 画像の幅
+ /// ピクセル単位でのフィルムにじみ半値半幅
+ static public void GaussianBlurFast(ref double[] pixels, int width, double hwhm)
+ {
+ GaussianBlurFast(pixels, width, hwhm, pixels);
+ }
+
///
/// GaussianBlurを施す。横方向と縦方向を分離するため、高速。
///
@@ -79,21 +105,22 @@ unsafe static public double[] GaussianBlur(double[] pixels, int width, double ra
/// 画像の幅
/// ピクセル単位でのフィルムにじみ半値半幅
///
- unsafe static public double[] GaussianBlurFast(double[] pixels, int width, double hwhm)
+ static private void GaussianBlurFast(double[] pixels, int width, double hwhm, double[] results)
{
+ int height = pixels.Length / width;
int limit = (int)(hwhm * 3) * 2 + 1;
- int height = pixels.Length / width;
int center = limit / 2;
-
- var results = new double[width * height];
+
if (limit == 1)
{
- Array.Copy(pixels, results, pixels.Length);
- return results;
+ if (pixels != results)
+ Array.Copy(pixels, results, pixels.Length);
+ return;
}
- var tmpPixels = ArrayPool.Shared.Rent(width * height);
+ double[] tmpPixels = ArrayPool.Shared.Rent(width * height);
+ double[] blurSumH = GC.AllocateUninitializedArray(height), blurSumW = GC.AllocateUninitializedArray(width);
try
{
var blur = new double[limit];
@@ -101,28 +128,49 @@ unsafe static public double[] GaussianBlurFast(double[] pixels, int width, doubl
blur[h] = Math.Exp(-(h - center) * (h - center) / hwhm / hwhm * Math.Log(2));
blur = Statistics.Normarize(blur);
+ Parallel.For(0, height, h =>
+ {
+ if (h < center)
+ blurSumH[h] = blur[(center - h)..].Sum();
+ else if (h >= height - center)
+ blurSumH[h] = blur[..(height - h + center)].Sum();
+ else
+ blurSumH[h] = 1;
+ });
Parallel.For(0, width, w =>
{
for (int h = 0; h < height; h++)
{
tmpPixels[h * width + w] = 0;
for (int n = Math.Max(0, center - h); n < Math.Min(blur.Length, height - h + center); n++)
- tmpPixels[h * width + w] += blur[n] * pixels[(h - center + n) * width + w];
+ tmpPixels[h * width + w] += blur[n] / blurSumH[h] * pixels[(h - center + n) * width + w];
}
});
+ Parallel.For(0, width, w =>
+ {
+ if (w < center)
+ blurSumW[w] = blur[(center - w)..].Sum();
+ else if (w >= width - center)
+ blurSumW[w] = blur[..(width - w + center)].Sum();
+ else
+ blurSumW[w] = 1;
+ });
Parallel.For(0, height, h =>
{
for (int w = 0; w < width; w++)
+ {
+ results[h * width + w] = 0;
for (int n = Math.Max(0, center - w); n < Math.Min(blur.Length, width - w + center); n++)
- results[h * width + w] += blur[n] * tmpPixels[h * width + w - center + n];
+ results[h * width + w] += blur[n] / blurSumW[w] * tmpPixels[h * width + w - center + n];
+ }
});
}
finally { ArrayPool.Shared.Return(tmpPixels); }
-
- return results;
}
+ #endregion
+
///
/// 周囲のピクセルと比べて、標準偏差 × threshold以上外れたピクセルは、周囲のピクセルの平均強度にする。
///
diff --git a/Crystallography/NativeWrapper.cs b/Crystallography/NativeWrapper.cs
index 9d5735a..ee428b0 100644
--- a/Crystallography/NativeWrapper.cs
+++ b/Crystallography/NativeWrapper.cs
@@ -1,12 +1,13 @@
-using MathNet.Numerics.LinearAlgebra.Complex;
+#region using, namespace
+using MathNet.Numerics.LinearAlgebra.Complex;
using System;
using System.Buffers;
using System.Collections.Generic;
using System.Linq;
using System.Numerics;
using System.Runtime.InteropServices;
-
namespace Crystallography;
+#endregion
public static partial class NativeWrapper
{
@@ -145,10 +146,15 @@ private static unsafe partial void _Histogram(
);
[LibraryImport("Crystallography.Native.dll")]
- private static unsafe partial void _GenerateTC(int dim, double thickness, double* _kg_z, double* _val, double* _vec, double* _result);
+ private static unsafe partial void _GenerateTC1(int dim, double thickness, double* _kg_z, double* _val, double* _vec, double* _tc_k);
+
+ [LibraryImport("Crystallography.Native.dll")]
+ private static unsafe partial void _GenerateTC2(int dim, double thickness, double* _kg_z, double* _val, double* _vec, double* _tc_k, double* _tc_kq);
[LibraryImport("Crystallography.Native.dll")]
private static unsafe partial void _RowVec_SqMat_ColVec(int dim, double* _rowVec, double* _sqMat, double* _colVec, double* _result);
+ [LibraryImport("Crystallography.Native.dll")]
+ private static unsafe partial void _STEM_INEL1(int dim, double* rowVec, int* n, double* r, double* sqMat, double* colVec, double* _result);
#endregion
@@ -383,32 +389,32 @@ unsafe static public void Divide(int dim, in Complex[] v1, in Complex[] v2, ref
#region Blend関数
unsafe static public void Blend(in int dim, in Complex[] c0, in Complex[] c1, in Complex[] c2, in Complex[] c3, in double r0, in double r1, in double r2, in double r3, ref Complex[] result)
{
- fixed (Complex* p0 = c0)
- fixed (Complex* p1 = c1)
- fixed (Complex* p2 = c2)
- fixed (Complex* p3 = c3)
- fixed (Complex* res = result)
+ fixed (Complex* p0 = c0, p1 = c1, p2 = c2, p3 = c3, res = result)
_Blend(dim * 2, (double*)p0, (double*)p1, (double*)p2, (double*)p3, r0, r1, r2, r3, (double*)res);
}
+
+ unsafe static public void Blend(in int dim, Complex* c0, Complex* c1, Complex* c2, Complex* c3, in double r0, in double r1, in double r2, in double r3, Complex* res)
+ {
+ //fixed (Complex* p0 = c0, p1 = c1, p2 = c2, p3 = c3, res = result)
+ _Blend(dim * 2, (double*)c0, (double*)c1, (double*)c2, (double*)c3, r0, r1, r2, r3, (double*)res);
+ }
+
unsafe static public void Blend(in int dim, in double[] c0, in double[] c1, in double[] c2, in double[] c3, in double r0, in double r1, in double r2, in double r3, ref double[] result)
{
- fixed (double* p0 = c0)
- fixed (double* p1 = c1)
- fixed (double* p2 = c2)
- fixed (double* p3 = c3)
- fixed (double* res = result)
+ fixed (double* p0 = c0, p1 = c1, p2 = c2, p3 = c3, res = result)
_Blend(dim, p0, p1, p2, p3, r0, r1, r2, r3, res);
}
unsafe static public void BlendAndConjugate(in int dim, in Complex[] c0, in Complex[] c1, in Complex[] c2, in Complex[] c3, in double r0, in double r1, in double r2, in double r3, ref Complex[] result)
{
- fixed (Complex* p0 = c0)
- fixed (Complex* p1 = c1)
- fixed (Complex* p2 = c2)
- fixed (Complex* p3 = c3)
- fixed (Complex* res = result)
+ fixed (Complex* p0 = c0, p1 = c1, p2 = c2, p3 = c3, res = result)
_BlendAndConjugate(dim, (double*)p0, (double*)p1, (double*)p2, (double*)p3, r0, r1, r2, r3, (double*)res);
}
+
+ unsafe static public void BlendAndConjugate(in int dim, Complex* c0, in Complex* c1, in Complex* c2, in Complex* c3, in double r0, in double r1, in double r2, in double r3, Complex* res)
+ {
+ _BlendAndConjugate(dim, (double*)c0, (double*)c1, (double*)c2, (double*)c3, r0, r1, r2, r3, (double*)res);
+ }
#endregion
#region Eigenライブラリーを利用して、PartialPivLuSolveを求める
@@ -469,14 +475,6 @@ unsafe static public void AdjointAndMultiply(int dim, Complex[] matrix1, Complex
#region 複素共役、転置
- unsafe static public void Conjugate(int dim, Complex[] matrix1, Complex[] matrix2, ref Complex[] result)
- {
- fixed (Complex* mtx1 = matrix1)
- fixed (Complex* mtx2 = matrix2)
- fixed (Complex* res = result)
- _AdjointAndMultiply(dim, (double*)mtx1, (double*)mtx2, (double*)res);
- }
-
unsafe static public void Adjoint(int dim, Complex[] matrix1, Complex[] matrix2, ref Complex[] result)
{
fixed (Complex* mtx1 = matrix1)
@@ -578,32 +576,20 @@ static unsafe public void MatrixExponential(in int dim, Complex[] mat, ref Compl
#region STEMの非弾性散乱電子強度の計算用の特殊関数
unsafe static public void AdjointMul_Mul_Mul(in int dim, in Complex[] mat1, in Complex[] mat2, in Complex[] mat3, ref Complex[] result)
{
- fixed (Complex* _mat1 = mat1)
- fixed (Complex* _mat2 = mat2)
- fixed (Complex* _mat3 = mat3)
- fixed (Complex* res = result)
+ fixed (Complex* _mat1 = mat1, _mat2 = mat2, _mat3 = mat3, res = result)
_AdJointMul_Mul_Mul(dim, (double*)_mat1, (double*)_mat2, (double*)_mat3, (double*)res);
}
unsafe static public void BlendAdjointMul_Mul_Mul(in int dim, in Complex[] c0, in Complex[] c1, in Complex[] c2, in Complex[] c3, double r0, double r1, double r2, double r3,
in Complex[] mat2, in Complex[] mat3, ref Complex[] result)
{
- fixed (Complex* p0 = c0)
- fixed (Complex* p1 = c1)
- fixed (Complex* p2 = c2)
- fixed (Complex* p3 = c3)
- fixed (Complex* _mat2 = mat2)
- fixed (Complex* _mat3 = mat3)
- fixed (Complex* res = result)
+ fixed (Complex* p0 = c0, p1 = c1, p2 = c2, p3 = c3, _mat2 = mat2,_mat3 = mat3, res = result)
_BlendAdJointMul_Mul_Mul(dim, (double*)p0, (double*)p1, (double*)p2, (double*)p3, r0, r1, r2, r3, (double*)_mat2, (double*)_mat3, (double*)res);
}
unsafe static public void TDS(in int dim, in Complex[] mat1, in Complex[] mat2, in Complex[] mat3, ref Complex[] result)
{
- fixed (Complex* _mat1 = mat1)
- fixed (Complex* _mat2 = mat2)
- fixed (Complex* _mat3 = mat3)
- fixed (Complex* res = result)
+ fixed (Complex* _mat1 = mat1, _mat2 = mat2, _mat3 = mat3, res = result)
_AdJointMul_Mul_Mul(dim, (double*)_mat1, (double*)_mat2, (double*)_mat3, (double*)res);
}
@@ -616,11 +602,21 @@ unsafe static public void TDS(in int dim, in Complex[] mat1, in Complex[] mat2,
///
///
///
- unsafe static public void GenerateTC(in int dim,in double thickness, in double[] kg_z, in Complex[] val, in Complex[] vec, ref Complex[] result)
+ //unsafe static public void GenerateTC(in int dim,in double thickness, in double[] kg_z, in Complex[] val, in Complex[] vec, ref Complex[] result)
+ //{
+ // fixed (double* _kg_z = kg_z)
+ // fixed (Complex* _val = val, _vec = vec, _result = result)
+ // _GenerateTC(dim, thickness, _kg_z, (double*)_val, (double*)_vec, (double*)_result);
+ //}
+
+ unsafe static public void GenerateTC1(in int dim, in double thickness, double* _kg_z, Complex* _val, Complex* _vec, Complex* _tc_k)
{
- fixed (double* _kg_z = kg_z)
- fixed (Complex* _val = val, _vec = vec, _result = result)
- _GenerateTC(dim, thickness, _kg_z, (double*)_val, (double*)_vec, (double*)_result);
+ _GenerateTC1(dim, thickness, _kg_z, (double*)_val, (double*)_vec, (double*)_tc_k);
+ }
+
+ unsafe static public void GenerateTC2(in int dim, in double thickness, double* _kg_z, Complex* _val, Complex* _vec, Complex* _tc_k, Complex* _tc_kq)
+ {
+ _GenerateTC2(dim, thickness, _kg_z, (double*)_val, (double*)_vec, (double*)_tc_k, (double*)_tc_kq);
}
///
@@ -633,16 +629,31 @@ unsafe static public void GenerateTC(in int dim,in double thickness, in double[]
///
unsafe static public Complex RowVec_SqMat_ColVec(in int dim, Complex[] rowVec, Complex[] sqMtx, Complex[] colVec)
{
- var result = new double[2];
+ var result = new Complex();
fixed (Complex* _rowVec = rowVec, _sqMtx = sqMtx, _colVec = colVec)
- fixed (double* _res = result)
- _RowVec_SqMat_ColVec(dim, (double*)_rowVec, (double*)_sqMtx, (double*)_colVec, _res);
+ _RowVec_SqMat_ColVec(dim, (double*)_rowVec, (double*)_sqMtx, (double*)_colVec, (double*)&result);
- return new Complex(result[0], result[1]);
+ return result;
+ }
+
+ unsafe static public Complex RowVec_SqMat_ColVec(in int dim, Complex* _rowVec, Complex* _sqMtx, Complex* _colVec)
+ {
+ var result = new Complex();
+ _RowVec_SqMat_ColVec(dim, (double*)_rowVec, (double*)_sqMtx, (double*)_colVec, (double*)&result);
+ return result;
+ }
+
+ unsafe static public Complex STEM_INEL1(in int dim, Complex* _rowVec, int[] n, double[] r, Complex* _sqMtx, Complex* _colVec)
+ {
+ var result = new Complex();
+ fixed (int* _n = n)
+ fixed (double* _r = r)
+ _STEM_INEL1(dim, (double*)_rowVec, _n, _r, (double*)_sqMtx, (double*)_colVec, (double*)&result);
+ return result;
}
#endregion
-
+
#region CBED
///
/// Eigenライブラリーを利用して固有値解を求めて、CBEDの解を求める
diff --git a/Crystallography/Profile.cs b/Crystallography/Profile.cs
index dbe4312..6578373 100644
--- a/Crystallography/Profile.cs
+++ b/Crystallography/Profile.cs
@@ -637,7 +637,7 @@ public object Clone()
public int ImageWidth = 0;
public int ImageHeight = 0;
- public override string ToString() => Name.ToString();
+ public override string ToString() => Name;
public DiffractionProfile()
{
diff --git a/IPAnalyzer/FormMain.cs b/IPAnalyzer/FormMain.cs
index be8d018..07d9f15 100644
--- a/IPAnalyzer/FormMain.cs
+++ b/IPAnalyzer/FormMain.cs
@@ -3363,9 +3363,9 @@ public void GetProfile(string fileName = "")
};
var targets = new int[] { -1 };
- //シーケンシャルイメージモードの時の処理
+
if (toolStripButtonImageSequence.Enabled)
- {
+ {//シーケンシャルイメージモードの時の処理
if (toolStripMenuItemAllSequentialImages.Checked)
targets = Enumerable.Range(0, FormSequentialImage.MaximumNumber).ToArray();
else if (toolStripMenuItemSelectedSequentialImages.Checked)
@@ -3374,17 +3374,20 @@ public void GetProfile(string fileName = "")
for (int i = 0; i < targets.Length; i++)
{
- //シーケンシャルイメージモードの時の処理
+
if (targets[0] != -1)
- {
+ {//シーケンシャルイメージモードの時の処理
FormSequentialImage.SkipCalcFreq = i != targets.Length - 1;
FormSequentialImage.AverageMode = false;
FormSequentialImage.SelectedIndex = targets[i];
if (Ring.ImageType == Ring.ImageTypeEnum.HDF5)
diffractionProfile.SrcWaveLength = UniversalConstants.Convert.EnergyToXrayWaveLength(Ring.SequentialImageEnergy[targets[i]]);
}
+
if (toolStripButtonImageSequence.Enabled)
diffractionProfile.Name = FileName + " " + FileNameSub;
+ else
+ diffractionProfile.Name = FileName;
//プロファイルを作成
diffractionProfile.OriginalProfile = Ring.GetProfile(IP);
diff --git a/IPAnalyzer/IPAnalyzer.csproj b/IPAnalyzer/IPAnalyzer.csproj
index a8df945..4405011 100644
--- a/IPAnalyzer/IPAnalyzer.csproj
+++ b/IPAnalyzer/IPAnalyzer.csproj
@@ -4,8 +4,8 @@
WinExe
net7.0-windows
true
- 2023.3.8.0528
- 2023.3.8.0528
+ 2023.3.18.0710
+ 2023.3.18.0710
App.ico
IPAnalyzer.Program
embedded
diff --git a/IPAnalyzer/Version.cs b/IPAnalyzer/Version.cs
index d5640ba..cbea1bb 100644
--- a/IPAnalyzer/Version.cs
+++ b/IPAnalyzer/Version.cs
@@ -9,6 +9,7 @@ static class Version
static public string History =
"History" +
+ "\r\n ver3.945(2023/03/18) Fixed a minor bug." +
"\r\n ver3.944(2023/03/07) Fixed a minor bug." +
"\r\n ver3.943(2023/03/07) Fixed bugs on the 'Macro function'." +
"\r\n ver3.942(2023/02/06) Improved 'Auto Procedure'" +