rotate matrix from 2 vectors

本文探讨了从两个已归一化的向量计算旋转矩阵的过程。通过使用点积计算角度及叉积确定轴,进而创建旋转矩阵。作者寻求更简洁高效的方法。

综合自己想的和能找到的方法大致如下,但是还是觉得比较繁琐:

 

如果见到有快速优雅的方法,请麻烦留个言告知下,谢谢。

import numpy as np from scipy.interpolate import griddata import os from tqdm import tqdm import glob import math # 配置参数 PHYS_SIZE = 2.4 GRID_RES = 256 CELL_SIZE = PHYS_SIZE / GRID_RES INPUT_DIR = './ASCII' # 45度数据目录 OUTPUT_DIR = './probe_data' # 对齐后的输出目录 # 创建输出目录 os.makedirs(OUTPUT_DIR, exist_ok=True) # 旋转矩阵常数 ROT_MATRIX = np.array([ [math.cos(math.radians(45)), -math.sin(math.radians(45))], [math.sin(math.radians(45)), math.cos(math.radians(45))] ]) # 风速旋转矩阵(逆旋转) WIND_ROT_MATRIX = np.array([ [math.cos(math.radians(-45)), -math.sin(math.radians(-45))], [math.sin(math.radians(-45)), math.cos(math.radians(-45))] ]) def rotate_coordinates(data): """将45度坐标旋转回0度坐标系""" # 提取坐标并旋转 coords = data[:, :2] rotated_coords = np.dot(coords, ROT_MATRIX.T) # 更新坐标数据 data[:, 0] = rotated_coords[:, 0] data[:, 1] = rotated_coords[:, 1] return data def rotate_wind_vectors(data): """将45度风速矢量旋转回0度坐标系""" # 提取风速分量 u = data[:, 3] v = data[:, 4] wind_vectors = np.column_stack((u, v)) # 旋转风速矢量 rotated_wind = np.dot(wind_vectors, WIND_ROT_MATRIX.T) # 更新风速数据 data[:, 3] = rotated_wind[:, 0] data[:, 4] = rotated_wind[:, 1] return data def create_centered_grid(): """创建中心对齐的网格坐标系""" start = -PHYS_SIZE / 2 + CELL_SIZE / 2 end = PHYS_SIZE / 2 - CELL_SIZE / 2 x = np.linspace(start, end, GRID_RES) y = np.linspace(start, end, GRID_RES) return np.meshgrid(x, y, indexing='xy') def generate_probe_matrix(raw_data, grid_x, grid_y): """为45度数据生成对齐的探针矩阵""" # 1. 旋转坐标和风速到0度基准 rotated_data = rotate_coordinates(raw_data.copy()) rotated_data = rotate_wind_vectors(rotated_data) # 2. 筛选在场域内的点 mask = (np.abs(rotated_data[:, 0]) <= PHYS_SIZE / 2) & \ (np.abs(rotated_data[:, 1]) <= PHYS_SIZE / 2) data = rotated_data[mask] if data.size == 0: return np.zeros((GRID_RES, GRID_RES, 3), dtype=np.float32) # 3. 准备插值 points = data[:, :2] values = data[:, [3, 4, 2]] # [Ux, Uy, P] # 4. 线性插值 interp = griddata(points, values, (grid_x, grid_y), method='linear') # 5. 处理NaN值 nan_mask = np.isnan(interp) if np.any(nan_mask): nearest = griddata(points, values, (grid_x, grid_y), method='nearest') interp[nan_mask] = nearest[nan_mask] return interp.astype(np.float32) # 主处理函数(更新版本) def process_45deg_cases(): """处理45度风向角数据并对齐到0度坐标系""" # 预计算网格坐标 grid_x, grid_y = create_centered_grid() # 获取所有输入文件 input_files = sorted(glob.glob(os.path.join(INPUT_DIR, "caseE45n-*"))) if not input_files: raise FileNotFoundError(f"在 {INPUT_DIR} 中未找到任何匹配 caseE45n-* 的文件") # 创建进度条 progress = tqdm(total=len(input_files), desc="处理45度案例") for i, input_path in enumerate(input_files): try: # 提取文件编号 file_num = os.path.basename(input_path).split("-")[1] # 设置输出路径 output_path = os.path.join(OUTPUT_DIR, f"probe_45_{file_num}.npz") # 读取原始数据 raw_data = np.loadtxt( input_path, dtype=np.float32, delimiter=None, skiprows=1, usecols=(1, 2, 4, 5, 6) # [x, y, P, Ux, Uy] ) # 生成对齐的探针矩阵 probe_matrix = generate_probe_matrix(raw_data, grid_x, grid_y) # 保存结果 np.savez_compressed( output_path, x=grid_x, y=grid_y, fields=probe_matrix ) # 更新进度 progress.update(1) except Exception as e: print(f"\n处理文件 {input_path} 时出错: {e}") continue progress.close() print(f"\n对齐处理完成! 共生成 {len(input_files)} 个对齐的探针矩阵到 {OUTPUT_DIR}") # 执行处理 if __name__ == "__main__": process_45deg_cases() 解析上述代码
最新发布
10-20
using OpenTK.GLControl; using OpenTK.Graphics.OpenGL4; using OpenTK.Mathematics; using Sof.ClientApplication.Framework.ViewModels; using System.Runtime.InteropServices; using Waferscan.GUI.Result.OpenGL; using Waferscan.GUI.Result.View; using Waferscan.GUI.Result.ViewModel; using Waferscan.Logics.Algorithm.Model; using Waferscan.Recipe.Enums; using Waferscan.Result.Models; namespace Waferscan.GUI.Result.WaferMap { public class Haze3dDraw : BaseMapDraw { public Channel Channel => _wafermap.Channel; public float MinThreshold { get; set; } public float MaxThreshold { get; set; } #region const const float SHOW_RATE = 0.999f; //显示区域缩小,防止部分像素点落在视野外 #endregion #region filed private Shader _shader; private int _hazeBufferObject; private int _hazeArrayObject; private int _meshBufferObject; private float[] _hazeDatas; private MeshData _meshData; private float[] _colors; #endregion public override void Init(ResultDisplayViewModel resultDisplay,WaferMapViewModel waferMap) { base.Init(resultDisplay, waferMap); ((Camera3D)_camera).UpdateViewport(_glControl.Width, _glControl.Height); _glControl.MakeCurrent(); _shader = new Shader("Shaders/shader_haze3d.vert", "Shaders/shader_haze3d.frag"); _hazeBufferObject = GL.GenBuffer(); GL.BindBuffer(BufferTarget.ArrayBuffer, _hazeBufferObject); _hazeArrayObject = GL.GenVertexArray(); GL.BindVertexArray(_hazeArrayObject); GL.VertexAttribPointer(0, 3, VertexAttribPointerType.Float, false, 3 * sizeof(float), 0); GL.EnableVertexAttribArray(0); //mesh _meshBufferObject = GL.GenBuffer(); GL.BindBuffer(BufferTarget.ElementArrayBuffer, _meshBufferObject); } protected override ICamera CreateCamera() { IsCtrlDown= true; var camera= new Camera3D(); return camera; } public override void KeyUp(object sender, System.Windows.Input.KeyEventArgs e) { } public override void Release() { _glControl.MakeCurrent(); // Unbind all the resources by binding the targets to 0/null. GL.BindBuffer(BufferTarget.ArrayBuffer, 0); GL.BindVertexArray(0); GL.UseProgram(0); // Delete all the resources. GL.DeleteBuffer(_hazeBufferObject); GL.DeleteBuffer(_meshBufferObject); GL.DeleteVertexArray(_hazeArrayObject); GL.DeleteProgram(_shader.Handle); } public override void Draw(ResultDisplayModel result) { _glControl.MakeCurrent(); // var hazeDatas = result.HazeDatas?[(int)Channel - 1]; if (hazeDatas == null || hazeDatas.Length < (int)Channel) { _hazeDatas = new float[0]; _glControl.Invalidate(); return; } //着色 var colorBar = _wafermap.ColorBar.BinsOrHazeItems; _colors = new float[colorBar.Count * 3]; for (int i = 0; i < colorBar.Count; i++) { _colors[i * 3] = (float)colorBar[i].R / 255; _colors[i * 3 + 1] = (float)colorBar[i].G / 255; _colors[i * 3 + 2] = (float)colorBar[i].B / 255; } _shader.SetVector3("colors[0]", _colors); // _hazeDatas = Normalize((int)result.WaferSize, hazeDatas); _meshData = PointCloud.BuildMesh(_hazeDatas, _hazeDatas.Length); //数据绑定 GL.BindBuffer(BufferTarget.ArrayBuffer, _hazeBufferObject); GL.BufferData(BufferTarget.ArrayBuffer, _hazeDatas.Length * sizeof(float), _hazeDatas, BufferUsageHint.StaticDraw); GL.BindBuffer(BufferTarget.ElementArrayBuffer, _meshBufferObject); GL.BufferData(BufferTarget.ElementArrayBuffer, _meshData.size * sizeof(uint), _meshData.ptr, BufferUsageHint.StaticDraw); //重绘 _glControl.Invalidate(); } public override void ReDraw() { } public float[] Normalize(int wafersize, HazeData[] hazeDatas) { var rSize = wafersize / 2 * 1000; MinThreshold = (float)_wafermap.ColorBar.HazeMinThreshold; MaxThreshold = (float)_wafermap.ColorBar.HazeMaxThreshold; if (MinThreshold == 0 && MaxThreshold == 0) { var haze2 = hazeDatas.OrderBy(x => x.Z).ToList(); MinThreshold = haze2[haze2.Count / 10].Z; MaxThreshold = haze2[haze2.Count / 10 * 9].Z; } var _hazeDatas = new float[hazeDatas.Length * 3]; for (int i = 0; i < hazeDatas.Length; i++) { _hazeDatas[i * 3] = (hazeDatas[i].X / rSize - 1) * SHOW_RATE; _hazeDatas[i * 3 + 1] = (hazeDatas[i].Y / rSize - 1) * SHOW_RATE; var signal = hazeDatas[i].Z; if (signal <= MinThreshold) _hazeDatas[i * 3 + 2] = 0; else if (signal >= MaxThreshold) _hazeDatas[i * 3 + 2] = 0.5f; else _hazeDatas[i * 3 + 2] = (signal - MinThreshold) / (MaxThreshold - MinThreshold) *0.5f; } return _hazeDatas; } private Matrix4 _modelMatrix = Matrix4.Identity; public override void Paint(object? sender, PaintEventArgs e) { lock (this) { _glControl.MakeCurrent(); } GL.Clear(ClearBufferMask.ColorBufferBit); if (_hazeDatas == null || _hazeDatas.Length == 0) return; _shader.Use(); _shader.SetMatrix4("haze3dModel",_modelMatrix); _shader.SetMatrix4("haze3dView", _camera.GetViewMatrix()); _shader.SetMatrix4("haze3dProjection",((Camera3D) _camera).GetProjectionMatrix()); GL.BindVertexArray(_hazeArrayObject); GL.DrawElements(PrimitiveType.Triangles, _meshData.size, DrawElementsType.UnsignedInt, 0); _glControl.SwapBuffers(); } } class PointCloud { [DllImport("Waferscan.Result.PointCloud.dll", CharSet = CharSet.Unicode, SetLastError = true)] public static extern MeshData BuildMesh(float[] datas, int size); } public struct MeshData { public IntPtr ptr; public int size; } } 这是Draw的代码,#version 330 core layout(location = 0) in vec3 aPosition; uniform mat4 haze3dView; uniform mat4 haze3dModel; uniform mat4 haze3dProjection; out vec3 vPos; void main() { gl_Position = haze3dProjection * haze3dView * haze3dModel * vec4(aPosition, 1.0); vPos = aPosition; }这是.vert的代码,using OpenTK.Graphics.OpenGL4; using OpenTK; using OpenTK.Mathematics; namespace Waferscan.GUI.Result.OpenGL { public class Camera3D : ICamera { // Default camera values const float YAW = -90.0f; const float PITCH = 0.0f; const float SPEED = 2.5f; const float SENSITIVITY = 0.1f; const float DEFAULT_ZOOM = 45.0f; // camera Attributes public Vector3 Position { get; set; }=new Vector3(0.0f,0.0f,5.0f); Vector3 Front; Vector3 Up; Vector3 Right; Vector3 WorldUp=new Vector3(0.0f,1.0f,0.0f); // euler Angles private float Yaw=YAW; private float Pitch=PITCH; private float _viewportWidth = 800; private float _viewportHeight = 600; public void UpdateViewport(int Height, int Width) { _viewportHeight = Height; _viewportWidth = Width; } // camera options float MovementSpeed= SPEED; float MouseSensitivity = SENSITIVITY; public float _fieldOfView = MathHelper.DegreesToRadians(DEFAULT_ZOOM); //Projection parameters public float FieldOfView { get=>_fieldOfView; set=>_fieldOfView=value; } public float AspectRatio { get; set; } = 16.0f / 9.0f; public float NearPlane { get; set; } = 0.1f; public float FarPlane { get; set; } = 1000.0f; //Model matrix public Matrix4 ModelMatrix { get; private set; } = Matrix4.Identity; //Fog effect parameters public float FogStart { get; set; } = 0.0f; public float FogEnd { get; set; } = 10.0f; public float Scale { get; set; } = 1.0f; public Camera3D() { UpdateModelMatrix(Vector3.Zero, Vector3.Zero, Vector3.One); updateCameraVectors(); } // returns the view matrix calculated using Euler Angles and the LookAt Matrix public Matrix4 GetViewMatrix() { return Matrix4.LookAt(Position, Position + Front, Up); } //returns the projection matrix calculated using the Euler angle and the LookAt matrix public Matrix4 GetProjectionMatrix() { float fov = _fieldOfView; float aspectRatio = _viewportWidth / _viewportHeight; float near = NearPlane; float far = FarPlane; return Matrix4.CreatePerspectiveFieldOfView(fov, aspectRatio, near, far ); } /// <summary> /// Update the model fov, aspectRatio, near, far) (position, rotation, scale) and save the current state /// </summary> public void UpdateModelMatrix(Vector3 position, Vector3 rotation, Vector3 scale) { //Matrix multiplication order: first scale, then rotate, and finally translate ModelMatrix = Matrix4.Identity; ModelMatrix *= Matrix4.CreateScale(scale); ModelMatrix *= Matrix4.CreateRotationX(MathHelper.DegreesToRadians(rotation.X)); ModelMatrix *= Matrix4.CreateRotationY(MathHelper.DegreesToRadians(rotation.Y)); ModelMatrix *= Matrix4.CreateRotationZ(MathHelper.DegreesToRadians(rotation.Z)); ModelMatrix *= Matrix4.CreateTranslation(position); } //Pass parameters to the shader public void SetShaderParameters(int shaderProgram) { if (shaderProgram == 0) { Console.WriteLine("Error: Shader program is not initialized!"); return; } //Pass model matrix var tempModel = ModelMatrix; int modelLoc = GL.GetUniformLocation(shaderProgram, "haze3dModel"); Console.WriteLine($"modelLoc: {modelLoc}"); Console.WriteLine($"ModelMatrix: {ModelMatrix}"); if (modelLoc != -1) GL.UniformMatrix4(modelLoc, false, ref tempModel); else Console.WriteLine("Warning: Uniform 'haze3dModel' not found in shader!"); //Pass view matrix Matrix4 viewMatrix =GetViewMatrix(); int viewLoc = GL.GetUniformLocation(shaderProgram, "haze3dView"); if (viewLoc != -1) GL.UniformMatrix4(viewLoc, false, ref viewMatrix); else Console.WriteLine("Warning: Uniform 'haze3dView' not found in shader!"); //Pass projection matrix Matrix4 projectionMatrix =GetProjectionMatrix(); int projectionLoc = GL.GetUniformLocation(shaderProgram, "haze3dProjection"); // 确保名称严格匹配着色器 if (projectionLoc == -1) { Console.WriteLine("ERROR: Uniform 'haze3dProjection' not found in shader!"); } if (projectionLoc != -1) GL.UniformMatrix4(projectionLoc, false, ref projectionMatrix); else Console.WriteLine("Warning: Uniform 'haze3dProjection' not found in shader!"); //Pass the fog effect parameters int fogStartLoc = GL.GetUniformLocation(shaderProgram, "fogStart"); if (fogStartLoc != -1) GL.Uniform1(fogStartLoc, FogStart); else Console.WriteLine("Warning: Uniform 'fogStart' not found in shader!"); int fogEndLoc = GL.GetUniformLocation(shaderProgram, "fogEnd"); if (fogEndLoc != -1) GL.Uniform1(fogEndLoc, FogEnd); else Console.WriteLine("Warning: Uniform 'fogEnd' not found in shader!"); } public Vector3 GetMouseViewPos(Vector2 mousePos) { throw new NotImplementedException(); } // processes input received from a mouse input system. Expects the offset value in both the x and y direction. public void ProcessMouseMovement(float xoffset, float yoffset) { xoffset *= MouseSensitivity; yoffset *= MouseSensitivity; Yaw += xoffset; Pitch += yoffset; // make sure that when pitch is out of bounds, screen doesn't get flipped Pitch = MathHelper.Clamp(Pitch, -89.0f, 89.0f); // Normalize yaw to 0-360 degrees Yaw = Yaw % 360.0f; // update Front, Right and Up Vectors using the updated Euler angles updateCameraVectors(); } // processes input received from a mouse scroll-wheel event. Only requires input on the vertical wheel-axis public void ProcessMouseScroll(float x, float y, float delta) { Scale = MathHelper.Clamp(Scale-delta,1.0f,45.0f); } // Update camera direction vectors private void updateCameraVectors() { // calculate the new Front vector Vector3 front = new Vector3 { X = (float)(Math.Cos(MathHelper.DegreesToRadians(Yaw)) * Math.Cos(MathHelper.DegreesToRadians(Pitch))), Y = (float)Math.Sin(MathHelper.DegreesToRadians(Pitch)), Z = (float)(Math.Sin(MathHelper.DegreesToRadians(Yaw)) * Math.Cos(MathHelper.DegreesToRadians(Pitch))), }; Front = Vector3.Normalize(front); // also re-calculate the Right and Up vector Right = Vector3.Normalize(Vector3.Cross(Front, WorldUp)); // normalize the vectors, because their length gets closer to 0 the more you look up or down which results in slower movement. Up = Vector3.Normalize(Vector3.Cross(Right, Front)); } public void ReSet() { //Reset the position Position = new Vector3(0.0f, 0.0f, 5.0f); //Reset Euler angle Yaw = YAW; Pitch = PITCH; //Reset the field of view angle FieldOfView = MathHelper.DegreesToRadians(DEFAULT_ZOOM); ModelMatrix = Matrix4.Identity; //update Front, Right and Up Vectors using the updated Euler angles updateCameraVectors(); } } } 这是Camera的代码 请问对于实现三个矩阵有什么问题呢
08-12
评论 5
成就一亿技术人!
拼手气红包6.0元
还能输入1000个字符
 
红包 添加红包
表情包 插入表情
 条评论被折叠 查看
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值