diff --git a/.github/workflows/mega-linter.yml b/.github/workflows/mega-linter.yml index 209b737..dab429a 100644 --- a/.github/workflows/mega-linter.yml +++ b/.github/workflows/mega-linter.yml @@ -29,7 +29,7 @@ jobs: fetch-depth: 0 - name: Run MegaLinter - uses: oxsecurity/megalinter/flavors/python@v9.4.0 + uses: oxsecurity/megalinter/flavors/python@v10.0.0 env: GITHUB_TOKEN: ${{ secrets.GITHUB_TOKEN }} @@ -40,4 +40,4 @@ jobs: name: megalinter-reports path: | megalinter-reports - mega-linter.log \ No newline at end of file + mega-linter.log diff --git a/.github/workflows/test.yml b/.github/workflows/test.yml index f1d8412..27ade51 100644 --- a/.github/workflows/test.yml +++ b/.github/workflows/test.yml @@ -1,8 +1,10 @@ +--- name: Test Runner PoseData -on: [push] +on: # yamllint disable-line rule:truthy + push: jobs: test: runs-on: self-hosted steps: - - run: echo "Runner is working!" \ No newline at end of file + - run: echo "Runner is working!" diff --git a/.mega-linter.yml b/.mega-linter.yml index b6199a8..9353d01 100644 --- a/.mega-linter.yml +++ b/.mega-linter.yml @@ -1,6 +1,6 @@ --- -# Keep MegaLinter focused on checks that are useful across these Python repos -# before strict typing, spelling, and security-scan baselines are established. +# Keep MegaLinter focused on project-owned code. SKEL is vendored upstream +# material and is therefore excluded from repository-wide lint enforcement. ENABLE_LINTERS: - ACTION_ACTIONLINT - MARKDOWN_MARKDOWNLINT @@ -13,3 +13,4 @@ VALIDATE_ALL_CODEBASE: true PRINT_ALPACA: false SHOW_ELAPSED_TIME: true FLAVOR_SUGGESTIONS: false +FILTER_REGEX_EXCLUDE: '(^SKEL/)' diff --git a/AMASS POSE Data Filtering/check_flagged_frames_against_pose_limits.py b/AMASS POSE Data Filtering/check_flagged_frames_against_pose_limits.py index cce76bd..8f9f3c6 100644 --- a/AMASS POSE Data Filtering/check_flagged_frames_against_pose_limits.py +++ b/AMASS POSE Data Filtering/check_flagged_frames_against_pose_limits.py @@ -134,7 +134,7 @@ def check_all_joint_limit_violations(skel_poses, n_frames): joint_names = [SMPL_JOINT_NAMES.get(j, f'joint_{j}') for j in joints_here] print(f"\nFrame {f}:") print(f" Raw-jump flagged SMPL joints: {joint_names}") - print(f" SKEL limit violations:") + print(" SKEL limit violations:") for name, val, limits in skel_violations[f]: print(f" {name}: {val:.3f} rad, limit {limits}") @@ -145,4 +145,4 @@ def check_all_joint_limit_violations(skel_poses, n_frames): print("\n--- All SKEL limit-violation frames ---") for f, vlist in sorted(skel_violations.items()): - print(f" Frame {f}: {[v[0] for v in vlist]}") \ No newline at end of file + print(f" Frame {f}: {[v[0] for v in vlist]}") diff --git a/AMASS POSE Data Filtering/tempCodeRunnerFile.py b/AMASS POSE Data Filtering/tempCodeRunnerFile.py deleted file mode 100644 index b62efce..0000000 --- a/AMASS POSE Data Filtering/tempCodeRunnerFile.py +++ /dev/null @@ -1,535 +0,0 @@ -import numpy as np -from scipy.spatial.transform import Rotation as R -from pathlib import Path -# import quaternion -import os -import matplotlib -import matplotlib.pyplot as plt -import matplotlib.patches as mpatches -import pyrecest -from pyrecest.filters import * -from pyrecest.distributions import * -import requests -import io -import gdown - -input_path = r"C:\Users\ragha\Desktop\important ids and documents\ml research prof.florian\KIT_Quaternions\3\912_3_01_poses_quaternions.npz" - - -file_id = "1F-XL8Tf59lbakNkMhEPjvkr3wT4O0MEW" # File Id of Motion Clip on Google Drive (.npz file) - -def process_sequence(file_id): - buffer = io.BytesIO() - gdown.download(id=file_id, output=buffer, quiet=False) - buffer.seek(0) - data = np.load(buffer, allow_pickle=True) - print(data.files) - poses = data['poses_quat'] - trans= data['trans'] - betas= data['betas'] - gender= data['gender'] - dmpls= data['dmpls'] - mocap_framerate= data['mocap_framerate'] - print(len(poses)) - return trans, betas, gender, dmpls, mocap_framerate, poses - -trans, betas, gender, dmpls, mocap_framerate, poses = process_sequence(file_id) - -# data = np.load(input_path, allow_pickle=True) - -joint_idx = 16 # joint used below for diagnostic plots -num_frames = len(poses) -num_joints = poses.shape[1] -n_particles = 1000 -process_noise_std = 0.01 -high_measurement_kappa = 500.0 -low_measurement_kappa = 50.0 -jump_threshold_std_multiplier = 3.0 -MAX_SAFE_KAPPA = 500.0 - - - -class OnlineKappaEstimator: - """ - Maintains a per-joint running estimate of mean jump and std, - then recomputes kappa bounds each frame via EMA. - """ - def __init__(self, n_joints, alpha=0.02, - initial_mean=None, initial_std=None): - """ - alpha: EMA smoothing factor. - ~0.02 → memory of ~50 frames (slow adaptation, stable). - ~0.1 → memory of ~10 frames (fast adaptation, responsive). - """ - self.alpha = alpha - self.ema_mean = initial_mean.copy() if initial_mean is not None \ - else np.full(n_joints, 0.1) - self.ema_var = (initial_std ** 2).copy() if initial_std is not None \ - else np.full(n_joints, 0.01) - - # def update(self, jumps): - # """jumps: shape (n_joints,) — geodesic distances for current frame.""" - # self.ema_mean += self.alpha * (jumps - self.ema_mean) - # self.ema_var += self.alpha * ( - # (jumps - self.ema_mean) ** 2 - self.ema_var - # ) - - def update(self, jumps): - """jumps: shape (n_joints,) — geodesic distances for current frame.""" - delta_pre = jumps - self.ema_mean # residual w.r.t. OLD mean - - self.ema_mean += self.alpha * delta_pre # update mean - - delta_post = jumps - self.ema_mean # residual w.r.t. NEW mean - - # Welford-style: cross product of pre/post deltas gives unbiased variance update - self.ema_var = (1.0 - self.alpha) * self.ema_var + \ - self.alpha * delta_pre * delta_post # ← correct - - - - def kappa_bounds(self, std_multiplier=3.0): - ema_std = np.sqrt(np.maximum(self.ema_var, 1e-8)) - threshold = self.ema_mean + std_multiplier * ema_std - - high_kappa = np.clip(1.0 / (2.0 * self.ema_mean ** 2), 10.0, 500.0) # ← 500 not 2000 - low_kappa = np.clip(1.0 / (2.0 * threshold ** 2), 1.0, 100.0) - - return high_kappa, low_kappa, threshold - - - - - - - -def normalize_quat(q): - q = q / np.linalg.norm(q, axis=-1, keepdims=True) - q = np.array(q, copy=True) - q *= np.where(q[..., -1:] < 0, -1.0, 1.0) - return q - - -def amass_to_filter_quat(q_amass): - """Convert AMASS quaternion format (w, x, y, z) to PyRecEst format (x, y, z, w).""" - q = np.array([q_amass[1], q_amass[2], q_amass[3], q_amass[0]]) - return normalize_quat(q) - - -def initialize_particles(q, n_particles, noise_std): - particles = q + np.random.randn(n_particles, 4) * noise_std - return normalize_quat(particles) - - -def quat_geodesic_distance(q1, q2): - """Angular distance between two quaternions in radians.""" - q1 = normalize_quat(q1) - q2 = normalize_quat(q2) - dot = np.clip(np.abs(np.dot(q1, q2)), 0, 1) # abs handles double cover - return 2 * np.arccos(dot) - - -def compute_observed_jump_statistics(poses): - observed_jumps = np.zeros((num_frames - 1, num_joints)) - - for frame_idx in range(1, num_frames): - for current_joint_idx in range(num_joints): - q_prev = amass_to_filter_quat(poses[frame_idx - 1, current_joint_idx, :]) - q_curr = amass_to_filter_quat(poses[frame_idx, current_joint_idx, :]) - observed_jumps[frame_idx - 1, current_joint_idx] = quat_geodesic_distance(q_prev, q_curr) - - mean_jumps = np.mean(observed_jumps, axis=0) - std_jumps = np.std(observed_jumps, axis=0) - thresholds = mean_jumps + jump_threshold_std_multiplier * std_jumps - - return observed_jumps, mean_jumps, std_jumps, thresholds - -def compute_kappa_bounds(mean_jump_per_joint, std_jump_per_joint, jump_thresholds): - """ - high_kappa: tight trust — derived from typical (mean) jump size. - low_kappa: loose trust — derived from the anomaly threshold. - Clipped to sane physical limits to avoid numerical issues. - """ - # At normal motion: spread ≈ mean_jump → kappa_high ≈ 1/(2*mean²) - high_kappa = np.clip( - 1.0 / (2.0 * mean_jump_per_joint ** 2 + 1e-8), - 10.0, 2000.0 - ) - # At anomalous motion: spread ≈ threshold → kappa_low ≈ 1/(2*threshold²) - low_kappa = np.clip( - 1.0 / (2.0 * jump_thresholds ** 2 + 1e-8), - 1.0, 200.0 - ) - return high_kappa, low_kappa - - - -# def adaptive_measurement_kappa(current_jump, mean_jump, threshold): -# if current_jump <= mean_jump: -# return high_measurement_kappa - -# if current_jump >= threshold or np.isclose(threshold, mean_jump): -# return low_measurement_kappa - -# jump_ratio = (current_jump - mean_jump) / (threshold - mean_jump) -# return high_measurement_kappa - jump_ratio * (high_measurement_kappa - low_measurement_kappa) - - -observed_jump_magnitudes, mean_jump_per_joint, std_jump_per_joint, jump_thresholds = ( - compute_observed_jump_statistics(poses) -) - -high_kappa_per_joint, low_kappa_per_joint = compute_kappa_bounds( - mean_jump_per_joint, std_jump_per_joint, jump_thresholds -) - - - -# EMA estimator — one object, tracks ALL joints simultaneously -# warm-started so frame 1 already has sensible bounds -kappa_estimator = OnlineKappaEstimator( - n_joints=num_joints, - alpha=0.02, - initial_mean=mean_jump_per_joint, - initial_std=std_jump_per_joint, -) - - - - -def adaptive_measurement_kappa(current_jump, mean_jump, threshold, - high_kappa, low_kappa, smoothness=3.0): - """ - Smooth exponential blend between high_kappa and low_kappa. - All inputs are scalars (called per-joint inside the loop). - smoothness: controls how sharply kappa drops as jump approaches threshold. - """ - jump_ratio = np.clip( - (current_jump - mean_jump) / (threshold - mean_jump + 1e-9), - 0.0, 1.0 - ) - # Exponential decay: stays near high_kappa until ratio rises, then drops fast - blend = 1.0 - np.exp(-smoothness * (1.0 - jump_ratio)) - # Remap so blend=0 → high_kappa, blend=1 → low_kappa - weight = np.exp(-smoothness * jump_ratio) - return high_kappa * weight + low_kappa * (1.0 - weight) - - - - - -def precompute_kappa_table(poses, num_frames, num_joints, - std_multiplier=3.0, smoothness=3.0): - """ - PASS 1: Single loop over all frames and joints. - Returns kappa_table shape (num_frames, num_joints) — - exact per-joint per-frame kappa, no EMA approximation. - """ - - # ── Step 1: compute all geodesic jumps ────────────────────────────────── - jump_magnitudes = np.zeros((num_frames - 1, num_joints)) - - for frame_idx in range(1, num_frames): - for joint_idx in range(num_joints): - q_prev = amass_to_filter_quat(poses[frame_idx - 1, joint_idx, :]) - q_curr = amass_to_filter_quat(poses[frame_idx, joint_idx, :]) - jump_magnitudes[frame_idx - 1, joint_idx] = quat_geodesic_distance(q_prev, q_curr) - - # ── Step 2: per-joint statistics (axis=0 → over frames) ───────────────── - mean_jump = np.mean(jump_magnitudes, axis=0) # shape (n_joints,) - std_jump = np.std(jump_magnitudes, axis=0) # shape (n_joints,) - threshold = mean_jump + std_multiplier * std_jump # shape (n_joints,) - - # ── Step 3: per-joint kappa bounds from statistics ─────────────────────── - high_kappa = np.clip(1.0 / (2.0 * mean_jump ** 2 + 1e-8), 10.0, 500.0) - low_kappa = np.clip(1.0 / (2.0 * threshold ** 2 + 1e-8), 1.0, 100.0) - - # ── Step 4: per-frame per-joint kappa via smooth blend ─────────────────── - # jump_magnitudes shape: (num_frames-1, num_joints) - # mean_jump, threshold shape: (num_joints,) → broadcast over frames - - jump_ratio = np.clip( - (jump_magnitudes - mean_jump) / (threshold - mean_jump + 1e-9), - 0.0, 1.0 - ) # shape (num_frames-1, num_joints) - - weight = np.exp(-smoothness * jump_ratio) # shape (num_frames-1, num_joints) - - kappa_table = np.clip( - high_kappa * weight + low_kappa * (1.0 - weight), - 1.0, 500.0 - ) # shape (num_frames-1, num_joints) - - return kappa_table, jump_magnitudes, mean_jump, std_jump, threshold - - -kappa_table, observed_jump_magnitudes, mean_jump_per_joint, \ - std_jump_per_joint, jump_thresholds = precompute_kappa_table( - poses, num_frames, num_joints - ) - -# One independent particle filter is maintained for each joint orientation. -particle_filters = [ - HyperhemisphericalParticleFilter(n_particles=n_particles, dim=3) - for _ in range(num_joints) -] - -# Store estimated orientations for every frame and every joint in (x, y, z, w) format. -estimates = np.zeros((num_frames, num_joints, 4)) -measurement_kappas = np.full((num_frames, num_joints), high_measurement_kappa) - -for current_joint_idx, pf in enumerate(particle_filters): - q0 = amass_to_filter_quat(poses[0, current_joint_idx, :]) - particles = initialize_particles(q0, n_particles, process_noise_std) - pf.set_state(HyperhemisphericalDiracDistribution(particles)) - estimates[0, current_joint_idx, :] = normalize_quat(pf.filter_state.mean()) - -for frame_idx in range(1, num_frames): - # Called ONCE per frame — returns shape (n_joints,) arrays - # EMA state at this point reflects all frames seen so far - # hk, lk, thresh_online = kappa_estimator.kappa_bounds() - - current_frame_jumps = observed_jump_magnitudes[frame_idx - 1] - - - for current_joint_idx, pf in enumerate(particle_filters): - # STEP 1 - PREDICT (random walk with noise) - particles = pf.filter_state.d - particles = particles + np.random.randn(len(particles), 4) * process_noise_std - particles = normalize_quat(particles) - pf.filter_state.d = particles - - # STEP 2 - OBSERVE current joint quaternion - q_obs = amass_to_filter_quat(poses[frame_idx, current_joint_idx, :]) - - # STEP 3 - UPDATE (reweight particles against observation) - current_jump = observed_jump_magnitudes[frame_idx - 1, current_joint_idx] - # measurement_kappa = adaptive_measurement_kappa( - # current_jump, - # mean_jump_per_joint[current_joint_idx], - # jump_thresholds[current_joint_idx], - # ) - # measurement_kappa = adaptive_measurement_kappa( - # current_jump = current_frame_jumps[current_joint_idx], - # mean_jump = kappa_estimator.ema_mean[current_joint_idx], - # threshold = thresh_online[current_joint_idx], - # high_kappa = hk[current_joint_idx], - # low_kappa = lk[current_joint_idx], - # ) - measurement_kappa = float(kappa_table[frame_idx - 1, current_joint_idx]) - measurement_kappas[frame_idx, current_joint_idx] = measurement_kappa - # measurement_kappa = float(np.clip(measurement_kappa, 1.0, MAX_SAFE_KAPPA)) - # measurement_kappas[frame_idx, current_joint_idx] = measurement_kappa - meas_noise = HyperhemisphericalWatsonDistribution(q_obs, kappa=measurement_kappa) - pf.update_nonlinear_using_likelihood(meas_noise.pdf) - - # STEP 4 - GET ESTIMATE - estimates[frame_idx, current_joint_idx, :] = normalize_quat(pf.filter_state.mean()) - # EMA update — AFTER all joints processed for this frame - # uses current_frame_jumps shape (n_joints,) — updates all joints at once - # kappa_estimator.update(current_frame_jumps) - -# Optional copy in AMASS/numpy-quaternion order (w, x, y, z), useful when saving -# estimates alongside the original AMASS pose data. -estimates_amass_order = estimates[:, :, [3, 0, 1, 2]] - -#Calculate differnce in observed ad estimated orientations for every frame - -# Calculate jump magnitude between consecutive observed frames -jump_magnitudes = [] -for frame_idx in range(1, num_frames): - q1 = estimates[frame_idx - 1, joint_idx, :] # previous frame - q2 = estimates[frame_idx, joint_idx, :] # current frame - q1 = normalize_quat(q1) - q2 = normalize_quat(q2) - - dist = quat_geodesic_distance(q1, q2) - jump_magnitudes.append(dist) - -jump_magnitudes = np.array(jump_magnitudes) - -# # Adaptive threshold: mean + 3*std -mean_jump = np.mean(jump_magnitudes) -std_jump = np.std(jump_magnitudes) -threshold = mean_jump + 3 * std_jump - -# fig, ax = plt.subplots(figsize=(14, 5)) -# frames = np.arange(1, num_frames) - -# # --- split into normal / anomaly series for cleaner legend --- -# anomaly_mask = jump_magnitudes > threshold -# normal_mask = ~anomaly_mask - -# # main jump line -# ax.plot(frames, jump_magnitudes, color='#378ADD', linewidth=1.2, -# alpha=0.85, zorder=2, label='jump magnitude') - -# # shade under the line -# ax.fill_between(frames, jump_magnitudes, alpha=0.08, color='#378ADD', zorder=1) - -# # threshold + mean lines -# ax.axhline(threshold, color='#E24B4A', linewidth=1.4, linestyle='--', -# zorder=3, label=f'threshold μ+3σ ({threshold:.3f} rad)') -# ax.axhline(mean_jump, color='#888780', linewidth=1.0, linestyle=':', -# zorder=3, label=f'mean ({mean_jump:.3f} rad)') - -# # anomaly scatter -# ax.scatter(frames[anomaly_mask], jump_magnitudes[anomaly_mask], -# color='#E24B4A', s=55, zorder=5, label=f'anomaly (n={anomaly_mask.sum()})') - -# # vertical drop-lines from anomaly dots to x-axis (optional, aids reading) -# for f, v in zip(frames[anomaly_mask], jump_magnitudes[anomaly_mask]): -# ax.vlines(f, 0, v, color='#E24B4A', linewidth=0.6, alpha=0.35, zorder=4) - -# # --- shaded band: mean ± 1σ --- -# ax.axhspan(mean_jump - std_jump, mean_jump + std_jump, -# color='#888780', alpha=0.07, zorder=0, label='±1σ band') - -# # labels & formatting -# ax.set_xlabel('Frame index', fontsize=11) -# ax.set_ylabel('Geodesic distance (rad)', fontsize=11) -# ax.set_title(f'Orientation jump magnitudes — joint {joint_idx}', fontsize=13, fontweight='normal') -# ax.set_xlim(frames[0], frames[-1]) -# ax.set_ylim(bottom=0) -# ax.grid(True, linewidth=0.4, alpha=0.5, linestyle='--') -# ax.spines[['top', 'right']].set_visible(False) -# ax.legend(fontsize=9, framealpha=0.85, loc='upper right') - -# # annotate anomaly frame indices -# for f, v in zip(frames[anomaly_mask], jump_magnitudes[anomaly_mask]): -# ax.annotate(f'f{f}', xy=(f, v), xytext=(4, 6), -# textcoords='offset points', fontsize=8, -# color='#E24B4A', fontweight='bold') - -# plt.tight_layout() -# plt.savefig(f'jump_magnitudes_joint{joint_idx}.png', dpi=150, bbox_inches='tight') -# plt.show() - -#2. Sudden jumps compared to estimated orientations -# Distance between filter estimate and observation each frame -filter_residuals = [] -for frame_idx in range(1, num_frames): - q_est = estimates[frame_idx, joint_idx, :] - - q_obs_raw = poses[frame_idx, joint_idx, :] - q_obs = np.array([q_obs_raw[1], q_obs_raw[2], q_obs_raw[3], q_obs_raw[0]]) - if q_obs[-1] < 0: q_obs = -q_obs - - dist = quat_geodesic_distance(q_est, q_obs) - filter_residuals.append(dist) - -filter_residuals = np.array(filter_residuals) - -# Calculate residuals: estimate vs observation for every frame -residuals = np.zeros(num_frames) - -for frame_idx in range(1, num_frames): - # filter estimate - q_est = estimates[frame_idx, joint_idx, :] - - # observed quaternion - q_obs_raw = poses[frame_idx, joint_idx, :] - q_obs = np.array([q_obs_raw[1], q_obs_raw[2], q_obs_raw[3], q_obs_raw[0]]) - if q_obs[-1] < 0: - q_obs = -q_obs - - residuals[frame_idx] = quat_geodesic_distance(q_est, q_obs) - -# Threshold -# mean_res = np.mean(residuals[1:]) -# std_res = np.std(residuals[1:]) -# threshold = mean_res + 3 * std_res - -# # Incorrect frames -# incorrect_frames = np.where(residuals > threshold)[0] - -# print(f"Mean residual : {np.degrees(mean_res):.2f} degrees") -# print(f"Std residual : {np.degrees(std_res):.2f} degrees") -# print(f"Threshold : {np.degrees(threshold):.2f} degrees") -# print(f"Incorrect frames: {incorrect_frames}") -# print(f"Residuals at incorrect frames (degrees):") -# for f in incorrect_frames: -# print(f" frame {f:4d}: {np.degrees(residuals[f]):.2f}°") - - -# plt.figure(figsize=(12, 4)) -# plt.plot(np.degrees(residuals), label='filter residual', color='steelblue') -# plt.axhline(np.degrees(threshold), color='red', linestyle='--', label=f'threshold ({np.degrees(threshold):.1f}°)') -# plt.scatter(incorrect_frames, np.degrees(residuals[incorrect_frames]), -# color='red', zorder=5, s=50, label=f'incorrect ({len(incorrect_frames)} frames)') -# plt.xlabel('Frame') -# plt.ylabel('Angular error (degrees)') -# plt.title('Incorrect observed orientations — Joint 0') -# plt.legend() -# plt.tight_layout() -# plt.show() - - -fig, axes = plt.subplots(2, 1, figsize=(14, 9), sharex=True) -fig.subplots_adjust(hspace=0.08) # tight gap since x-axis is shared - -frames = np.arange(1, num_frames) - -# ── shared helpers ──────────────────────────────────────────────────────────── -def plot_jump_panel(ax, data, label_y, title, color='#378ADD'): - mean_v = np.mean(data) - std_v = np.std(data) - thresh = mean_v + 3 * std_v - mask = data > thresh - - ax.plot(frames, data, color=color, linewidth=1.2, alpha=0.85, zorder=2, - label=label_y) - ax.fill_between(frames, data, alpha=0.08, color=color, zorder=1) - - ax.axhline(thresh, color='#E24B4A', linewidth=1.4, linestyle='--', zorder=3, - label=f'threshold μ+3σ ({thresh:.3f} rad)') - ax.axhline(mean_v, color='#888780', linewidth=1.0, linestyle=':', zorder=3, - label=f'mean ({mean_v:.3f} rad)') - ax.axhspan(mean_v - std_v, mean_v + std_v, - color='#888780', alpha=0.07, zorder=0, label='±1σ band') - - ax.scatter(frames[mask], data[mask], - color='#E24B4A', s=55, zorder=5, - label=f'anomaly (n={mask.sum()})') - for f, v in zip(frames[mask], data[mask]): - ax.vlines(f, 0, v, color='#E24B4A', linewidth=0.6, alpha=0.35, zorder=4) - ax.annotate(f'f{f}', xy=(f, v), xytext=(4, 6), - textcoords='offset points', fontsize=8, - color='#E24B4A', fontweight='bold') - - ax.set_ylabel('Geodesic distance (rad)', fontsize=10) - ax.set_title(title, fontsize=12, fontweight='normal', pad=6) - ax.set_ylim(bottom=0) - ax.grid(True, linewidth=0.4, alpha=0.5, linestyle='--') - ax.spines[['top', 'right']].set_visible(False) - ax.legend(fontsize=8.5, framealpha=0.85, loc='upper right') - - return thresh, mask # caller can use if needed - -# ── panel 1 : consecutive-frame jumps ──────────────────────────────────────── -plot_jump_panel( - axes[0], jump_magnitudes, - label_y='jump magnitude', - title=f'Consecutive-frame orientation jumps — joint {joint_idx}', -) - -# ── panel 2 : filter residuals (estimate vs observation) ───────────────────── -plot_jump_panel( - axes[1], filter_residuals, - label_y='filter residual', - title=f'Filter residual (estimate vs observation) — joint {joint_idx}', - color='#1D9E75', # teal to distinguish from panel 1 -) - -axes[1].set_xlabel('Frame index', fontsize=10) -axes[0].set_xlim(frames[0], frames[-1]) # shared x propagates automatically - -plt.suptitle(f'Joint {joint_idx} — orientation diagnostics', fontsize=13, - y=1.01, fontweight='normal') - -plt.tight_layout() -output_dir = Path("outputs") -output_dir.mkdir(exist_ok=True) -plt.savefig(output_dir / f'orientation_diagnostics_joint{joint_idx}.png', dpi=150, - bbox_inches='tight') -plt.show() diff --git a/AMASS POSE Data Filtering/temp_check_skel_anomalies.py b/AMASS POSE Data Filtering/temp_check_skel_anomalies.py deleted file mode 100644 index 4638be4..0000000 --- a/AMASS POSE Data Filtering/temp_check_skel_anomalies.py +++ /dev/null @@ -1,47 +0,0 @@ -import pickle -import numpy as np -from skel.kin_skel import pose_param_names, pose_limits -# skel_pkl_path = r"C:\Users\ragha\Desktop\important ids and documents\ml research prof.florian\PoseDataParticleFilter\SKEL\SKEL\data\skel\skel_male.pkl" - -skel_pkl_path = r"C:\Users\ragha\Desktop\important ids and documents\ml research prof.florian\PoseDataParticleFilter\SKEL\SKEL\output\Subject_1_F_1_poses\Subject_1_F_1_poses_skel.pkl" -with open(skel_pkl_path, "rb") as f: - skel_data = pickle.load(f) -skel_poses = skel_data['poses'] # (num_frames, 46) -num_frames = skel_poses.shape[0] - - - -def fitted_pose_reference_ranges(raw_limits): - """Return ranges in the signs used by fitted SKEL pose vectors. - - ``kin_skel.pose_limits`` stores the right scapula-elevation range with - the opposite sign. SKEL applies the right parameter directly, but - negates the left parameter in ``left_scapula``. Therefore a symmetric - pose has positive right elevation and negative left elevation. - """ - ranges = {name: tuple(sorted(bounds)) for name, bounds in raw_limits.items()} - ranges.update({ - 'scapula_elevation_r': (0.1, 0.4), - 'scapula_elevation_l': (-0.4, -0.1), - }) - return ranges - - -# These are diagnostic reference ranges. The SKEL fitter does not enforce -# them as hard constraints during optimization, so an out-of-range value is -# a fit-quality signal, not proof that the source AMASS frame is invalid. -pose_limits = fitted_pose_reference_ranges(pose_limits) - -print(f"{'Parameter':<25} {'ViolationRate':>13} {'DataMin':>10} {'DataMax':>10} {'LimitLo':>10} {'LimitHi':>10} {'Overlap?':>9}") -print("-" * 95) - - -for i, name in enumerate(pose_param_names): - if name not in pose_limits: - continue - lo, hi = pose_limits[name] - vals = skel_poses[:, i] - violations = (vals < lo) | (vals > hi) - rate = violations.mean() - overlaps = (vals.min() <= hi) and (vals.max() >= lo) - print(f"{name:<25} {rate*100:>12.1f}% {vals.min():>10.3f} {vals.max():>10.3f} {lo:>10.3f} {hi:>10.3f} {str(overlaps):>9}") diff --git a/Detect_ground_penetration/tempCodeRunnerFile.py b/Detect_ground_penetration/tempCodeRunnerFile.py deleted file mode 100644 index ffbb373..0000000 --- a/Detect_ground_penetration/tempCodeRunnerFile.py +++ /dev/null @@ -1 +0,0 @@ -you \ No newline at end of file diff --git a/Quaternion Representation/amaas_data_quaternions_representaion.py b/Quaternion Representation/amaas_data_quaternions_representaion.py index 9ff566e..c896dda 100644 --- a/Quaternion Representation/amaas_data_quaternions_representaion.py +++ b/Quaternion Representation/amaas_data_quaternions_representaion.py @@ -1,8 +1,5 @@ import numpy as np -from scipy.spatial.transform import Rotation as R -from pathlib import Path import quaternion -import os input_path= r"C:\Users\ragha\Desktop\important ids and documents\ml research prof.florian\PoseDataParticleFilter\Subject_1_F_1_poses_corrupted.npz" output_path= r"C:\Users\ragha\Desktop\important ids and documents\ml research prof.florian\PoseDataParticleFilter\Subject_1_F_1_poses_corrupted_quaternions.npz" @@ -98,4 +95,4 @@ # print("root_orient_quat :", root_orient_quat.shape) # (T, 1, 4) # print("body_pose_quat :", body_pose_quat.shape) # (T, 21, 4) # print("hand_pose_quat :", hand_pose_quat.shape) # (T, 30, 4) -# print("trans :", trans.shape) # (T, 3) \ No newline at end of file +# print("trans :", trans.shape) # (T, 3) diff --git a/Quaternion Representation/tempCodeRunnerFile.py b/Quaternion Representation/tempCodeRunnerFile.py deleted file mode 100644 index 9ff566e..0000000 --- a/Quaternion Representation/tempCodeRunnerFile.py +++ /dev/null @@ -1,101 +0,0 @@ -import numpy as np -from scipy.spatial.transform import Rotation as R -from pathlib import Path -import quaternion -import os - -input_path= r"C:\Users\ragha\Desktop\important ids and documents\ml research prof.florian\PoseDataParticleFilter\Subject_1_F_1_poses_corrupted.npz" -output_path= r"C:\Users\ragha\Desktop\important ids and documents\ml research prof.florian\PoseDataParticleFilter\Subject_1_F_1_poses_corrupted_quaternions.npz" - -# Create output folder if it doesn't exist -# os.makedirs(OUTPUT_FOLDER, exist_ok=True) - -#Get All files in the input folder -# all_files = [f for f in os.listdir(INPUT_FOLDER) if f.endswith('.npz')] -# print(f"Found {len(all_files)} files to process") - -# Load your .npz file -# for idx, filename in enumerate(all_files): -# input_path = os.path.join(INPUT_FOLDER, filename) -# output_path = os.path.join(OUTPUT_FOLDER, filename.replace('.npz', '_quaternions.npz')) - -# print(f"\n[{idx+1}/{len(all_files)}] Processing: {filename}") - -try: - data = np.load(input_path, allow_pickle=True) - - trans = data['trans'] - gender = str(data['gender']) - framerate = float(data['mocap_framerate']) - betas = data['betas'] # [16] - dmpls = data['dmpls'] # [T, 8] - poses = data['poses'] # [T, 156] - - T = poses.shape[0] - - # ── SMPL-H joint layout (156 = 52 joints × 3) ────────────────────────────── - # poses[:, 0:3] → root orientation (1 joint) - # poses[:, 3:66] → body joints (21 joints) - # poses[:, 66:156] → hand joints (30 joints, 15 per hand) - - NUM_JOINTS = 52 # SMPL-H total - - - # ── Reshape to [T, 52, 3] axis-angle ─────────────────────────────────────── - poses_aa = poses.reshape(T, NUM_JOINTS, 3) # axis-angle per joint - Quaternion_arr = quaternion.from_rotation_vector(poses_aa) - quaternion_poses_arr = quaternion.as_float_array(Quaternion_arr) - - # Verify - norms = np.linalg.norm(quaternion_poses_arr, axis=-1) - assert np.allclose(norms, 1.0), "Quaternion norms are not 1.0!" - - np.savez( - output_path, - trans = data['trans'], - betas = data['betas'], - gender = data['gender'], - dmpls = data['dmpls'], - mocap_framerate = data['mocap_framerate'], - poses_quat = quaternion_poses_arr, # (T, 52, 4) [w, x, y, z] - ) - print(f"Saved | frames={T} | shape={quaternion_poses_arr.shape} | path={output_path}") - - -except Exception as e: - print(f"FAILED: Subject_1_F_1_poses_quaternion | Error: {e}") - raise - - - -# convert axis angle to quaternion for one frame [Not in Use] - - # # print(poses.shape) - # poses_aa = poses.reshape(T, NUM_JOINTS, 3) # axis-angle per joint - # # print(poses_aa.shape) - # # temp_vect= poses_aa[0] - # # frame_a= [temp_vect[0][0], temp_vect[0][1], temp_vect[0][2]] - # Quaternion_arr = quaternion.from_rotation_vector(poses_aa) - # quaternion_poses_arr = quaternion.as_float_array(Quaternion_arr) - # # print(poses_quat_arr.shape) - - -# print(betas) -# # print(poses_aa.shape) -# # ── Convert axis-angle → quaternion [T, 52, 4] ───────────────────────────── -# # scipy Rotation expects shape [N, 3] so flatten time & joints, then reshape back -# poses_flat = poses_aa.reshape(-1, 3) # [T*52, 3] -# rots = R.from_rotvec(poses_flat) # axis-angle = rotvec -# quats = rots.as_quat() # [T*52, 4] (x, y, z, w) -# poses_quat = quats.reshape(T, NUM_JOINTS, 4) # [T, 52, 4] -# print(quats.shape) -# print(poses_quat) -# # ── Split into semantic groups ────────────────────────────────────────────── -# root_orient_quat = poses_quat[:, 0:1, :] # [T, 1, 4] global orientation -# body_pose_quat = poses_quat[:, 1:22, :] # [T, 21, 4] body joints -# hand_pose_quat = poses_quat[:, 22:, :] # [T, 30, 4] hand joints - -# print("root_orient_quat :", root_orient_quat.shape) # (T, 1, 4) -# print("body_pose_quat :", body_pose_quat.shape) # (T, 21, 4) -# print("hand_pose_quat :", hand_pose_quat.shape) # (T, 30, 4) -# print("trans :", trans.shape) # (T, 3) \ No newline at end of file diff --git a/inject_synthetic_anomalies.py b/inject_synthetic_anomalies.py index 8bd88c7..3c92158 100644 --- a/inject_synthetic_anomalies.py +++ b/inject_synthetic_anomalies.py @@ -290,11 +290,11 @@ def window_free(start, end): print(f"\nInjected {len(events)} events covering " f"{int(bool_mask.any(axis=1).sum())} / {n_frames} frames " f"({100 * bool_mask.any(axis=1).sum() / n_frames:.1f}%) as positives.") -print(f"\nBy type:") +print("\nBy type:") for t in TYPES: n = sum(1 for e in events if e["type"] == t) print(f" {t:10s}: {n} events") -print(f"\nBy severity:") +print("\nBy severity:") for tier in SEVERITY_DEG: n = sum(1 for e in events if e["severity_tier"] == tier) print(f" {tier:10s}: {n} events")