diff --git a/.git-blame-ignore-revs b/.git-blame-ignore-revs new file mode 100644 index 0000000..8a6937b --- /dev/null +++ b/.git-blame-ignore-revs @@ -0,0 +1 @@ +ef006dbd1e31cc7cae2fae978401a818ee2025d1 diff --git a/.gitignore b/.gitignore index 68bc17f..6119f53 100644 --- a/.gitignore +++ b/.gitignore @@ -158,3 +158,6 @@ cython_debug/ # and can be added to the global gitignore or merged into this file. For a more nuclear # option (not recommended) you can uncomment the following to ignore the entire idea folder. #.idea/ + +test/ +.vscode/ diff --git a/K-ShakeTune/K-SnT_axis.cfg b/K-ShakeTune/K-SnT_axis.cfg index e51d691..9cdb0a2 100644 --- a/K-ShakeTune/K-SnT_axis.cfg +++ b/K-ShakeTune/K-SnT_axis.cfg @@ -13,7 +13,7 @@ gcode: {% set scv = params.SCV|default(None) %} {% set max_sm = params.MAX_SMOOTHING|default(None) %} {% set keep_results = params.KEEP_N_RESULTS|default(3)|int %} - {% set keep_csv = params.KEEP_CSV|default(True) %} + {% set keep_csv = params.KEEP_CSV|default(0)|int %} {% set X, Y = False, False %} @@ -27,17 +27,21 @@ gcode: { action_raise_error("AXIS selection invalid. Should be either all, x or y!") } {% endif %} - {% if scv is none %} + {% if scv is none or scv == "" %} {% set scv = printer.toolhead.square_corner_velocity %} {% endif %} + {% if max_sm == "" %} + {% set max_sm = none %} + {% endif %} + {% if X %} TEST_RESONANCES AXIS=X OUTPUT=raw_data NAME=x FREQ_START={min_freq} FREQ_END={max_freq} HZ_PER_SEC={hz_per_sec} M400 RESPOND MSG="X axis frequency profile generation..." RESPOND MSG="This may take some time (1-3min)" - RUN_SHELL_COMMAND CMD=shaketune PARAMS="--type shaper --scv {scv} {% if max_sm is not none %}--max_smoothing {max_sm}{% endif %} {% if keep_csv %}--keep_csv{% endif %}" + RUN_SHELL_COMMAND CMD=shaketune PARAMS="--type shaper --scv {scv} {% if max_sm is not none %}--max_smoothing {max_sm}{% endif %} {% if keep_csv %}--keep_csv{% endif %} --keep_results {keep_results}" {% endif %} {% if Y %} @@ -46,8 +50,5 @@ gcode: RESPOND MSG="Y axis frequency profile generation..." RESPOND MSG="This may take some time (1-3min)" - RUN_SHELL_COMMAND CMD=shaketune PARAMS="--type shaper --scv {scv} {% if max_sm is not none %}--max_smoothing {max_sm}{% endif %} {% if keep_csv %}--keep_csv{% endif %}" + RUN_SHELL_COMMAND CMD=shaketune PARAMS="--type shaper --scv {scv} {% if max_sm is not none %}--max_smoothing {max_sm}{% endif %} {% if keep_csv %}--keep_csv{% endif %} --keep_results {keep_results}" {% endif %} - - M400 - RUN_SHELL_COMMAND CMD=shaketune PARAMS="--type clean --keep_results {keep_results}" diff --git a/K-ShakeTune/K-SnT_belts.cfg b/K-ShakeTune/K-SnT_belts.cfg index 47bdedc..cd4987c 100644 --- a/K-ShakeTune/K-SnT_belts.cfg +++ b/K-ShakeTune/K-SnT_belts.cfg @@ -3,14 +3,14 @@ ################################################ # Written by Frix_x#0161 # -[gcode_macro BELTS_SHAPER_CALIBRATION] +[gcode_macro COMPARE_BELTS_RESPONSES] description: Perform a custom half-axis test to analyze and compare the frequency profiles of individual belts on CoreXY printers gcode: {% set min_freq = params.FREQ_START|default(5)|float %} {% set max_freq = params.FREQ_END|default(133.33)|float %} {% set hz_per_sec = params.HZ_PER_SEC|default(1)|float %} {% set keep_results = params.KEEP_N_RESULTS|default(3)|int %} - {% set keep_csv = params.KEEP_CSV|default(True) %} + {% set keep_csv = params.KEEP_CSV|default(0)|int %} TEST_RESONANCES AXIS=1,1 OUTPUT=raw_data NAME=b FREQ_START={min_freq} FREQ_END={max_freq} HZ_PER_SEC={hz_per_sec} M400 @@ -20,6 +20,4 @@ gcode: RESPOND MSG="Belts comparative frequency profile generation..." RESPOND MSG="This may take some time (3-5min)" - RUN_SHELL_COMMAND CMD=shaketune PARAMS="--type belts {% if keep_csv %}--keep_csv{% endif %}" - M400 - RUN_SHELL_COMMAND CMD=shaketune PARAMS="--type clean --keep_results {keep_results}" + RUN_SHELL_COMMAND CMD=shaketune PARAMS="--type belts {% if keep_csv %}--keep_csv{% endif %} --keep_results {keep_results}" diff --git a/K-ShakeTune/K-SnT_vibrations.cfg b/K-ShakeTune/K-SnT_vibrations.cfg index 351974b..a0a9ddd 100644 --- a/K-ShakeTune/K-SnT_vibrations.cfg +++ b/K-ShakeTune/K-SnT_vibrations.cfg @@ -1,106 +1,33 @@ -################################################ -###### VIBRATIONS AND SPEED OPTIMIZATIONS ###### -################################################ +######################################### +###### MACHINE VIBRATIONS ANALYSIS ###### +######################################### # Written by Frix_x#0161 # -[gcode_macro VIBRATIONS_CALIBRATION] +[gcode_macro CREATE_VIBRATIONS_PROFILE] gcode: - {% set size = params.SIZE|default(60)|int %} # size of the area where the movements are done - {% set direction = params.DIRECTION|default('XY') %} # can be set to either XY, AB, ABXY, A, B, X, Y, Z + {% set size = params.SIZE|default(100)|int %} # size of the circle where the angled lines are done {% set z_height = params.Z_HEIGHT|default(20)|int %} # z height to put the toolhead before starting the movements - - {% set min_speed = params.MIN_SPEED|default(20)|float * 60 %} # minimum feedrate for the movements {% set max_speed = params.MAX_SPEED|default(200)|float * 60 %} # maximum feedrate for the movements {% set speed_increment = params.SPEED_INCREMENT|default(2)|float * 60 %} # feedrate increment between each move + {% set feedrate_travel = params.TRAVEL_SPEED|default(200)|int * 60 %} # travel feedrate between moves {% set accel = params.ACCEL|default(3000)|int %} # accel value used to move on the pattern {% set accel_chip = params.ACCEL_CHIP|default("adxl345") %} # ADXL chip name in the config {% set keep_results = params.KEEP_N_RESULTS|default(3)|int %} - {% set keep_csv = params.KEEP_CSV|default(True) %} + {% set keep_csv = params.KEEP_CSV|default(0)|int %} {% set mid_x = printer.toolhead.axis_maximum.x|float / 2 %} {% set mid_y = printer.toolhead.axis_maximum.y|float / 2 %} - {% set nb_samples = ((max_speed - min_speed) / speed_increment + 1) | int %} + {% set min_speed = 2 * 60 %} # minimum feedrate for the movements is set to 2mm/s + {% set nb_speed_samples = ((max_speed - min_speed) / speed_increment + 1) | int %} {% set accel = [accel, printer.configfile.settings.printer.max_accel]|min %} {% set old_accel = printer.toolhead.max_accel %} {% set old_cruise_ratio = printer.toolhead.minimum_cruise_ratio %} {% set old_sqv = printer.toolhead.square_corner_velocity %} - {% set direction_factor = { - 'XY' : { - 'start' : {'x': -0.5, 'y': -0.5 }, - 'move_factors' : { - '0' : {'x': 0.5, 'y': -0.5, 'z': 0.0 }, - '1' : {'x': 0.5, 'y': 0.5, 'z': 0.0 }, - '2' : {'x': -0.5, 'y': 0.5, 'z': 0.0 }, - '3' : {'x': -0.5, 'y': -0.5, 'z': 0.0 } - } - }, - 'AB' : { - 'start' : {'x': 0.0, 'y': 0.0 }, - 'move_factors' : { - '0' : {'x': 0.5, 'y': -0.5, 'z': 0.0 }, - '1' : {'x': -0.5, 'y': 0.5, 'z': 0.0 }, - '2' : {'x': 0.0, 'y': 0.0, 'z': 0.0 }, - '3' : {'x': 0.5, 'y': 0.5, 'z': 0.0 }, - '4' : {'x': -0.5, 'y': -0.5, 'z': 0.0 }, - '5' : {'x': 0.0, 'y': 0.0, 'z': 0.0 } - } - }, - 'ABXY' : { - 'start' : {'x': -0.5, 'y': 0.5 }, - 'move_factors' : { - '0' : {'x': -0.5, 'y': -0.5, 'z': 0.0 }, - '1' : {'x': 0.5, 'y': -0.5, 'z': 0.0 }, - '2' : {'x': -0.5, 'y': 0.5, 'z': 0.0 }, - '3' : {'x': 0.5, 'y': 0.5, 'z': 0.0 }, - '4' : {'x': -0.5, 'y': -0.5, 'z': 0.0 }, - '5' : {'x': -0.5, 'y': 0.5, 'z': 0.0 } - } - }, - 'B' : { - 'start' : {'x': 0.5, 'y': 0.5 }, - 'move_factors' : { - '0' : {'x': -0.5, 'y': -0.5, 'z': 0.0 }, - '1' : {'x': 0.5, 'y': 0.5, 'z': 0.0 } - } - }, - 'A' : { - 'start' : {'x': -0.5, 'y': 0.5 }, - 'move_factors' : { - '0' : {'x': 0.5, 'y': -0.5, 'z': 0.0 }, - '1' : {'x': -0.5, 'y': 0.5, 'z': 0.0 } - } - }, - 'X' : { - 'start' : {'x': -0.5, 'y': 0.0 }, - 'move_factors' : { - '0' : {'x': 0.5, 'y': 0.0, 'z': 0.0 }, - '1' : {'x': -0.5, 'y': 0.0, 'z': 0.0 } - } - }, - 'Y' : { - 'start' : {'x': 0.0, 'y': 0.5 }, - 'move_factors' : { - '0' : {'x': 0.0, 'y': -0.5, 'z': 0.0 }, - '1' : {'x': 0.0, 'y': 0.5, 'z': 0.0 } - } - }, - 'Z' : { - 'start' : {'x': 0.0, 'y': 0.0 }, - 'move_factors' : { - '0' : {'x': 0.0, 'y': 0.0, 'z': 1.0 }, - '1' : {'x': 0.0, 'y': 0.0, 'z': 0.0 } - } - }, - 'E' : { - 'start' : {'x': 0.0, 'y': 0.0 }, - 'move_factor' : 0.05 - } - } - %} + {% set kinematics = printer.configfile.settings.printer.kinematics %} {% if not 'xyz' in printer.toolhead.homed_axes %} @@ -111,56 +38,177 @@ gcode: { action_raise_error("Only 2 decimal digits are allowed for SPEED_INCREMENT") } {% endif %} - {% if (size / (max_speed / 60)) < 0.25 and direction != 'E' %} + {% if (size / (max_speed / 60)) < 0.25 %} { action_raise_error("SIZE is too small for this MAX_SPEED. Increase SIZE or decrease MAX_SPEED!") } {% endif %} - {% if not (direction in direction_factor) %} - { action_raise_error("DIRECTION is not valid. Only XY, AB, ABXY, A, B, X, Y, Z or E is allowed!") } - {% endif %} - {action_respond_info("")} - {action_respond_info("Starting speed and vibration calibration")} + {action_respond_info("Starting machine vibrations profile measurement")} {action_respond_info("This operation can not be interrupted by normal means. Hit the \"emergency stop\" button to stop it if needed")} {action_respond_info("")} - SAVE_GCODE_STATE NAME=STATE_VIBRATIONS_CALIBRATION + SAVE_GCODE_STATE NAME=CREATE_VIBRATIONS_PROFILE G90 # Set the wanted acceleration values (not too high to avoid oscillation, not too low to be able to reach constant speed on each segments) SET_VELOCITY_LIMIT ACCEL={accel} MINIMUM_CRUISE_RATIO=0 SQUARE_CORNER_VELOCITY={[(accel / 1000), 5.0]|max} - + # Going to the start position G1 Z{z_height} F{feedrate_travel / 10} - G1 X{mid_x + (size * direction_factor[direction].start.x) } Y{mid_y + (size * direction_factor[direction].start.y)} F{feedrate_travel} - - # vibration pattern for each frequency - {% for curr_sample in range(0, nb_samples) %} - {% set curr_speed = min_speed + curr_sample * speed_increment %} - RESPOND MSG="{"Current speed: %.2f mm/s" % (curr_speed / 60)|float}" + G1 X{mid_x } Y{mid_y} F{feedrate_travel} + + + {% if kinematics == "cartesian" %} + # Cartesian motors are on X and Y axis directly + RESPOND MSG="Cartesian kinematics mode" + {% set main_angles = [0, 90] %} + {% elif kinematics == "corexy" %} + # CoreXY motors are on A and B axis (45 and 135 degrees) + RESPOND MSG="CoreXY kinematics mode" + {% set main_angles = [45, 135] %} + {% else %} + { action_raise_error("Only Cartesian and CoreXY kinematics are supported at the moment for the vibrations measurement tool!") } + {% endif %} - ACCELEROMETER_MEASURE CHIP={accel_chip} - {% if direction == 'E' %} - G0 E{curr_speed*direction_factor[direction].move_factor} F{curr_speed} - {% else %} - {% for key, factor in direction_factor[direction].move_factors|dictsort %} - G1 X{mid_x + (size * factor.x) } Y{mid_y + (size * factor.y)} Z{z_height + (size * factor.z)} F{curr_speed} + {% set pi = (3.141592653589793) | float %} + {% set tau = (pi * 2) | float %} + + + {% for curr_angle in main_angles %} + {% for curr_speed_sample in range(0, nb_speed_samples) %} + {% set curr_speed = min_speed + curr_speed_sample * speed_increment %} + {% set rad_angle_full = (curr_angle|float * pi / 180) %} + + # ----------------------------------------------------------------------------------------------------------- + # Here are some maths to approximate the sin and cos values of rad_angle in Jinja + # Thanks a lot to Aubey! for sharing the idea of using hardcoded Taylor series and + # the associated bit of code to do it easily! This is pure madness! + {% set rad_angle = ((rad_angle_full % tau) - (tau / 2)) | float %} + + {% if rad_angle < (-(tau / 4)) %} + {% set rad_angle = (rad_angle + (tau / 2)) | float %} + {% set final_mult = (-1) %} + {% elif rad_angle > (tau / 4) %} + {% set rad_angle = (rad_angle - (tau / 2)) | float %} + {% set final_mult = (-1) %} + {% else %} + {% set final_mult = (1) %} + {% endif %} + + {% set sin0 = (rad_angle) %} + {% set sin1 = ((rad_angle ** 3) / 6) | float %} + {% set sin2 = ((rad_angle ** 5) / 120) | float %} + {% set sin3 = ((rad_angle ** 7) / 5040) | float %} + {% set sin4 = ((rad_angle ** 9) / 362880) | float %} + {% set sin5 = ((rad_angle ** 11) / 39916800) | float %} + {% set sin6 = ((rad_angle ** 13) / 6227020800) | float %} + {% set sin7 = ((rad_angle ** 15) / 1307674368000) | float %} + {% set sin = (-(sin0 - sin1 + sin2 - sin3 + sin4 - sin5 + sin6 - sin7) * final_mult) | float %} + + {% set cos0 = (1) | float %} + {% set cos1 = ((rad_angle ** 2) / 2) | float %} + {% set cos2 = ((rad_angle ** 4) / 24) | float %} + {% set cos3 = ((rad_angle ** 6) / 720) | float %} + {% set cos4 = ((rad_angle ** 8) / 40320) | float %} + {% set cos5 = ((rad_angle ** 10) / 3628800) | float %} + {% set cos6 = ((rad_angle ** 12) / 479001600) | float %} + {% set cos7 = ((rad_angle ** 14) / 87178291200) | float %} + {% set cos = (-(cos0 - cos1 + cos2 - cos3 + cos4 - cos5 + cos6 - cos7) * final_mult) | float %} + # ----------------------------------------------------------------------------------------------------------- + + # Reduce the segments length for the lower speed range (0-100mm/s). The minimum length is 1/3 of the SIZE and is gradually increased + # to the nominal SIZE at 100mm/s. No further size changes are made above this speed. The goal is to ensure that the print head moves + # enough to collect enough data for vibration analysis, without doing unnecessary distance to save time. At higher speeds, the full + # segments lengths are used because the head moves faster and travels more distance in the same amount of time and we want enough data + {% if curr_speed < (100 * 60) %} + {% set segment_length_multiplier = 1/5 + 4/5 * (curr_speed / 60) / 100 %} + {% else %} + {% set segment_length_multiplier = 1 %} + {% endif %} + + # Calculate angle coordinates using trigonometry and length multiplier and move to start point + {% set dx = (size / 2) * cos * segment_length_multiplier %} + {% set dy = (size / 2) * sin * segment_length_multiplier %} + G1 X{mid_x - dx} Y{mid_y - dy} F{feedrate_travel} + + # Adjust the number of back and forth movements based on speed to also save time on lower speed range + # 3 movements are done by default, reduced to 2 between 150-250mm/s and to 1 under 150mm/s. + {% set movements = 3 %} + {% if curr_speed < (150 * 60) %} + {% set movements = 1 %} + {% elif curr_speed < (250 * 60) %} + {% set movements = 2 %} + {% endif %} + + ACCELEROMETER_MEASURE CHIP={accel_chip} + + # Back and forth movements to record the vibrations at constant speed in both direction + {% for n in range(movements) %} + G1 X{mid_x + dx} Y{mid_y + dy} F{curr_speed} + G1 X{mid_x - dx} Y{mid_y - dy} F{curr_speed} {% endfor %} - {% endif %} - ACCELEROMETER_MEASURE CHIP={accel_chip} NAME=sp{("%.2f" % (curr_speed / 60)|float)|replace('.','_')}n1 - G4 P300 + + ACCELEROMETER_MEASURE CHIP={accel_chip} NAME=an{("%.2f" % curr_angle|float)|replace('.','_')}sp{("%.2f" % (curr_speed / 60)|float)|replace('.','_')} + G4 P300 - M400 + M400 + {% endfor %} {% endfor %} - RESPOND MSG="Machine and motors vibration graph generation..." - RESPOND MSG="This may take some time (3-5min)" - RUN_SHELL_COMMAND CMD=shaketune PARAMS="--type vibrations --axis_name {direction} --accel {accel|int} --chip_name {accel_chip} {% if keep_csv %}--keep_csv{% endif %}" - M400 - RUN_SHELL_COMMAND CMD=shaketune PARAMS="--type clean --keep_results {keep_results}" - # Restore the previous acceleration values SET_VELOCITY_LIMIT ACCEL={old_accel} MINIMUM_CRUISE_RATIO={old_cruise_ratio} SQUARE_CORNER_VELOCITY={old_sqv} - RESTORE_GCODE_STATE NAME=STATE_VIBRATIONS_CALIBRATION + # Extract the TMC names and configuration + {% set ns_x = namespace(path='') %} + {% set ns_y = namespace(path='') %} + + {% for item in printer %} + {% set parts = item.split() %} + {% if parts|length == 2 and parts[0].startswith('tmc') and parts[0][3:].isdigit() %} + {% if parts[1] == 'stepper_x' %} + {% set ns_x.path = parts[0] %} + {% elif parts[1] == 'stepper_y' %} + {% set ns_y.path = parts[0] %} + {% endif %} + {% endif %} + {% endfor %} + + {% if ns_x.path and ns_y.path %} + {% set metadata = + "stepper_x_tmc:" ~ ns_x.path ~ "|" + "stepper_x_run_current:" ~ (printer[ns_x.path + ' stepper_x'].run_current | round(2) | string) ~ "|" + "stepper_x_hold_current:" ~ (printer[ns_x.path + ' stepper_x'].hold_current | round(2) | string) ~ "|" + "stepper_y_tmc:" ~ ns_y.path ~ "|" + "stepper_y_run_current:" ~ (printer[ns_y.path + ' stepper_y'].run_current | round(2) | string) ~ "|" + "stepper_y_hold_current:" ~ (printer[ns_y.path + ' stepper_y'].hold_current | round(2) | string) ~ "|" + %} + + {% set autotune_x = printer.configfile.config['autotune_tmc stepper_x'] if 'autotune_tmc stepper_x' in printer.configfile.config else none %} + {% set autotune_y = printer.configfile.config['autotune_tmc stepper_y'] if 'autotune_tmc stepper_y' in printer.configfile.config else none %} + {% if autotune_x and autotune_y %} + {% set stepper_x_voltage = autotune_x.voltage if autotune_x.voltage else '24.0' %} + {% set stepper_y_voltage = autotune_y.voltage if autotune_y.voltage else '24.0' %} + {% set metadata = metadata ~ + "autotune_enabled:True|" + "stepper_x_motor:" ~ autotune_x.motor ~ "|" + "stepper_x_voltage:" ~ stepper_x_voltage ~ "|" + "stepper_y_motor:" ~ autotune_y.motor ~ "|" + "stepper_y_voltage:" ~ stepper_y_voltage ~ "|" + %} + {% else %} + {% set metadata = metadata ~ "autotune_enabled:False|" %} + {% endif %} + + DUMP_TMC STEPPER=stepper_x + DUMP_TMC STEPPER=stepper_y + + {% else %} + { action_respond_info("No TMC drivers found for X and Y steppers") } + {% endif %} + + RESPOND MSG="Machine vibrations profile generation..." + RESPOND MSG="This may take some time (3-5min)" + RUN_SHELL_COMMAND CMD=shaketune PARAMS="--type vibrations --accel {accel|int} --kinematics {kinematics} {% if metadata %}--metadata {metadata}{% endif %} --chip_name {accel_chip} {% if keep_csv %}--keep_csv{% endif %} --keep_results {keep_results}" + + RESTORE_GCODE_STATE NAME=CREATE_VIBRATIONS_PROFILE diff --git a/K-ShakeTune/scripts/common_func.py b/K-ShakeTune/scripts/common_func.py deleted file mode 100644 index 203054b..0000000 --- a/K-ShakeTune/scripts/common_func.py +++ /dev/null @@ -1,124 +0,0 @@ -#!/usr/bin/env python3 - -# Common functions for the Shake&Tune package -# Written by Frix_x#0161 # - -import math -import os, sys -from importlib import import_module -from pathlib import Path -import numpy as np -from scipy.signal import spectrogram -from git import GitCommandError, Repo - - -def parse_log(logname): - with open(logname) as f: - for header in f: - if not header.startswith('#'): - break - if not header.startswith('freq,psd_x,psd_y,psd_z,psd_xyz'): - # Raw accelerometer data - return np.loadtxt(logname, comments='#', delimiter=',') - # Power spectral density data or shaper calibration data - raise ValueError("File %s does not contain raw accelerometer data and therefore " - "is not supported by Shake&Tune. Please use the official Klipper " - "script to process it instead." % (logname,)) - - -def setup_klipper_import(kdir): - kdir = os.path.expanduser(kdir) - sys.path.append(os.path.join(kdir, 'klippy')) - return import_module('.shaper_calibrate', 'extras') - - -# This is used to print the current S&T version on top of the png graph file -def get_git_version(): - try: - # Get the absolute path of the script, resolving any symlinks - # Then get 2 times to parent dir to be at the git root folder - script_path = Path(__file__).resolve() - repo_path = script_path.parents[2] - repo = Repo(repo_path) - - try: - version = repo.git.describe('--tags') - except GitCommandError: - # If no tag is found, use the simplified commit SHA instead - version = repo.head.commit.hexsha[:7] - return version - - except Exception as e: - return None - - -# This is Klipper's spectrogram generation function adapted to use Scipy -def compute_spectrogram(data): - N = data.shape[0] - Fs = N / (data[-1, 0] - data[0, 0]) - # Round up to a power of 2 for faster FFT - M = 1 << int(.5 * Fs - 1).bit_length() - window = np.kaiser(M, 6.) - - def _specgram(x): - return spectrogram(x, fs=Fs, window=window, nperseg=M, noverlap=M//2, - detrend='constant', scaling='density', mode='psd') - - d = {'x': data[:, 1], 'y': data[:, 2], 'z': data[:, 3]} - f, t, pdata = _specgram(d['x']) - for axis in 'yz': - pdata += _specgram(d[axis])[2] - return pdata, t, f - - -# Compute natural resonant frequency and damping ratio by using the half power bandwidth method with interpolated frequencies -def compute_mechanical_parameters(psd, freqs): - max_power_index = np.argmax(psd) - fr = freqs[max_power_index] - max_power = psd[max_power_index] - - half_power = max_power / math.sqrt(2) - idx_below = np.where(psd[:max_power_index] <= half_power)[0][-1] - idx_above = np.where(psd[max_power_index:] <= half_power)[0][0] + max_power_index - freq_below_half_power = freqs[idx_below] + (half_power - psd[idx_below]) * (freqs[idx_below + 1] - freqs[idx_below]) / (psd[idx_below + 1] - psd[idx_below]) - freq_above_half_power = freqs[idx_above - 1] + (half_power - psd[idx_above - 1]) * (freqs[idx_above] - freqs[idx_above - 1]) / (psd[idx_above] - psd[idx_above - 1]) - - bandwidth = freq_above_half_power - freq_below_half_power - bw1 = math.pow(bandwidth/fr,2) - bw2 = math.pow(bandwidth/fr,4) - - zeta = math.sqrt(0.5-math.sqrt(1/(4+4*bw1-bw2))) - - return fr, zeta, max_power_index - -# This find all the peaks in a curve by looking at when the derivative term goes from positive to negative -# Then only the peaks found above a threshold are kept to avoid capturing peaks in the low amplitude noise of a signal -def detect_peaks(data, indices, detection_threshold, relative_height_threshold=None, window_size=5, vicinity=3): - # Smooth the curve using a moving average to avoid catching peaks everywhere in noisy signals - kernel = np.ones(window_size) / window_size - smoothed_data = np.convolve(data, kernel, mode='valid') - mean_pad = [np.mean(data[:window_size])] * (window_size // 2) - smoothed_data = np.concatenate((mean_pad, smoothed_data)) - - # Find peaks on the smoothed curve - smoothed_peaks = np.where((smoothed_data[:-2] < smoothed_data[1:-1]) & (smoothed_data[1:-1] > smoothed_data[2:]))[0] + 1 - smoothed_peaks = smoothed_peaks[smoothed_data[smoothed_peaks] > detection_threshold] - - # Additional validation for peaks based on relative height - valid_peaks = smoothed_peaks - if relative_height_threshold is not None: - valid_peaks = [] - for peak in smoothed_peaks: - peak_height = smoothed_data[peak] - np.min(smoothed_data[max(0, peak-vicinity):min(len(smoothed_data), peak+vicinity+1)]) - if peak_height > relative_height_threshold * smoothed_data[peak]: - valid_peaks.append(peak) - - # Refine peak positions on the original curve - refined_peaks = [] - for peak in valid_peaks: - local_max = peak + np.argmax(data[max(0, peak-vicinity):min(len(data), peak+vicinity+1)]) - vicinity - refined_peaks.append(local_max) - - num_peaks = len(refined_peaks) - - return num_peaks, np.array(refined_peaks), indices[refined_peaks] diff --git a/K-ShakeTune/scripts/graph_vibrations.py b/K-ShakeTune/scripts/graph_vibrations.py deleted file mode 100755 index 4a7515d..0000000 --- a/K-ShakeTune/scripts/graph_vibrations.py +++ /dev/null @@ -1,466 +0,0 @@ -#!/usr/bin/env python3 - -################################################## -###### SPEED AND VIBRATIONS PLOTTING SCRIPT ###### -################################################## -# Written by Frix_x#0161 # - -# Be sure to make this script executable using SSH: type 'chmod +x ./graph_vibrations.py' when in the folder ! - -##################################################################### -################ !!! DO NOT EDIT BELOW THIS LINE !!! ################ -##################################################################### - -import optparse, matplotlib, re, os, operator -from datetime import datetime -from collections import OrderedDict -import numpy as np -import matplotlib.pyplot as plt -import matplotlib.font_manager, matplotlib.ticker, matplotlib.gridspec - -matplotlib.use('Agg') - -from locale_utils import set_locale, print_with_c_locale -from common_func import compute_mechanical_parameters, detect_peaks, get_git_version, parse_log, setup_klipper_import - - -PEAKS_DETECTION_THRESHOLD = 0.05 -PEAKS_RELATIVE_HEIGHT_THRESHOLD = 0.04 -VALLEY_DETECTION_THRESHOLD = 0.1 # Lower is more sensitive - -KLIPPAIN_COLORS = { - "purple": "#70088C", - "orange": "#FF8D32", - "dark_purple": "#150140", - "dark_orange": "#F24130", - "red_pink": "#F2055C" -} - - -###################################################################### -# Computation -###################################################################### - -# Call to the official Klipper input shaper object to do the PSD computation -def calc_freq_response(data): - helper = shaper_calibrate.ShaperCalibrate(printer=None) - return helper.process_accelerometer_data(data) - - -def compute_vibration_spectrogram(datas, group, max_freq): - psd_list = [] - first_freqs = None - signal_axes = ['x', 'y', 'z', 'all'] - - for i in range(0, len(datas), group): - # Round up to the nearest power of 2 for faster FFT - N = datas[i].shape[0] - T = datas[i][-1,0] - datas[i][0,0] - M = 1 << int((N/T) * 0.5 - 1).bit_length() - if N <= M: - # If there is not enough lines in the array to be able to round up to the - # nearest power of 2, we need to pad some zeros at the end of the array to - # avoid entering a blocking state from Klipper shaper_calibrate.py - datas[i] = np.pad(datas[i], [(0, (M-N)+1), (0, 0)], mode='constant', constant_values=0) - - freqrsp = calc_freq_response(datas[i]) - for n in range(group - 1): - data = datas[i + n + 1] - - # Round up to the nearest power of 2 for faster FFT - N = data.shape[0] - T = data[-1,0] - data[0,0] - M = 1 << int((N/T) * 0.5 - 1).bit_length() - if N <= M: - # If there is not enough lines in the array to be able to round up to the - # nearest power of 2, we need to pad some zeros at the end of the array to - # avoid entering a blocking state from Klipper shaper_calibrate.py - data = np.pad(data, [(0, (M-N)+1), (0, 0)], mode='constant', constant_values=0) - - freqrsp.add_data(calc_freq_response(data)) - - if not psd_list: - # First group, just put it in the result list - first_freqs = freqrsp.freq_bins - psd = freqrsp.psd_sum[first_freqs <= max_freq] - px = freqrsp.psd_x[first_freqs <= max_freq] - py = freqrsp.psd_y[first_freqs <= max_freq] - pz = freqrsp.psd_z[first_freqs <= max_freq] - psd_list.append([psd, px, py, pz]) - else: - # Not the first group, we need to interpolate every new signals - # to the first one to equalize the frequency_bins between them - signal_normalized = dict() - freqs = freqrsp.freq_bins - for axe in signal_axes: - signal = freqrsp.get_psd(axe) - signal_normalized[axe] = np.interp(first_freqs, freqs, signal) - - # Remove data above max_freq on all axes and add to the result list - psd = signal_normalized['all'][first_freqs <= max_freq] - px = signal_normalized['x'][first_freqs <= max_freq] - py = signal_normalized['y'][first_freqs <= max_freq] - pz = signal_normalized['z'][first_freqs <= max_freq] - psd_list.append([psd, px, py, pz]) - - return np.array(first_freqs[first_freqs <= max_freq]), np.array(psd_list) - - -def compute_speed_profile(speeds, freqs, psd_list): - # Preallocate arrays as psd_list is known and consistent - pwrtot_sum = np.zeros(len(psd_list)) - pwrtot_x = np.zeros(len(psd_list)) - pwrtot_y = np.zeros(len(psd_list)) - pwrtot_z = np.zeros(len(psd_list)) - - for i, psd in enumerate(psd_list): - pwrtot_sum[i] = np.trapz(psd[0], freqs) - pwrtot_x[i] = np.trapz(psd[1], freqs) - pwrtot_y[i] = np.trapz(psd[2], freqs) - pwrtot_z[i] = np.trapz(psd[3], freqs) - - # Resample the signals to get a better detection of the valleys of low energy - # and avoid getting limited by the speed increment defined by the user - resampled_speeds, resampled_power_sum = resample_signal(speeds, pwrtot_sum) - _, resampled_pwrtot_x = resample_signal(speeds, pwrtot_x) - _, resampled_pwrtot_y = resample_signal(speeds, pwrtot_y) - _, resampled_pwrtot_z = resample_signal(speeds, pwrtot_z) - - return resampled_speeds, [resampled_power_sum, resampled_pwrtot_x, resampled_pwrtot_y, resampled_pwrtot_z] - - -def compute_motor_profile(power_spectral_densities): - # Sum the PSD across all speeds for each frequency of the spectrogram. Basically this - # is equivalent to sum up all the spectrogram column by column to plot the total on the right - motor_total_vibration = np.sum([psd[0] for psd in power_spectral_densities], axis=0) - - # Then a very little smoothing of the signal is applied to avoid too much noise and sharp peaks on it and simplify - # the resonance frequency and damping ratio estimation later on. Also, too much smoothing is bad and would alter the results - smoothed_motor_total_vibration = np.convolve(motor_total_vibration, np.ones(10)/10, mode='same') - - return smoothed_motor_total_vibration - - -# The goal is to find zone outside of peaks (flat low energy zones) to advise them as good speeds range to use in the slicer -def identify_low_energy_zones(power_total): - valleys = [] - - # Calculate the mean and standard deviation of the entire power_total - mean_energy = np.mean(power_total) - std_energy = np.std(power_total) - - # Define a threshold value as mean minus a certain number of standard deviations - threshold_value = mean_energy - VALLEY_DETECTION_THRESHOLD * std_energy - - # Find valleys in power_total based on the threshold - in_valley = False - start_idx = 0 - for i, value in enumerate(power_total): - if not in_valley and value < threshold_value: - in_valley = True - start_idx = i - elif in_valley and value >= threshold_value: - in_valley = False - valleys.append((start_idx, i)) - - # If the last point is still in a valley, close the valley - if in_valley: - valleys.append((start_idx, len(power_total) - 1)) - - max_signal = np.max(power_total) - - # Calculate mean energy for each valley as a percentage of the maximum of the signal - valley_means_percentage = [] - for start, end in valleys: - if not np.isnan(np.mean(power_total[start:end])): - valley_means_percentage.append((start, end, (np.mean(power_total[start:end]) / max_signal) * 100)) - - # Sort valleys based on mean percentage values - sorted_valleys = sorted(valley_means_percentage, key=lambda x: x[2]) - - return sorted_valleys - - -# Resample the signal to achieve denser data points in order to get more precise valley placing and -# avoid having to use the original sampling of the signal (that is equal to the speed increment used for the test) -def resample_signal(speeds, power_total, new_spacing=0.1): - new_speeds = np.arange(speeds[0], speeds[-1] + new_spacing, new_spacing) - new_power_total = np.interp(new_speeds, speeds, power_total) - return np.array(new_speeds), np.array(new_power_total) - - -###################################################################### -# Graphing -###################################################################### - -def plot_speed_profile(ax, speeds, power_total, num_peaks, peaks, low_energy_zones): - # For this function, we have two array for the speeds. Indeed, since the power total sum was resampled to better detect - # the valleys of low energy later on, we also need the resampled speed array to plot it. For the rest - ax.set_title("Machine speed profile", fontsize=14, color=KLIPPAIN_COLORS['dark_orange'], weight='bold') - ax.set_xlabel('Speed (mm/s)') - ax.set_ylabel('Energy') - - ax2 = ax.twinx() - ax2.yaxis.set_visible(False) - - max_y = power_total[0].max() + power_total[0].max() * 0.05 - ax.set_xlim([speeds.min(), speeds.max()]) - ax.set_ylim([0, max_y]) - ax2.set_ylim([0, max_y]) - - ax.plot(speeds, power_total[0], label="X+Y+Z", color='purple', zorder=5) - ax.plot(speeds, power_total[1], label="X", color='red') - ax.plot(speeds, power_total[2], label="Y", color='green') - ax.plot(speeds, power_total[3], label="Z", color='blue') - - if peaks.size: - ax.plot(speeds[peaks], power_total[0][peaks], "x", color='black', markersize=8) - for idx, peak in enumerate(peaks): - fontcolor = 'red' - fontweight = 'bold' - ax.annotate(f"{idx+1}", (speeds[peak], power_total[0][peak]), - textcoords="offset points", xytext=(8, 5), - ha='left', fontsize=13, color=fontcolor, weight=fontweight) - ax2.plot([], [], ' ', label=f'Number of peaks: {num_peaks}') - else: - ax2.plot([], [], ' ', label=f'No peaks detected') - - for idx, (start, end, energy) in enumerate(low_energy_zones): - ax.axvline(speeds[start], color='red', linestyle='dotted', linewidth=1.5) - ax.axvline(speeds[end], color='red', linestyle='dotted', linewidth=1.5) - ax2.fill_between(speeds[start:end], 0, power_total[0][start:end], color='green', alpha=0.2, label=f'Zone {idx+1}: {speeds[start]:.1f} to {speeds[end]:.1f} mm/s (mean energy: {energy:.2f}%)') - - ax.xaxis.set_minor_locator(matplotlib.ticker.AutoMinorLocator()) - ax.yaxis.set_minor_locator(matplotlib.ticker.AutoMinorLocator()) - ax.grid(which='major', color='grey') - ax.grid(which='minor', color='lightgrey') - fontP = matplotlib.font_manager.FontProperties() - fontP.set_size('small') - ax.legend(loc='upper left', prop=fontP) - ax2.legend(loc='upper right', prop=fontP) - - return - - -def plot_vibration_spectrogram(ax, speeds, freqs, power_spectral_densities, peaks, fr, max_freq): - # Prepare the spectrum data - spectrum = np.empty([len(freqs), len(speeds)]) - for i in range(len(speeds)): - for j in range(len(freqs)): - spectrum[j, i] = power_spectral_densities[i][0][j] - - ax.set_title("Vibrations spectrogram", fontsize=14, color=KLIPPAIN_COLORS['dark_orange'], weight='bold') - # ax.pcolormesh(speeds, freqs, spectrum, norm=matplotlib.colors.LogNorm(), - # cmap='inferno', shading='gouraud') - - ax.imshow(spectrum, norm=matplotlib.colors.LogNorm(), cmap='inferno', - aspect='auto', extent=[speeds[0], speeds[-1], freqs[0], freqs[-1]], - origin='lower', interpolation='antialiased') - - # Add peaks lines in the spectrogram to get hint from peaks found in the first graph - if peaks is not None: - for idx, peak in enumerate(peaks): - ax.axvline(peak, color='cyan', linestyle='dotted', linewidth=0.75) - ax.annotate(f"Peak {idx+1}", (peak, freqs[-1]*0.9), - textcoords="data", color='cyan', rotation=90, fontsize=10, - verticalalignment='top', horizontalalignment='right') - - # Add motor resonance line - if fr is not None and fr > 25: - ax.axhline(fr, color='cyan', linestyle='dotted', linewidth=1) - ax.annotate(f"Motor resonance", (speeds[-1]*0.95, fr+2), - textcoords="data", color='cyan', fontsize=10, - verticalalignment='bottom', horizontalalignment='right') - - ax.set_ylim([0., max_freq]) - ax.set_ylabel('Frequency (hz)') - ax.set_xlabel('Speed (mm/s)') - - return - - -def plot_motor_profile(ax, freqs, motor_vibration_power, motor_fr, motor_zeta, motor_max_power_index): - ax.set_title("Motors frequency profile", fontsize=14, color=KLIPPAIN_COLORS['dark_orange'], weight='bold') - ax.set_xlabel('Energy') - ax.set_ylabel('Frequency (hz)') - - ax2 = ax.twinx() - ax2.yaxis.set_visible(False) - - ax.set_ylim([freqs.min(), freqs.max()]) - ax.set_xlim([0, motor_vibration_power.max() + motor_vibration_power.max() * 0.1]) - - # Plot the profile curve - ax.plot(motor_vibration_power, freqs, color=KLIPPAIN_COLORS['orange']) - - # Tag the resonance peak - ax.plot(motor_vibration_power[motor_max_power_index], freqs[motor_max_power_index], "x", color='black', markersize=8) - fontcolor = KLIPPAIN_COLORS['purple'] - fontweight = 'bold' - ax.annotate(f"R", (motor_vibration_power[motor_max_power_index], freqs[motor_max_power_index]), - textcoords="offset points", xytext=(8, 8), - ha='right', fontsize=13, color=fontcolor, weight=fontweight) - - # Add the legend - ax2.plot([], [], ' ', label="Motor resonant frequency (ω0): %.1fHz" % (motor_fr)) - ax2.plot([], [], ' ', label="Motor damping ratio (ζ): %.3f" % (motor_zeta)) - - ax.xaxis.set_minor_locator(matplotlib.ticker.AutoMinorLocator()) - ax.yaxis.set_minor_locator(matplotlib.ticker.AutoMinorLocator()) - ax.grid(which='major', color='grey') - ax.grid(which='minor', color='lightgrey') - fontP = matplotlib.font_manager.FontProperties() - fontP.set_size('small') - ax2.legend(loc='upper right', prop=fontP) - - return - - -###################################################################### -# Startup and main routines -###################################################################### - -def extract_speed(logname): - try: - speed = re.search('sp(.+?)n', os.path.basename(logname)).group(1).replace('_','.') - except AttributeError: - raise ValueError("File %s does not contain speed in its name and therefore " - "is not supported by this script." % (logname,)) - return float(speed) - - -def sort_and_slice(raw_speeds, raw_datas, remove): - # Sort to get the speeds and their datas aligned and in ascending order - raw_speeds, raw_datas = zip(*sorted(zip(raw_speeds, raw_datas), key=operator.itemgetter(0))) - - # Optionally remove the beginning and end of each data file to get only - # the constant speed part of the segments and remove the start/stop phase - sliced_datas = [] - for data in raw_datas: - sliced = round((len(data) * remove / 100) / 2) - sliced_datas.append(data[sliced:len(data)-sliced]) - - return raw_speeds, sliced_datas - - -def vibrations_calibration(lognames, klipperdir="~/klipper", axisname=None, accel=None, max_freq=1000., remove=0): - set_locale() - global shaper_calibrate - shaper_calibrate = setup_klipper_import(klipperdir) - - # Parse the raw data and get them ready for analysis - raw_datas = [parse_log(filename) for filename in lognames] - raw_speeds = [extract_speed(filename) for filename in lognames] - speeds, datas = sort_and_slice(raw_speeds, raw_datas, remove) - del raw_datas, raw_speeds - - # As we assume that we have the same number of file for each speed increment, we can group - # the PSD results by this number (to combine all the segments of the pattern at a constant speed) - group_by = speeds.count(speeds[0]) - - # Remove speeds duplicates and graph the processed datas - speeds = list(OrderedDict((x, True) for x in speeds).keys()) - - # Compute speed profile, vibration spectrogram and motor resonance profile - freqs, psd = compute_vibration_spectrogram(datas, group_by, max_freq) - upsampled_speeds, speeds_powers = compute_speed_profile(speeds, freqs, psd) - motor_vibration_power = compute_motor_profile(psd) - - # Peak detection and low energy valleys (good speeds) identification between the peaks - num_peaks, vibration_peaks, peaks_speeds = detect_peaks( - speeds_powers[0], upsampled_speeds, - PEAKS_DETECTION_THRESHOLD * speeds_powers[0].max(), - PEAKS_RELATIVE_HEIGHT_THRESHOLD, 10, 10 - ) - low_energy_zones = identify_low_energy_zones(speeds_powers[0]) - - # Print the vibration peaks info in the console - formated_peaks_speeds = ["{:.1f}".format(pspeed) for pspeed in peaks_speeds] - print_with_c_locale("Vibrations peaks detected: %d @ %s mm/s (avoid setting a speed near these values in your slicer print profile)" % (num_peaks, ", ".join(map(str, formated_peaks_speeds)))) - - # Motor resonance estimation - motor_fr, motor_zeta, motor_max_power_index = compute_mechanical_parameters(motor_vibration_power, freqs) - if motor_fr > 25: - print_with_c_locale("Motors have a main resonant frequency at %.1fHz with an estimated damping ratio of %.3f" % (motor_fr, motor_zeta)) - else: - print_with_c_locale("The detected resonance frequency of the motors is too low (%.1fHz). This is probably due to the test run with too high acceleration!" % motor_fr) - print_with_c_locale("Try lowering the ACCEL value before restarting the macro to ensure that only constant speeds are recorded and that the dynamic behavior in the corners is not impacting the measurements.") - - # Create graph layout - fig, ((ax1, ax2), (ax3, ax4)) = plt.subplots(2, 2, gridspec_kw={ - 'height_ratios':[4, 3], - 'width_ratios':[5, 3], - 'bottom':0.050, - 'top':0.890, - 'left':0.057, - 'right':0.985, - 'hspace':0.166, - 'wspace':0.138 - }) - ax2.remove() # top right graph is not used and left blank for now... - fig.set_size_inches(14, 11.6) - - # Add title - title_line1 = "VIBRATIONS MEASUREMENT TOOL" - fig.text(0.075, 0.965, title_line1, ha='left', va='bottom', fontsize=20, color=KLIPPAIN_COLORS['purple'], weight='bold') - try: - filename_parts = (lognames[0].split('/')[-1]).split('_') - dt = datetime.strptime(f"{filename_parts[1]} {filename_parts[2].split('-')[0]}", "%Y%m%d %H%M%S") - title_line2 = dt.strftime('%x %X') - if axisname is not None: - title_line2 += ' -- ' + str(axisname).upper() + ' axis' - if accel is not None: - title_line2 += ' at ' + str(accel) + ' mm/s²' - except: - print_with_c_locale("Warning: CSV filename look to be different than expected (%s)" % (lognames[0])) - title_line2 = lognames[0].split('/')[-1] - fig.text(0.075, 0.957, title_line2, ha='left', va='top', fontsize=16, color=KLIPPAIN_COLORS['dark_purple']) - - # Plot the graphs - plot_speed_profile(ax1, upsampled_speeds, speeds_powers, num_peaks, vibration_peaks, low_energy_zones) - plot_motor_profile(ax4, freqs, motor_vibration_power, motor_fr, motor_zeta, motor_max_power_index) - plot_vibration_spectrogram(ax3, speeds, freqs, psd, peaks_speeds, motor_fr, max_freq) - - # Adding a small Klippain logo to the top left corner of the figure - ax_logo = fig.add_axes([0.001, 0.924, 0.075, 0.075], anchor='NW') - ax_logo.imshow(plt.imread(os.path.join(os.path.dirname(os.path.abspath(__file__)), 'klippain.png'))) - ax_logo.axis('off') - - # Adding Shake&Tune version in the top right corner - st_version = get_git_version() - if st_version is not None: - fig.text(0.995, 0.985, st_version, ha='right', va='bottom', fontsize=8, color=KLIPPAIN_COLORS['purple']) - - return fig - - -def main(): - # Parse command-line arguments - usage = "%prog [options] " - opts = optparse.OptionParser(usage) - opts.add_option("-o", "--output", type="string", dest="output", - default=None, help="filename of output graph") - opts.add_option("-a", "--axis", type="string", dest="axisname", - default=None, help="axis name to be printed on the graph") - opts.add_option("-c", "--accel", type="int", dest="accel", - default=None, help="accel value to be printed on the graph") - opts.add_option("-f", "--max_freq", type="float", default=1000., - help="maximum frequency to graph") - opts.add_option("-r", "--remove", type="int", default=0, - help="percentage of data removed at start/end of each CSV files") - opts.add_option("-k", "--klipper_dir", type="string", dest="klipperdir", - default="~/klipper", help="main klipper directory") - options, args = opts.parse_args() - if len(args) < 1: - opts.error("No CSV file(s) to analyse") - if options.output is None: - opts.error("You must specify an output file.png to use the script (option -o)") - if options.remove > 50 or options.remove < 0: - opts.error("You must specify a correct percentage (option -r) in the 0-50 range") - - fig = vibrations_calibration(args, options.klipperdir, options.axisname, options.accel, options.max_freq, options.remove) - fig.savefig(options.output, dpi=150) - - -if __name__ == '__main__': - main() diff --git a/K-ShakeTune/scripts/is_workflow.py b/K-ShakeTune/scripts/is_workflow.py deleted file mode 100755 index f072c63..0000000 --- a/K-ShakeTune/scripts/is_workflow.py +++ /dev/null @@ -1,303 +0,0 @@ -#!/usr/bin/env python3 - -############################################ -###### INPUT SHAPER KLIPPAIN WORKFLOW ###### -############################################ -# Written by Frix_x#0161 # - -# This script is designed to be used with gcode_shell_commands directly from Klipper -# Use the provided Shake&Tune macros instead! - - -import optparse -import os -import time -import glob -import sys -import shutil -import tarfile -from datetime import datetime - -################################################################################################################# -RESULTS_FOLDER = os.path.expanduser('~/printer_data/config/K-ShakeTune_results') -KLIPPER_FOLDER = os.path.expanduser('~/klipper') -################################################################################################################# - -from graph_belts import belts_calibration -from graph_shaper import shaper_calibration -from graph_vibrations import vibrations_calibration -from analyze_axesmap import axesmap_calibration - -RESULTS_SUBFOLDERS = ['belts', 'inputshaper', 'vibrations'] - - -def is_file_open(filepath): - for proc in os.listdir('/proc'): - if proc.isdigit(): - for fd in glob.glob(f'/proc/{proc}/fd/*'): - try: - if os.path.samefile(fd, filepath): - return True - except FileNotFoundError: - # Klipper has already released the CSV file - pass - except PermissionError: - # Unable to check for this particular process due to permissions - pass - return False - - -def create_belts_graph(keep_csv): - current_date = datetime.now().strftime('%Y%m%d_%H%M%S') - lognames = [] - - globbed_files = glob.glob('/tmp/raw_data_axis*.csv') - if not globbed_files: - print("No CSV files found in the /tmp folder to create the belt graphs!") - sys.exit(1) - if len(globbed_files) < 2: - print("Not enough CSV files found in the /tmp folder. Two files are required for the belt graphs!") - sys.exit(1) - - sorted_files = sorted(globbed_files, key=os.path.getmtime, reverse=True) - - for filename in sorted_files[:2]: - # Wait for the file handler to be released by Klipper - while is_file_open(filename): - time.sleep(2) - - # Extract the tested belt from the filename and rename/move the CSV file to the result folder - belt = os.path.basename(filename).split('_')[3].split('.')[0].upper() - new_file = os.path.join(RESULTS_FOLDER, RESULTS_SUBFOLDERS[0], f'belt_{current_date}_{belt}.csv') - shutil.move(filename, new_file) - os.sync() # Sync filesystem to avoid problems - - # Save the file path for later - lognames.append(new_file) - - # Wait for the file handler to be released by the move command - while is_file_open(new_file): - time.sleep(2) - - # Generate the belts graph and its name - fig = belts_calibration(lognames, KLIPPER_FOLDER) - png_filename = os.path.join(RESULTS_FOLDER, RESULTS_SUBFOLDERS[0], f'belts_{current_date}.png') - fig.savefig(png_filename, dpi=150) - - # Remove the CSV files if the user don't want to keep them - if not keep_csv: - for csv in lognames: - if os.path.exists(csv): - os.remove(csv) - - return - - -def create_shaper_graph(keep_csv, max_smoothing, scv): - current_date = datetime.now().strftime('%Y%m%d_%H%M%S') - - # Get all the files and sort them based on last modified time to select the most recent one - globbed_files = glob.glob('/tmp/raw_data*.csv') - if not globbed_files: - print("No CSV files found in the /tmp folder to create the input shaper graphs!") - sys.exit(1) - - sorted_files = sorted(globbed_files, key=os.path.getmtime, reverse=True) - filename = sorted_files[0] - - # Wait for the file handler to be released by Klipper - while is_file_open(filename): - time.sleep(2) - - # Extract the tested axis from the filename and rename/move the CSV file to the result folder - axis = os.path.basename(filename).split('_')[3].split('.')[0].upper() - new_file = os.path.join(RESULTS_FOLDER, RESULTS_SUBFOLDERS[1], f'resonances_{current_date}_{axis}.csv') - shutil.move(filename, new_file) - os.sync() # Sync filesystem to avoid problems - - # Wait for the file handler to be released by the move command - while is_file_open(new_file): - time.sleep(2) - - # Generate the shaper graph and its name - fig = shaper_calibration([new_file], KLIPPER_FOLDER, max_smoothing=max_smoothing, scv=scv) - png_filename = os.path.join(RESULTS_FOLDER, RESULTS_SUBFOLDERS[1], f'resonances_{current_date}_{axis}.png') - fig.savefig(png_filename, dpi=150) - - # Remove the CSV file if the user don't want to keep it - if not keep_csv: - if os.path.exists(new_file): - os.remove(new_file) - - return axis - - -def create_vibrations_graph(axis_name, accel, chip_name, keep_csv): - current_date = datetime.now().strftime('%Y%m%d_%H%M%S') - lognames = [] - - globbed_files = glob.glob(f'/tmp/{chip_name}-*.csv') - if not globbed_files: - print("No CSV files found in the /tmp folder to create the vibration graphs!") - sys.exit(1) - if len(globbed_files) < 3: - print("Not enough CSV files found in the /tmp folder. At least 3 files are required for the vibration graphs!") - sys.exit(1) - - for filename in globbed_files: - # Wait for the file handler to be released by Klipper - while is_file_open(filename): - time.sleep(2) - - # Cleanup of the filename and moving it in the result folder - cleanfilename = os.path.basename(filename).replace(chip_name, f'vibr_{current_date}') - new_file = os.path.join(RESULTS_FOLDER, RESULTS_SUBFOLDERS[2], cleanfilename) - shutil.move(filename, new_file) - - # Save the file path for later - lognames.append(new_file) - - # Sync filesystem to avoid problems as there is a lot of file copied - os.sync() - time.sleep(5) - - # Generate the vibration graph and its name - fig = vibrations_calibration(lognames, KLIPPER_FOLDER, axis_name, accel) - png_filename = os.path.join(RESULTS_FOLDER, RESULTS_SUBFOLDERS[2], f'vibrations_{current_date}_{axis_name}.png') - fig.savefig(png_filename, dpi=150) - - # Archive all the csv files in a tarball in case the user want to keep them - if keep_csv: - with tarfile.open(os.path.join(RESULTS_FOLDER, RESULTS_SUBFOLDERS[2], f'vibrations_{current_date}_{axis_name}.tar.gz'), 'w:gz') as tar: - for csv_file in lognames: - tar.add(csv_file, recursive=False) - - # Remove the remaining CSV files not needed anymore (tarball is safe if it was created) - for csv_file in lognames: - if os.path.exists(csv_file): - os.remove(csv_file) - - return - - -def find_axesmap(accel, chip_name): - current_date = datetime.now().strftime('%Y%m%d_%H%M%S') - result_filename = os.path.join(RESULTS_FOLDER, f'axes_map_{current_date}.txt') - lognames = [] - - globbed_files = glob.glob(f'/tmp/{chip_name}-*.csv') - if not globbed_files: - print("No CSV files found in the /tmp folder to analyze and find the axes_map!") - sys.exit(1) - - sorted_files = sorted(globbed_files, key=os.path.getmtime, reverse=True) - filename = sorted_files[0] - - # Wait for the file handler to be released by Klipper - while is_file_open(filename): - time.sleep(2) - - # Analyze the CSV to find the axes_map parameter - lognames.append(filename) - results = axesmap_calibration(lognames, accel) - - with open(result_filename, 'w') as f: - f.write(results) - - return - - -# Utility function to get old files based on their modification time -def get_old_files(folder, extension, limit): - files = [os.path.join(folder, f) for f in os.listdir(folder) if f.endswith(extension)] - files.sort(key=lambda x: os.path.getmtime(x), reverse=True) - return files[limit:] - -def clean_files(keep_results): - # Define limits based on STORE_RESULTS - keep1 = keep_results + 1 - keep2 = 2 * keep_results + 1 - - # Find old files in each directory - old_belts_files = get_old_files(os.path.join(RESULTS_FOLDER, RESULTS_SUBFOLDERS[0]), '.png', keep1) - old_inputshaper_files = get_old_files(os.path.join(RESULTS_FOLDER, RESULTS_SUBFOLDERS[1]), '.png', keep2) - old_vibrations_files = get_old_files(os.path.join(RESULTS_FOLDER, RESULTS_SUBFOLDERS[2]), '.png', keep1) - - # Remove the old belt files - for old_file in old_belts_files: - file_date = "_".join(os.path.splitext(os.path.basename(old_file))[0].split('_')[1:3]) - for suffix in ['A', 'B']: - csv_file = os.path.join(RESULTS_FOLDER, RESULTS_SUBFOLDERS[0], f'belt_{file_date}_{suffix}.csv') - if os.path.exists(csv_file): - os.remove(csv_file) - os.remove(old_file) - - # Remove the old shaper files - for old_file in old_inputshaper_files: - csv_file = os.path.join(RESULTS_FOLDER, RESULTS_SUBFOLDERS[1], os.path.splitext(os.path.basename(old_file))[0] + ".csv") - if os.path.exists(csv_file): - os.remove(csv_file) - os.remove(old_file) - - # Remove the old vibrations files - for old_file in old_vibrations_files: - os.remove(old_file) - tar_file = os.path.join(RESULTS_FOLDER, RESULTS_SUBFOLDERS[2], os.path.splitext(os.path.basename(old_file))[0] + ".tar.gz") - if os.path.exists(tar_file): - os.remove(tar_file) - - -def main(): - # Parse command-line arguments - usage = "%prog [options] " - opts = optparse.OptionParser(usage) - opts.add_option("-t", "--type", type="string", dest="type", - default=None, help="type of output graph to produce") - opts.add_option("--accel", type="int", default=None, dest="accel_used", - help="acceleration used during the vibration macro or axesmap macro") - opts.add_option("--axis_name", type="string", default=None, dest="axis_name", - help="axis tested during the vibration macro") - opts.add_option("--chip_name", type="string", default="adxl345", dest="chip_name", - help="accelerometer chip name in klipper used during the vibration macro or the axesmap macro") - opts.add_option("-n", "--keep_results", type="int", default=3, dest="keep_results", - help="number of results to keep in the result folder after each run of the script") - opts.add_option("-c", "--keep_csv", action="store_true", default=False, dest="keep_csv", - help="weither or not to keep the CSV files alongside the PNG graphs image results") - opts.add_option("--scv", "--square_corner_velocity", type="float", dest="scv", default=5., - help="square corner velocity used to compute max accel for axis shapers graphs") - opts.add_option("--max_smoothing", type="float", dest="max_smoothing", default=None, - help="maximum shaper smoothing to allow") - options, args = opts.parse_args() - - if options.type is None: - opts.error("You must specify the type of output graph you want to produce (option -t)") - elif options.type.lower() is None or options.type.lower() not in ['belts', 'shaper', 'vibrations', 'axesmap', 'clean']: - opts.error("Type of output graph need to be in the list of 'belts', 'shaper', 'vibrations', 'axesmap' or 'clean'") - else: - graph_mode = options.type - - # Check if results folders are there or create them before doing anything else - for result_subfolder in RESULTS_SUBFOLDERS: - folder = os.path.join(RESULTS_FOLDER, result_subfolder) - if not os.path.exists(folder): - os.makedirs(folder) - - if graph_mode.lower() == 'belts': - create_belts_graph(keep_csv=options.keep_csv) - print(f"Belt graph created. You will find the results in {RESULTS_FOLDER}/{RESULTS_SUBFOLDERS[0]}") - elif graph_mode.lower() == 'shaper': - axis = create_shaper_graph(keep_csv=options.keep_csv, max_smoothing=options.max_smoothing, scv=options.scv) - print(f"{axis} input shaper graph created. You will find the results in {RESULTS_FOLDER}/{RESULTS_SUBFOLDERS[1]}") - elif graph_mode.lower() == 'vibrations': - create_vibrations_graph(axis_name=options.axis_name, accel=options.accel_used, chip_name=options.chip_name, keep_csv=options.keep_csv) - print(f"{options.axis_name} vibration graph created. You will find the results in {RESULTS_FOLDER}/{RESULTS_SUBFOLDERS[2]}") - elif graph_mode.lower() == 'axesmap': - print(f"WARNING: AXES_MAP_CALIBRATION is currently very experimental and may produce incorrect results... Please validate the output!") - find_axesmap(accel=options.accel_used, chip_name=options.chip_name) - elif graph_mode.lower() == 'clean': - print(f"Cleaning output folder to keep only the last {options.keep_results} results...") - clean_files(keep_results=options.keep_results) - - -if __name__ == '__main__': - main() diff --git a/K-ShakeTune/scripts/shaketune.sh b/K-ShakeTune/scripts/shaketune.sh deleted file mode 100755 index 952a161..0000000 --- a/K-ShakeTune/scripts/shaketune.sh +++ /dev/null @@ -1,5 +0,0 @@ -#!/usr/bin/env bash - -source ~/klippain_shaketune-env/bin/activate -python ~/klippain_shaketune/K-ShakeTune/scripts/is_workflow.py "$@" -deactivate diff --git a/K-ShakeTune/shaketune.sh b/K-ShakeTune/shaketune.sh new file mode 100755 index 0000000..53af59f --- /dev/null +++ b/K-ShakeTune/shaketune.sh @@ -0,0 +1,10 @@ +#!/usr/bin/env bash + +# This script is used to run the Shake&Tune Python scripts as a module +# from the project root directory using its virtual environment +# Usage: ./shaketune.sh + +source ~/klippain_shaketune-env/bin/activate +cd ~/klippain_shaketune +python -m src.is_workflow "$@" +deactivate diff --git a/K-ShakeTune/shaketune_cmd.cfg b/K-ShakeTune/shaketune_cmd.cfg index e4d667c..8891eda 100644 --- a/K-ShakeTune/shaketune_cmd.cfg +++ b/K-ShakeTune/shaketune_cmd.cfg @@ -1,5 +1,5 @@ [gcode_shell_command shaketune] -command: ~/printer_data/config/K-ShakeTune/scripts/shaketune.sh +command: ~/printer_data/config/K-ShakeTune/shaketune.sh timeout: 600.0 verbose: True diff --git a/README.md b/README.md index bfdad50..0fa696c 100644 --- a/README.md +++ b/README.md @@ -13,9 +13,9 @@ It operates in two steps: Check out the **[detailed documentation of the Shake&Tune module here](./docs/README.md)**. You can also look at the documentation for each type of graph by directly clicking on them below to better understand your results and tune your machine! -| [Belts graph](./docs/macros/belts_tuning.md) | [Axis input shaper graphs](./docs/macros/axis_tuning.md) | [Vibrations graph](./docs/macros/vibrations_tuning.md) | +| [Belts graph](./docs/macros/belts_tuning.md) | [Axis input shaper graphs](./docs/macros/axis_tuning.md) | [Vibrations graph](./docs/macros/vibrations_profile.md) | |:----------------:|:------------:|:---------------------:| -| [](./docs/macros/belts_tuning.md) | [](./docs/macros/axis_tuning.md) | [](./docs/macros/vibrations_tuning.md) | +| [](./docs/macros/belts_tuning.md) | [](./docs/macros/axis_tuning.md) | [](./docs/macros/vibrations_profile.md) | > **Note**: > @@ -34,17 +34,13 @@ Follow these steps to install the Shake&Tune module in your printer: [include K-ShakeTune/*.cfg] ``` - > **Note**: - > - > Due to some breaking changes in the resonance testing code on the Klipper side, Shake&Tune has been modified to take advantage of this, and thus S&T v2.6+ will only support a Klipper version from Feb 17th 2024. If you are using an older version of Klipper, you must use S&T <=2.5.x - ## Usage Ensure your machine is homed, then invoke one of the following macros as needed: - - `AXES_MAP_CALIBRATION` to automatically find Klipper's `axes_map` parameter for your accelerometer orientation (be careful, this is experimental for now). - - `BELTS_SHAPER_CALIBRATION` for belt resonance graphs, useful for verifying belt tension and differential belt paths behavior. - - `AXES_SHAPER_CALIBRATION` for input shaper graphs to mitigate ringing/ghosting by tuning Klipper's input shaper system. - - `VIBRATIONS_CALIBRATION` for machine and motors vibration graphs, used to optimize your slicer speed profiles and TMC drivers parameters. - - `EXCITATE_AXIS_AT_FREQ` to sustain a specific excitation frequency, useful to let you inspect and find out what is resonating. + - `AXES_MAP_CALIBRATION` to automatically find Klipper's `axes_map` parameter for your accelerometer orientation (be careful, this is experimental for now and known to give bad results). + - `COMPARE_BELTS_RESPONSES` for a differential belt resonance graph, useful for checking relative belt tensions and belt path behaviors on a CoreXY printer. + - `AXES_SHAPER_CALIBRATION` for standard input shaper graphs, used to mitigate ringing/ghosting by tuning Klipper's input shaper filters. + - `CREATE_VIBRATIONS_PROFILE` for vibrations graphs as a function of toolhead direction and speed, used to find problematic ranges where the printer could be exposed to more VFAs and optimize your slicer speed profiles and TMC driver parameters. + - `EXCITATE_AXIS_AT_FREQ` to maintain a specific excitation frequency, useful to inspect and find out what is resonating. For further insights on the usage of these macros and the generated graphs, refer to the [K-Shake&Tune module documentation](./docs/README.md). diff --git a/docs/README.md b/docs/README.md index cb55236..d0f8f87 100644 --- a/docs/README.md +++ b/docs/README.md @@ -8,14 +8,14 @@ First, check out the **[input shaping and tuning generalities](./is_tuning_gener Then look at the documentation for each type of graph by clicking on them below tu run the tests and better understand your results to tune your machine! -| [Belts graph](./macros/belts_tuning.md) | [Axis input shaper graphs](./macros/axis_tuning.md) | [Vibrations graph](./macros/vibrations_tuning.md) | +| [Belt response comparison](./macros/belts_tuning.md) | [Axis input shaper graphs](./macros/axis_tuning.md) | [Vibrations profile](./macros/vibrations_profile.md) | |:----------------:|:------------:|:---------------------:| -| [](./macros/belts_tuning.md) | [](./macros/axis_tuning.md) | [](./macros/vibrations_tuning.md) | +| [](./macros/belts_tuning.md) | [](./macros/axis_tuning.md) | [](./macros/vibrations_profile.md) | ## Additional macros -### AXES_MAP_CALIBRATION +### AXES_MAP_CALIBRATION (experimental) All graphs generated by this package show plots based on accelerometer measurements, typically labeled with the X, Y, and Z axes. It's important to note that if the accelerometer is rotated, its axes may not align correctly with the machine axes, making the plots more difficult to interpret, analyze, and understand. The `AXES_MAP_CALIBRATION` is designed to automatically measure the alignement of the accelerometer in order to set it correctly. diff --git a/docs/images/vibrations_example.png b/docs/images/vibrations_example.png index 03c9406..ec9059a 100644 Binary files a/docs/images/vibrations_example.png and b/docs/images/vibrations_example.png differ diff --git a/docs/images/vibrations_graphs/angular_speed_energy_profile.png b/docs/images/vibrations_graphs/angular_speed_energy_profile.png new file mode 100644 index 0000000..df12173 Binary files /dev/null and b/docs/images/vibrations_graphs/angular_speed_energy_profile.png differ diff --git a/docs/images/vibrations_graphs/global_speed_energy_profile.png b/docs/images/vibrations_graphs/global_speed_energy_profile.png new file mode 100644 index 0000000..b865a94 Binary files /dev/null and b/docs/images/vibrations_graphs/global_speed_energy_profile.png differ diff --git a/docs/images/vibrations_graphs/motor_frequency_profile.png b/docs/images/vibrations_graphs/motor_frequency_profile.png new file mode 100644 index 0000000..5393a8e Binary files /dev/null and b/docs/images/vibrations_graphs/motor_frequency_profile.png differ diff --git a/docs/images/vibrations_graphs/polar_angle_energy_profile.png b/docs/images/vibrations_graphs/polar_angle_energy_profile.png new file mode 100644 index 0000000..0cac32a Binary files /dev/null and b/docs/images/vibrations_graphs/polar_angle_energy_profile.png differ diff --git a/docs/images/vibrations_graphs/sd2w_spectrogram.png b/docs/images/vibrations_graphs/sd2w_spectrogram.png deleted file mode 100644 index 519177a..0000000 Binary files a/docs/images/vibrations_graphs/sd2w_spectrogram.png and /dev/null differ diff --git a/docs/images/vibrations_graphs/vibration_graph_explanation.png b/docs/images/vibrations_graphs/vibration_graph_explanation.png deleted file mode 100644 index 17b9aee..0000000 Binary files a/docs/images/vibrations_graphs/vibration_graph_explanation.png and /dev/null differ diff --git a/docs/images/vibrations_graphs/vibration_graph_explanation2.png b/docs/images/vibrations_graphs/vibration_graph_explanation2.png deleted file mode 100644 index 0875c52..0000000 Binary files a/docs/images/vibrations_graphs/vibration_graph_explanation2.png and /dev/null differ diff --git a/docs/images/vibrations_graphs/vibrations_heatmaps.png b/docs/images/vibrations_graphs/vibrations_heatmaps.png new file mode 100644 index 0000000..7f58dd2 Binary files /dev/null and b/docs/images/vibrations_graphs/vibrations_heatmaps.png differ diff --git a/docs/macros/axis_tuning.md b/docs/macros/axis_tuning.md index 261feef..a37168b 100644 --- a/docs/macros/axis_tuning.md +++ b/docs/macros/axis_tuning.md @@ -18,7 +18,7 @@ Then, call the `AXES_SHAPER_CALIBRATION` macro and look for the graphs in the re |SCV|printer square corner velocity|Square corner velocity you want to use to calculate shaper recommendations. Using higher SCV values usually results in more smoothing and lower maximum accelerations| |MAX_SMOOTHING|None|Max smoothing allowed when calculating shaper recommendations| |KEEP_N_RESULTS|3|Total number of results to keep in the result folder after running the test. The older results are automatically cleaned up| -|KEEP_CSV|True|Weither or not to keep the CSV data file alonside the PNG graphs| +|KEEP_CSV|0|Weither or not to keep the CSV data file alonside the PNG graphs| ## Graphs description @@ -105,7 +105,7 @@ Here's how to troubleshoot the issue: Such graph patterns can arise from various factors, and there isn't a one-size-fits-all solution. To address them: 1. A wobbly table can be the cause. So first thing to do is to try with the printer directly on the floor. - 1. Ensure optimal belt tension using the [`BELTS_SHAPER_CALIBRATION` macro](./belts_tuning.md). + 1. Ensure optimal belt tension using the [`COMPARE_BELTS_RESPONSES` macro](./belts_tuning.md). 1. If problems persist, it might be due to an improperly squared gantry. For correction, refer to [Nero3D's de-racking video](https://youtu.be/cOn6u9kXvy0?si=ZCSdWU6br3Y9rGsy). 1. If it's still there... you will need to find out what is resonating to fix it. You can use the `EXCITATE_AXIS_AT_FREQ` macro to help you find it. diff --git a/docs/macros/belts_tuning.md b/docs/macros/belts_tuning.md index 6ee3d46..225ed00 100644 --- a/docs/macros/belts_tuning.md +++ b/docs/macros/belts_tuning.md @@ -1,13 +1,13 @@ # Belt relative difference measurements -The `BELTS_SHAPER_CALIBRATION` macro is dedicated for CoreXY machines where it can help you to diagnose belt path problems by measuring and plotting the differences between their behavior. It will also help you tension your belts at the same tension. +The `COMPARE_BELTS_RESPONSES` macro is dedicated for CoreXY machines where it can help you to diagnose belt path problems by measuring and plotting the differences between their behavior. It will also help you tension your belts at the same tension. ## Usage **Before starting, ensure that the belts are properly tensioned**. For example, you can follow the [Voron belt tensioning documentation](https://docs.vorondesign.com/tuning/secondary_printer_tuning.html#belt-tension). This is crucial: you need a good starting point to then iterate from it! -Then, call the `BELTS_SHAPER_CALIBRATION` macro and look for the graphs in the results folder. Here are the parameters available: +Then, call the `COMPARE_BELTS_RESPONSES` macro and look for the graphs in the results folder. Here are the parameters available: | parameters | default value | description | |-----------:|---------------|-------------| @@ -15,7 +15,7 @@ Then, call the `BELTS_SHAPER_CALIBRATION` macro and look for the graphs in the r |FREQ_END|133|Maximum excitation frequency| |HZ_PER_SEC|1|Number of Hz per seconds for the test| |KEEP_N_RESULTS|3|Total number of results to keep in the result folder after running the test. The older results are automatically cleaned up| -|KEEP_CSV|True|Weither or not to keep the CSV data files alonside the PNG graphs| +|KEEP_CSV|0|Weither or not to keep the CSV data files alonside the PNG graphs| ## Graphs description diff --git a/docs/macros/vibrations_profile.md b/docs/macros/vibrations_profile.md new file mode 100644 index 0000000..d11aeb0 --- /dev/null +++ b/docs/macros/vibrations_profile.md @@ -0,0 +1,99 @@ +# Machine vibrations profiles + +The `CREATE_VIBRATIONS_PROFILE` macro analyzes accelerometer data to plot the vibration profile of your 3D printer. The resulting graphs highlight optimal print speeds and angles that produce the least amount of vibration. It provides a technical basis for adjustments in your slicer profiles, but also in hardware setup and TMC driver parameters to improve print quality and reduce VFAs (vertical fines artifacts). + + > **Warning** + > + > You will need to calibrate the standard input shaper algorithms of Klipper using the other macros first! This test should be used as a last step to calibrate your printer with Shake&Tune. + + +## Usage + +Call the `CREATE_VIBRATIONS_PROFILE` macro with the speed range you want to measure. Here are the parameters available: + +| parameters | default value | description | +|-----------:|---------------|-------------| +|SIZE|100|maximum size in mm of the circle in which the recorded movements take place| +|Z_HEIGHT|20|z height to put the toolhead before starting the movements. Be careful, if your accelerometer is mounted under the nozzle, increase it to avoid crashing it on the bed of the machine| +|ACCEL|3000 (or max printer accel)|accel in mm/s^2 used for all moves. Try to keep it relatively low to avoid dynamic effects that alter the measurements, but high enough to achieve a constant speed for >~70% of the segments. 3000 is a reasonable default for most printers, unless you want to record at very high speed, in which case you will want to increase SIZE and decrease ACCEL a bit.| +|MAX_SPEED|200|maximum speed of the toolhead in mm/s to record for analysis| +|SPEED_INCREMENT|2|toolhead speed increments in mm/s between each movement| +|TRAVEL_SPEED|200|speed in mm/s used for all the travels moves| +|ACCEL_CHIP|"adxl345"|accelerometer chip name in the config| +|KEEP_N_RESULTS|3|Total number of results to keep in the result folder after running the test. The older results are automatically cleaned up| +|KEEP_CSV|0|Weither or not to keep the CSV data files alonside the PNG graphs (archived in a tarball)| + + +## Graphs description + +The `CREATE_VIBRATIONS_PROFILE` macro results are constituted of a set of 6 plots. At the top of the figure you can also see all the detected motor, current and TMC driver parameters. These notes are just for reference in case you want to tinker with them and don't forget what you changed between each run of the macro. + +![](../images/vibrations_example.png) + +### Global Speed Energy Profile + +| Example | description | +|:-----|-------------| +|![](../images/vibrations_graphs/global_speed_energy_profile.png)|This plot shows the relationship between toolhead speed (mm/s) and vibrational energy, providing a global view of how speed impacts vibration across all movements. By using speeds from the green zones, your printer will run more smoothly and you will minimize vibrations and related fine artifacts in prints| + +This graph is the most important one of this tool. You want to use it to adapt your slicer profile, especially by looking at the "vibration metric" curve, that will helps you find which speeds can be problematic for your printer. Here's the magic behind it, broken down into two key parts: + 1. **Spectrum Variance**: This is like the mood ring of your printer, showing how the vibes (a.k.a vibrations) change when printing from different angles. If the "vibration metric" is low, it means your printer is keeping its cool, staying consistent no matter the angle. But if it spikes, it's a sign that some angles are making your printer jitter more than a caffeinated squirrel. *Imagine it like this: You're looking for a chill party vibe where the music's good at every angle, not one where you turn a corner and suddenly it's too loud or too soft.* + 2. **Spectrum Max**: This one's about the max volume of the party, or how loud the strongest vibration is across all angles at any speed. We're aiming to avoid the speeds that crank up the volume too high, causing a resonance rave in the motors. *Think of it this way: You don't want the base so high that it feels like your heart's going to beat out of your chest. We're looking for a nice background level where everyone can chat and have a good time.* + +And why do we care so much about finding these speeds? Because during a print, the toolhead will move in all directions depending on the geometry, and we want a speed that's like a good friend, reliable no matter what the situation. Fortunately, since the motors in our printers share their vibes without non-linear mixing and just add up (think of it as each doing its own dance without bumping into each other), we can find those happy green zones on the graph: these are the speeds that keep the vibe cool and the energy just right, making them perfect for all your print jobs. + +### Polar Angle Energy Profile + +| Example | description | +|:-----|-------------| +|![](../images/vibrations_graphs/polar_angle_energy_profile.png)|Shows how vibrational energy varies with the direction where the toolhead is running. It helps in identifying angles that produce less vibration, and potentially detect assymetries in the belt paths for a CoreXY printer| + +This plot is like your go-to playlist for finding those angles where the vibe is just right. But here's the thing: when printing, your toolhead will groove in all directions and angles, depending on the geometry of your parts, so sticking to just one angle isn't possible. My tip to make the most of this chart for your prints: if you're working on something rectangular, try to align it so that most of the edges match the angles that's least likely to make your printer jitter. For those sleek CoreXY printers, aiming for 45/135 degrees is usually a hit, while the trusty Cartesian printers groove best at 0/90 degrees. And for everything else? Well, there's not much more to do here except rely on the [Global Speed Energy Profile chart](#global-speed-energy-profile) to tune your slicer profile speeds instead. + +Now, onto the symmetry indicator. Think of this tool as the dance coach for your printer, especially designed for those with a symmetrical setup like CoreXY models. It's all about using some pretty neat math (cross-correlation, to be exact) to check out the vibes from both sides of the dance floor. Picture it as a top-notch party dancer, scanning the room at every angle, judging each dancer, and only giving top marks when everyone is perfectly in sync. This tool is ace at catching any sneakiness in your motor control or belt path, highlighting any "butterfly" shapes or even the slightest variations in the motors' resonance patterns. It's like having a magnifying glass that points out exactly where the party fouls are, helping you to fix them and keep your prints rolling out smooth and stunning. + +### Angular Speed Energy Profiles + +| Example | description | +|:-----|-------------| +|![](../images/vibrations_graphs/angular_speed_energy_profile.png)|Provides a detailed view of how energy distribution changes with speed for specific angles. It's useful for fine-tuning speeds for different directions of motion, or for tracking and diagnosing your printer's behavior across the major axes| + +This chart is like a snapshot, capturing the vibe at certain angles of your printing party. But remember, it's just a glimpse into a few specific angles and doesn't fully reveal the whole dance floor where the toolhead moves in every direction, vibing with the unique geometry of your parts. So, think of it as a way to peek into how everyone's grooving in each corner of the party. It's great for a quick check-up to see how the vibe is holding up, but when it comes to setting the rhythm of your slicer speeds, you're going to want to use the [Global Speed Energy Profile chart](#global-speed-energy-profile) instead. + +### Vibrations Heatmaps + +| Example | description | +|:-----|-------------| +|![](../images/vibrations_graphs/vibrations_heatmaps.png)|Both plots provides a comprehensive overview of vibrational energy across speeds and angles. It visually identifies zones of high and low energy, aiding in the comprehensive understanding of the printer motors behavior. It's what is captured by the accelerometer and the base of all the other plots| + +Both heatmaps lay down the vibe of vibrational energy across all speeds and angles, painting a picture of how the beat spreads throughout your printer's dance floor. The polar heatmap gives you a 360-degree whirl of the action, while the regular one lays it out in a classic 2D groove, yet both are vibing to the same tune and showing you where the energy's hot and popping and where it's cool and mellow across your printer's operational range. Think of it as the unique fingerprint of your motor's behavior captured by the accelerometer, it's the raw rhythm of your printer in action. + +Because the scale is both normalized and logarithmic, you're looking for a heatmap (or spectrogram) that has a cool, consistent "orangish" vibe throughout, signaling not so much change over the spectrum with fairly low motor resonances. See areas in your heatmap that swing from deep purple/black to bright white/yellow? That's a sign that your printer motors are hitting high resonances at certain angles and speed combinations that are above the baseline vibrations outside of those areas. But remember, this is just the lay of the land, a snapshot of the scene: tweaking this vibe directly may not be easy, but you can still [play around with the TMC driver parameters](#improving-the-results) to adjust the beats and find a smoother rhythm. + +### Motor Frequency Profile + +| Example | description | +|:-----|-------------| +|![](../images/vibrations_graphs/motor_frequency_profile.png)|Identifies the resonant frequencies of the motors and their damping ratios. Informative for now, but will be used later| + +For now, this graph is purely informational and is a measurement of the motor's natural resonance profile. Think of this plot as a sneak peek at the inner workings of your printer's dance floor. It's not quite ready to hit the main stage for practical use, but just you wait... Keep an eye on this chart as it hints at future remixes where you'll get to play DJ and tweak and tune your printer's performance like never before. + + +## Improving the results + +These graphs essentially depict the behavior of the motor control on your machine. While there isn't much room for easy adjustments to enhance them, most of you should only utilize them to configure your slicer profile to avoid problematic speeds. + +However, if you want to go the rabbit hole, as the data in these graphs largely hinges on the type of motors, their physical characteristic and the way they are controlled by the TMC drivers black magic, there are opportunities for optimization. Tweaking TMC parameters allow to adjust the peaks, enhance machine performance, or diminish overall machine noise. For this process, I recommend to directly use the [Klipper TMC Autotune](https://github.com/andrewmcgr/klipper_tmc_autotune) plugin, which should simplify everything considerably. But keep in mind that it's still an experimental plugin and it's not perfect. + +For individuals inclined to reach the bottom of the rabbit hole and that want to handle this manually, the use of an oscilloscope is mandatory. Majority of the necessary resources are available directly on the Trinamics TMC website: + 1. You should first consult the datasheet specific to your TMC model for guidance on parameter names and their respective uses. + 2. Then to tune the parameters, have a look at the application notes available on their platform, especially [AN001](https://www.trinamic.com/fileadmin/assets/Support/AppNotes/AN001-SpreadCycle.pdf), [AN002](https://www.trinamic.com/fileadmin/assets/Support/AppNotes/AN002-StallGuard2.pdf), [AN003](https://www.trinamic.com/fileadmin/assets/Support/AppNotes/AN003_-_DcStep_Basics_and_Wizard.pdf) and [AN009](https://www.trinamic.com/fileadmin/assets/Support/AppNotes/AN009_Tuning_coolStep.pdf). + 3. For a more comprehensive understanding, you might also want to explore [AN015](https://www.trinamic.com/fileadmin/assets/Support/AppNotes/AN015-StealthChop_Performance.pdf) and [AN021](https://www.trinamic.com/fileadmin/assets/Support/AppNotes/AN021-StealthChop_Performance_comparison_V1.12.pdf ), although they are more geared towards enhancing comprehension than calibration, akin to the TMC datasheet. + +For reference, the default settings used in Klipper are: +``` +#driver_TBL: 2 +#driver_TOFF: 3 +#driver_HEND: 0 +#driver_HSTRT: 5 +``` diff --git a/docs/macros/vibrations_tuning.md b/docs/macros/vibrations_tuning.md deleted file mode 100644 index a414d6e..0000000 --- a/docs/macros/vibrations_tuning.md +++ /dev/null @@ -1,61 +0,0 @@ -# Vibrations measurements - -The `VIBRATIONS_CALIBRATION` macro helps you to identify the speed settings that exacerbate the vibrations of the machine (ie. where the frame and motors resonate badly). This will help you to find the clean speed ranges where the machine is more silent and less prone to vertical fine artifacts on the prints. - - > **Warning** - > - > You will first need to calibrate the standard input shaper algorithm of Klipper using the other macros! This test should not be used before as it would be useless and the results invalid. - - -## Usage - -Call the `VIBRATIONS_CALIBRATION` macro with the direction and speed range you want to measure. Here are the parameters available: - -| parameters | default value | description | -|-----------:|---------------|-------------| -|SIZE|60|size in mm of the area where the movements are done| -|DIRECTION|"XY"|direction vector where you want to do the measurements. Can be set to either "XY", "AB", "ABXY", "A", "B", "X", "Y", "Z", "E"| -|Z_HEIGHT|20|z height to put the toolhead before starting the movements. Be careful, if your accelerometer is mounted under the nozzle, increase it to avoid crashing it on the bed of the machine| -|ACCEL|3000 (or max printer accel)|accel in mm/s^2 used for all the moves. Try to keep it relatively low to avoid bad oscillations that affect the measurements, but but high enough to reach constant speed for >~70% of the segments| -|MIN_SPEED|20|minimum speed of the toolhead in mm/s for the movements| -|MAX_SPEED|200|maximum speed of the toolhead in mm/s for the movements| -|SPEED_INCREMENT|2|speed increments of the toolhead in mm/s between every movements| -|TRAVEL_SPEED|200|speed in mm/s used for all the travels moves| -|ACCEL_CHIP|"adxl345"|accelerometer chip name in the config| -|KEEP_N_RESULTS|3|Total number of results to keep in the result folder after running the test. The older results are automatically cleaned up| -|KEEP_CSV|True|Weither or not to keep the CSV data files alonside the PNG graphs (archived in a tarball)| - - -## Graphs description - -![](../images/vibrations_graphs/vibration_graph_explanation.png) -![](../images/vibrations_graphs/vibration_graph_explanation2.png) - -## Improving the results - -These graphs essentially depict the behavior of the motor control on your machine. While there isn't much room for easy adjustments to enhance them, most of you should only utilize them to configure your slicer profile to avoid problematic speeds. - -However, if you want to go the rabbit hole, as the data in these graphs largely hinges on the type of motors, their physical characteristic and the way they are controlled by the TMC drivers black magic, there are opportunities for optimization. Tweaking TMC parameters allow to adjust the peaks, enhance machine performance, or diminish overall machine noise. For this process, I recommend to directly use the [Klipper TMC Autotune](https://github.com/andrewmcgr/klipper_tmc_autotune) plugin, which should simplify everything considerably. But keep in mind that it's still an experimental plugin and it's not perfect. - -For individuals inclined to reach the bottom of the rabbit hole and that want to handle this manually, the use of an oscilloscope is mandatory. Majority of the necessary resources are available directly on the Trinamics TMC website: - 1. You should first consult the datasheet specific to your TMC model for guidance on parameter names and their respective uses. - 2. Then to tune the parameters, have a look at the application notes available on their platform, especially [AN001](https://www.trinamic.com/fileadmin/assets/Support/AppNotes/AN001-SpreadCycle.pdf), [AN002](https://www.trinamic.com/fileadmin/assets/Support/AppNotes/AN002-StallGuard2.pdf), [AN003](https://www.trinamic.com/fileadmin/assets/Support/AppNotes/AN003_-_DcStep_Basics_and_Wizard.pdf) and [AN009](https://www.trinamic.com/fileadmin/assets/Support/AppNotes/AN009_Tuning_coolStep.pdf). - 3. For a more comprehensive understanding, you might also want to explore [AN015](https://www.trinamic.com/fileadmin/assets/Support/AppNotes/AN015-StealthChop_Performance.pdf) and [AN021](https://www.trinamic.com/fileadmin/assets/Support/AppNotes/AN021-StealthChop_Performance_comparison_V1.12.pdf ), although they are more geared towards enhancing comprehension than calibration, akin to the TMC datasheet. - -For reference, the default settings used in Klipper are: -``` -#driver_TBL: 2 -#driver_TOFF: 3 -#driver_HEND: 0 -#driver_HSTRT: 5 -``` - -### Semi-blank spectrogram (LIS2DW) - -The integration of LIS2DW as a resonance measuring device in Klipper is becoming more and more common, especially because some manufacturers are promoting its superiority over the established ADXL345. It's indeed a new generation chip that should be better to measure traditional "accelerations". However, a detailed comparison of their datasheets and practical measurements paints a more complex picture: the LIS2DW boasts greater sensitivity, but it has a lower sampling rate and produce significant aliasing. - -This lower sampling rate is problematic for the vibration graph because it only records data up to 200 Hz, which is too low to produce an accurate graph. This will be seen as a small low frequency band on the spectrogram with a blank area for higher frequencies and incorrect data printed in the speed profile and motor frequency profile. - -| LIS2DW vibration measurement | -| --- | -| ![](../images/vibrations_graphs/sd2w_spectrogram.png) | diff --git a/pyproject.toml b/pyproject.toml new file mode 100644 index 0000000..42306ed --- /dev/null +++ b/pyproject.toml @@ -0,0 +1,29 @@ +[project] +name = "Shake&Tune" +description = "Klipper streamlined input shaper workflow and calibration tools" +readme = "README.md" +requires-python = ">= 3.9" +authors = [ + {name = "Félix Boisselier", email = "felix@fboisselier.fr"} +] +keywords = ["klipper", "input shaper", "calibration", "3d printer"] +license = {file = "LICENSE"} + +[project.urls] +Repository = "https://github.com/Frix-x/klippain-shaketune" +Documentation = "https://github.com/Frix-x/klippain-shaketune/tree/main/docs" +Issues = "https://github.com/Frix-x/klippain-shaketune/issues" +Changelog = "https://github.com/Frix-x/klippain-shaketune/releases" + +[tool.ruff] +line-length = 120 # We all have modern screens now and I believe this should be brought in line with current technology +indent-width = 4 +target-version = "py39" + +[tool.ruff.lint] +select = ["E4", "E7", "E9", "F", "B"] +unfixable = ["B"] + +[tool.ruff.format] +quote-style = "single" +skip-magic-trailing-comma = false diff --git a/src/graph_creators/__init.py__ b/src/graph_creators/__init.py__ new file mode 100644 index 0000000..e69de29 diff --git a/K-ShakeTune/scripts/analyze_axesmap.py b/src/graph_creators/analyze_axesmap.py old mode 100755 new mode 100644 similarity index 73% rename from K-ShakeTune/scripts/analyze_axesmap.py rename to src/graph_creators/analyze_axesmap.py index 0c38641..9376cfc --- a/K-ShakeTune/scripts/analyze_axesmap.py +++ b/src/graph_creators/analyze_axesmap.py @@ -5,17 +5,12 @@ ###################################### # Written by Frix_x#0161 # -# Be sure to make this script executable using SSH: type 'chmod +x ./analyze_axesmap.py' when in the folder ! - -##################################################################### -################ !!! DO NOT EDIT BELOW THIS LINE !!! ################ -##################################################################### - import optparse + import numpy as np -from locale_utils import print_with_c_locale from scipy.signal import butter, filtfilt +from ..helpers.locale_utils import print_with_c_locale NUM_POINTS = 500 @@ -24,6 +19,7 @@ # Computation ###################################################################### + def accel_signal_filter(data, cutoff=2, fs=100, order=5): nyq = 0.5 * fs normal_cutoff = cutoff / nyq @@ -32,10 +28,12 @@ def accel_signal_filter(data, cutoff=2, fs=100, order=5): filtered_data -= np.mean(filtered_data) return filtered_data + def find_first_spike(data): min_index, max_index = np.argmin(data), np.argmax(data) return ('-', min_index) if min_index < max_index else ('', max_index) + def get_movement_vector(data, start_idx, end_idx): if start_idx < end_idx: vector = [] @@ -45,21 +43,19 @@ def get_movement_vector(data, start_idx, end_idx): else: return np.zeros(3) + def angle_between(v1, v2): v1_u = v1 / np.linalg.norm(v1) v2_u = v2 / np.linalg.norm(v2) return np.arccos(np.clip(np.dot(v1_u, v2_u), -1.0, 1.0)) + def compute_errors(filtered_data, spikes_sorted, accel_value, num_points): # Get the movement start points in the correct order from the sorted bag of spikes movement_starts = [spike[0][1] for spike in spikes_sorted] # Theoretical unit vectors for X, Y, Z printer axes - printer_axes = { - 'x': np.array([1, 0, 0]), - 'y': np.array([0, 1, 0]), - 'z': np.array([0, 0, 1]) - } + printer_axes = {'x': np.array([1, 0, 0]), 'y': np.array([0, 1, 0]), 'z': np.array([0, 0, 1])} alignment_errors = {} sensitivity_errors = {} @@ -82,6 +78,7 @@ def compute_errors(filtered_data, spikes_sorted, accel_value, num_points): # Startup and main routines ###################################################################### + def parse_log(logname): with open(logname) as f: for header in f: @@ -91,26 +88,28 @@ def parse_log(logname): # Raw accelerometer data return np.loadtxt(logname, comments='#', delimiter=',') # Power spectral density data or shaper calibration data - raise ValueError("File %s does not contain raw accelerometer data and therefore " - "is not supported by this script. Please use the official Klipper " - "calibrate_shaper.py script to process it instead." % (logname,)) + raise ValueError( + 'File %s does not contain raw accelerometer data and therefore ' + 'is not supported by this script. Please use the official Klipper ' + 'calibrate_shaper.py script to process it instead.' % (logname,) + ) def axesmap_calibration(lognames, accel=None): # Parse the raw data and get them ready for analysis raw_datas = [parse_log(filename) for filename in lognames] if len(raw_datas) > 1: - raise ValueError("Analysis of multiple CSV files at once is not possible with this script") + raise ValueError('Analysis of multiple CSV files at once is not possible with this script') - filtered_data = [accel_signal_filter(raw_datas[0][:, i+1]) for i in range(3)] + filtered_data = [accel_signal_filter(raw_datas[0][:, i + 1]) for i in range(3)] spikes = [find_first_spike(filtered_data[i]) for i in range(3)] spikes_sorted = sorted([(spikes[0], 'x'), (spikes[1], 'y'), (spikes[2], 'z')], key=lambda x: x[0][1]) # Using the previous variables to get the axes_map and errors - axes_map = ','.join([f"{spike[0][0]}{spike[1]}" for spike in spikes_sorted]) + axes_map = ','.join([f'{spike[0][0]}{spike[1]}' for spike in spikes_sorted]) # alignment_error, sensitivity_error = compute_errors(filtered_data, spikes_sorted, accel, NUM_POINTS) - results = f"Detected axes_map:\n {axes_map}\n" + results = f'Detected axes_map:\n {axes_map}\n' # TODO: work on this function that is currently not giving good results... # results += "Accelerometer angle deviation:\n" @@ -127,21 +126,21 @@ def axesmap_calibration(lognames, accel=None): def main(): # Parse command-line arguments - usage = "%prog [options] " + usage = '%prog [options] ' opts = optparse.OptionParser(usage) - opts.add_option("-o", "--output", type="string", dest="output", - default=None, help="filename of output graph") - opts.add_option("-a", "--accel", type="string", dest="accel", - default=None, help="acceleration value used to do the movements") + opts.add_option('-o', '--output', type='string', dest='output', default=None, help='filename of output graph') + opts.add_option( + '-a', '--accel', type='string', dest='accel', default=None, help='acceleration value used to do the movements' + ) options, args = opts.parse_args() if len(args) < 1: - opts.error("No CSV file(s) to analyse") + opts.error('No CSV file(s) to analyse') if options.accel is None: - opts.error("You must specify the acceleration value used when generating the CSV file (option -a)") + opts.error('You must specify the acceleration value used when generating the CSV file (option -a)') try: accel_value = float(options.accel) except ValueError: - opts.error("Invalid acceleration value. It should be a numeric value.") + opts.error('Invalid acceleration value. It should be a numeric value.') results = axesmap_calibration(args, accel_value) print_with_c_locale(results) diff --git a/K-ShakeTune/scripts/graph_belts.py b/src/graph_creators/graph_belts.py old mode 100755 new mode 100644 similarity index 63% rename from K-ShakeTune/scripts/graph_belts.py rename to src/graph_creators/graph_belts.py index d30158e..edbb316 --- a/K-ShakeTune/scripts/graph_belts.py +++ b/src/graph_creators/graph_belts.py @@ -5,27 +5,31 @@ ################################################# # Written by Frix_x#0161 # -# Be sure to make this script executable using SSH: type 'chmod +x ./graph_belts.py' when in the folder! - -##################################################################### -################ !!! DO NOT EDIT BELOW THIS LINE !!! ################ -##################################################################### - -import optparse, matplotlib, sys, importlib, os -from datetime import datetime +import optparse +import os from collections import namedtuple +from datetime import datetime + +import matplotlib +import matplotlib.colors +import matplotlib.font_manager +import matplotlib.pyplot as plt +import matplotlib.ticker import numpy as np -import matplotlib.pyplot as plt -import matplotlib.font_manager, matplotlib.ticker, matplotlib.colors from scipy.interpolate import griddata matplotlib.use('Agg') -from locale_utils import set_locale, print_with_c_locale -from common_func import compute_spectrogram, detect_peaks, get_git_version, parse_log, setup_klipper_import - +from ..helpers.common_func import ( + compute_curve_similarity_factor, + compute_spectrogram, + detect_peaks, + parse_log, + setup_klipper_import, +) +from ..helpers.locale_utils import print_with_c_locale, set_locale -ALPHABET = "ABCDEFGHIJKLMNOPQRSTUVWXYZ" # For paired peaks names +ALPHABET = 'ABCDEFGHIJKLMNOPQRSTUVWXYZ' # For paired peaks names PEAKS_DETECTION_THRESHOLD = 0.20 CURVE_SIMILARITY_SIGMOID_K = 0.6 @@ -37,11 +41,11 @@ SignalData = namedtuple('CalibrationData', ['freqs', 'psd', 'peaks', 'paired_peaks', 'unpaired_peaks']) KLIPPAIN_COLORS = { - "purple": "#70088C", - "orange": "#FF8D32", - "dark_purple": "#150140", - "dark_orange": "#F24130", - "red_pink": "#F2055C" + 'purple': '#70088C', + 'orange': '#FF8D32', + 'dark_purple': '#150140', + 'dark_orange': '#F24130', + 'red_pink': '#F2055C', } @@ -49,27 +53,6 @@ # Computation of the PSD graph ###################################################################### -# Calculate or estimate a "similarity" factor between two PSD curves and scale it to a percentage. This is -# used here to quantify how close the two belts path behavior and responses are close together. -def compute_curve_similarity_factor(signal1, signal2): - freqs1 = signal1.freqs - psd1 = signal1.psd - freqs2 = signal2.freqs - psd2 = signal2.psd - - # Interpolate PSDs to match the same frequency bins and do a cross-correlation - psd2_interp = np.interp(freqs1, freqs2, psd2) - cross_corr = np.correlate(psd1, psd2_interp, mode='full') - - # Find the peak of the cross-correlation and compute a similarity normalized by the energy of the signals - peak_value = np.max(cross_corr) - similarity = peak_value / (np.sqrt(np.sum(psd1**2) * np.sum(psd2_interp**2))) - - # Apply sigmoid scaling to get better numbers and get a final percentage value - scaled_similarity = sigmoid_scale(-np.log(1 - similarity), CURVE_SIMILARITY_SIGMOID_K) - - return scaled_similarity - # This function create pairs of peaks that are close in frequency on two curves (that are known # to be resonances points and must be similar on both belts on a CoreXY kinematic) @@ -81,37 +64,37 @@ def pair_peaks(peaks1, freqs1, psd1, peaks2, freqs2, psd2): for p2 in peaks2: distances.append(abs(freqs1[p1] - freqs2[p2])) distances = np.array(distances) - + median_distance = np.median(distances) iqr = np.percentile(distances, 75) - np.percentile(distances, 25) - + threshold = median_distance + 1.5 * iqr threshold = min(threshold, 10) - + # Pair the peaks using the dynamic thresold paired_peaks = [] unpaired_peaks1 = list(peaks1) unpaired_peaks2 = list(peaks2) - + while unpaired_peaks1 and unpaired_peaks2: min_distance = threshold + 1 pair = None - + for p1 in unpaired_peaks1: for p2 in unpaired_peaks2: distance = abs(freqs1[p1] - freqs2[p2]) if distance < min_distance: min_distance = distance pair = (p1, p2) - - if pair is None: # No more pairs below the threshold + + if pair is None: # No more pairs below the threshold break - + p1, p2 = pair paired_peaks.append(((p1, freqs1[p1], psd1[p1]), (p2, freqs2[p2], psd2[p2]))) unpaired_peaks1.remove(p1) unpaired_peaks2.remove(p2) - + return paired_peaks, unpaired_peaks1, unpaired_peaks2 @@ -119,6 +102,7 @@ def pair_peaks(peaks1, freqs1, psd1, peaks2, freqs2, psd2): # Computation of the differential spectrogram ###################################################################### + # Interpolate source_data (2D) to match target_x and target_y in order to # get similar time and frequency dimensions for the differential spectrogram def interpolate_2d(target_x, target_y, source_x, source_y, source_data): @@ -146,7 +130,7 @@ def compute_combined_spectrogram(data1, data2): # Interpolate the spectrograms pdata2_interpolated = interpolate_2d(bins1, t1, bins2, t2, pdata2) - + # Combine them in two form: a summed diff for the MHI computation and a diverging diff for the spectrogram colors combined_sum = np.abs(pdata1 - pdata2_interpolated) combined_divergent = pdata1 - pdata2_interpolated @@ -168,58 +152,61 @@ def compute_mhi(combined_data, similarity_coefficient, num_unpaired_peaks): total_variability_metric = np.sum(np.abs(np.gradient(filtered_data))) # Scale the metric to a percentage using the threshold (found empirically on a large number of user data shared to me) base_percentage = (np.log1p(total_variability_metric) / np.log1p(DC_THRESHOLD_METRIC)) * 100 - + # Adjust the percentage based on the similarity_coefficient to add a grain of salt adjusted_percentage = base_percentage * (1 - DC_GRAIN_OF_SALT_FACTOR * (similarity_coefficient / 100)) # Adjust the percentage again based on the number of unpaired peaks to add a second grain of salt peak_confidence = num_unpaired_peaks / DC_MAX_UNPAIRED_PEAKS_ALLOWED final_percentage = (1 - peak_confidence) * adjusted_percentage + peak_confidence * 100 - + # Ensure the result lies between 0 and 100 by clipping the computed value final_percentage = np.clip(final_percentage, 0, 100) - + return final_percentage, mhi_lut(final_percentage) # LUT to transform the MHI into a textual value easy to understand for the users of the script -def mhi_lut(mhi): +def mhi_lut(mhi): ranges = [ - (0, 30, "Excellent mechanical health"), - (30, 45, "Good mechanical health"), - (45, 55, "Acceptable mechanical health"), - (55, 70, "Potential signs of a mechanical issue"), - (70, 85, "Likely a mechanical issue"), - (85, 100, "Mechanical issue detected") + (0, 30, 'Excellent mechanical health'), + (30, 45, 'Good mechanical health'), + (45, 55, 'Acceptable mechanical health'), + (55, 70, 'Potential signs of a mechanical issue'), + (70, 85, 'Likely a mechanical issue'), + (85, 100, 'Mechanical issue detected'), ] for lower, upper, message in ranges: if lower < mhi <= upper: return message - return "Error computing MHI value" + return 'Error computing MHI value' ###################################################################### # Graphing ###################################################################### + def plot_compare_frequency(ax, lognames, signal1, signal2, similarity_factor, max_freq): # Get the belt name for the legend to avoid putting the full file name signal1_belt = (lognames[0].split('/')[-1]).split('_')[-1][0] signal2_belt = (lognames[1].split('/')[-1]).split('_')[-1][0] if signal1_belt == 'A' and signal2_belt == 'B': - signal1_belt += " (axis 1,-1)" - signal2_belt += " (axis 1, 1)" + signal1_belt += ' (axis 1,-1)' + signal2_belt += ' (axis 1, 1)' elif signal1_belt == 'B' and signal2_belt == 'A': - signal1_belt += " (axis 1, 1)" - signal2_belt += " (axis 1,-1)" + signal1_belt += ' (axis 1, 1)' + signal2_belt += ' (axis 1,-1)' else: - print_with_c_locale("Warning: belts doesn't seem to have the correct name A and B (extracted from the filename.csv)") + print_with_c_locale( + "Warning: belts doesn't seem to have the correct name A and B (extracted from the filename.csv)" + ) # Plot the two belts PSD signals - ax.plot(signal1.freqs, signal1.psd, label="Belt " + signal1_belt, color=KLIPPAIN_COLORS['purple']) - ax.plot(signal2.freqs, signal2.psd, label="Belt " + signal2_belt, color=KLIPPAIN_COLORS['orange']) + ax.plot(signal1.freqs, signal1.psd, label='Belt ' + signal1_belt, color=KLIPPAIN_COLORS['purple']) + ax.plot(signal2.freqs, signal2.psd, label='Belt ' + signal2_belt, color=KLIPPAIN_COLORS['orange']) # Trace the "relax region" (also used as a threshold to filter and detect the peaks) psd_lowest_max = min(signal1.psd.max(), signal2.psd.max()) @@ -234,38 +221,71 @@ def plot_compare_frequency(ax, lognames, signal1, signal2, similarity_factor, ma for _, (peak1, peak2) in enumerate(signal1.paired_peaks): label = ALPHABET[paired_peak_count] - amplitude_offset = abs(((signal2.psd[peak2[0]] - signal1.psd[peak1[0]]) / max(signal1.psd[peak1[0]], signal2.psd[peak2[0]])) * 100) + amplitude_offset = abs( + ((signal2.psd[peak2[0]] - signal1.psd[peak1[0]]) / max(signal1.psd[peak1[0]], signal2.psd[peak2[0]])) * 100 + ) frequency_offset = abs(signal2.freqs[peak2[0]] - signal1.freqs[peak1[0]]) - offsets_table_data.append([f"Peaks {label}", f"{frequency_offset:.1f} Hz", f"{amplitude_offset:.1f} %"]) - - ax.plot(signal1.freqs[peak1[0]], signal1.psd[peak1[0]], "x", color='black') - ax.plot(signal2.freqs[peak2[0]], signal2.psd[peak2[0]], "x", color='black') - ax.plot([signal1.freqs[peak1[0]], signal2.freqs[peak2[0]]], [signal1.psd[peak1[0]], signal2.psd[peak2[0]]], ":", color='gray') - - ax.annotate(label + "1", (signal1.freqs[peak1[0]], signal1.psd[peak1[0]]), - textcoords="offset points", xytext=(8, 5), - ha='left', fontsize=13, color='black') - ax.annotate(label + "2", (signal2.freqs[peak2[0]], signal2.psd[peak2[0]]), - textcoords="offset points", xytext=(8, 5), - ha='left', fontsize=13, color='black') + offsets_table_data.append([f'Peaks {label}', f'{frequency_offset:.1f} Hz', f'{amplitude_offset:.1f} %']) + + ax.plot(signal1.freqs[peak1[0]], signal1.psd[peak1[0]], 'x', color='black') + ax.plot(signal2.freqs[peak2[0]], signal2.psd[peak2[0]], 'x', color='black') + ax.plot( + [signal1.freqs[peak1[0]], signal2.freqs[peak2[0]]], + [signal1.psd[peak1[0]], signal2.psd[peak2[0]]], + ':', + color='gray', + ) + + ax.annotate( + label + '1', + (signal1.freqs[peak1[0]], signal1.psd[peak1[0]]), + textcoords='offset points', + xytext=(8, 5), + ha='left', + fontsize=13, + color='black', + ) + ax.annotate( + label + '2', + (signal2.freqs[peak2[0]], signal2.psd[peak2[0]]), + textcoords='offset points', + xytext=(8, 5), + ha='left', + fontsize=13, + color='black', + ) paired_peak_count += 1 for peak in signal1.unpaired_peaks: - ax.plot(signal1.freqs[peak], signal1.psd[peak], "x", color='black') - ax.annotate(str(unpaired_peak_count + 1), (signal1.freqs[peak], signal1.psd[peak]), - textcoords="offset points", xytext=(8, 5), - ha='left', fontsize=13, color='red', weight='bold') + ax.plot(signal1.freqs[peak], signal1.psd[peak], 'x', color='black') + ax.annotate( + str(unpaired_peak_count + 1), + (signal1.freqs[peak], signal1.psd[peak]), + textcoords='offset points', + xytext=(8, 5), + ha='left', + fontsize=13, + color='red', + weight='bold', + ) unpaired_peak_count += 1 for peak in signal2.unpaired_peaks: - ax.plot(signal2.freqs[peak], signal2.psd[peak], "x", color='black') - ax.annotate(str(unpaired_peak_count + 1), (signal2.freqs[peak], signal2.psd[peak]), - textcoords="offset points", xytext=(8, 5), - ha='left', fontsize=13, color='red', weight='bold') + ax.plot(signal2.freqs[peak], signal2.psd[peak], 'x', color='black') + ax.annotate( + str(unpaired_peak_count + 1), + (signal2.freqs[peak], signal2.psd[peak]), + textcoords='offset points', + xytext=(8, 5), + ha='left', + fontsize=13, + color='red', + weight='bold', + ) unpaired_peak_count += 1 # Add estimated similarity to the graph - ax2 = ax.twinx() # To split the legends in two box + ax2 = ax.twinx() # To split the legends in two box ax2.yaxis.set_visible(False) ax2.plot([], [], ' ', label=f'Estimated similarity: {similarity_factor:.1f}%') ax2.plot([], [], ' ', label=f'Number of unpaired peaks: {unpaired_peak_count}') @@ -279,17 +299,32 @@ def plot_compare_frequency(ax, lognames, signal1, signal2, similarity_factor, ma ax.xaxis.set_minor_locator(matplotlib.ticker.AutoMinorLocator()) ax.yaxis.set_minor_locator(matplotlib.ticker.AutoMinorLocator()) - ax.ticklabel_format(axis='y', style='scientific', scilimits=(0,0)) + ax.ticklabel_format(axis='y', style='scientific', scilimits=(0, 0)) ax.grid(which='major', color='grey') ax.grid(which='minor', color='lightgrey') fontP = matplotlib.font_manager.FontProperties() fontP.set_size('small') - ax.set_title('Belts Frequency Profiles (estimated similarity: {:.1f}%)'.format(similarity_factor), fontsize=14, color=KLIPPAIN_COLORS['dark_orange'], weight='bold') + ax.set_title( + 'Belts Frequency Profiles (estimated similarity: {:.1f}%)'.format(similarity_factor), + fontsize=14, + color=KLIPPAIN_COLORS['dark_orange'], + weight='bold', + ) # Print the table of offsets ontop of the graph below the original legend (upper right) if len(offsets_table_data) > 0: - columns = ["", "Frequency delta", "Amplitude delta", ] - offset_table = ax.table(cellText=offsets_table_data, colLabels=columns, bbox=[0.66, 0.75, 0.33, 0.15], loc='upper right', cellLoc='center') + columns = [ + '', + 'Frequency delta', + 'Amplitude delta', + ] + offset_table = ax.table( + cellText=offsets_table_data, + colLabels=columns, + bbox=[0.66, 0.75, 0.33, 0.15], + loc='upper right', + cellLoc='center', + ) offset_table.auto_set_font_size(False) offset_table.set_fontsize(8) offset_table.auto_set_column_width([0, 1, 2]) @@ -306,19 +341,35 @@ def plot_compare_frequency(ax, lognames, signal1, signal2, similarity_factor, ma def plot_difference_spectrogram(ax, signal1, signal2, t, bins, combined_divergent, textual_mhi, max_freq): - ax.set_title(f"Differential Spectrogram", fontsize=14, color=KLIPPAIN_COLORS['dark_orange'], weight='bold') + ax.set_title('Differential Spectrogram', fontsize=14, color=KLIPPAIN_COLORS['dark_orange'], weight='bold') ax.plot([], [], ' ', label=f'{textual_mhi} (experimental)') - + # Draw the differential spectrogram with a specific custom norm to get orange or purple values where there is signal or white near zeros # imgshow is better suited here than pcolormesh since its result is already rasterized and we doesn't need to keep vector graphics # when saving to a final .png file. Using it also allow to save ~150-200MB of RAM during the "fig.savefig" operation. - colors = [KLIPPAIN_COLORS['dark_orange'], KLIPPAIN_COLORS['orange'], 'white', KLIPPAIN_COLORS['purple'], KLIPPAIN_COLORS['dark_purple']] - cm = matplotlib.colors.LinearSegmentedColormap.from_list('klippain_divergent', list(zip([0, 0.25, 0.5, 0.75, 1], colors))) + colors = [ + KLIPPAIN_COLORS['dark_orange'], + KLIPPAIN_COLORS['orange'], + 'white', + KLIPPAIN_COLORS['purple'], + KLIPPAIN_COLORS['dark_purple'], + ] + cm = matplotlib.colors.LinearSegmentedColormap.from_list( + 'klippain_divergent', list(zip([0, 0.25, 0.5, 0.75, 1], colors)) + ) norm = matplotlib.colors.TwoSlopeNorm(vmin=np.min(combined_divergent), vcenter=0, vmax=np.max(combined_divergent)) - ax.imshow(combined_divergent.T, cmap=cm, norm=norm, aspect='auto', extent=[t[0], t[-1], bins[0], bins[-1]], interpolation='bilinear', origin='lower') + ax.imshow( + combined_divergent.T, + cmap=cm, + norm=norm, + aspect='auto', + extent=[t[0], t[-1], bins[0], bins[-1]], + interpolation='bilinear', + origin='lower', + ) ax.set_xlabel('Frequency (hz)') - ax.set_xlim([0., max_freq]) + ax.set_xlim([0.0, max_freq]) ax.set_ylabel('Time (s)') ax.set_ylim([0, bins[-1]]) @@ -330,18 +381,32 @@ def plot_difference_spectrogram(ax, signal1, signal2, t, bins, combined_divergen unpaired_peak_count = 0 for _, peak in enumerate(signal1.unpaired_peaks): ax.axvline(signal1.freqs[peak], color=KLIPPAIN_COLORS['red_pink'], linestyle='dotted', linewidth=1.5) - ax.annotate(f"Peak {unpaired_peak_count + 1}", (signal1.freqs[peak], t[-1]*0.05), - textcoords="data", color=KLIPPAIN_COLORS['red_pink'], rotation=90, fontsize=10, - verticalalignment='bottom', horizontalalignment='right') - unpaired_peak_count +=1 + ax.annotate( + f'Peak {unpaired_peak_count + 1}', + (signal1.freqs[peak], t[-1] * 0.05), + textcoords='data', + color=KLIPPAIN_COLORS['red_pink'], + rotation=90, + fontsize=10, + verticalalignment='bottom', + horizontalalignment='right', + ) + unpaired_peak_count += 1 for _, peak in enumerate(signal2.unpaired_peaks): ax.axvline(signal2.freqs[peak], color=KLIPPAIN_COLORS['red_pink'], linestyle='dotted', linewidth=1.5) - ax.annotate(f"Peak {unpaired_peak_count + 1}", (signal2.freqs[peak], t[-1]*0.05), - textcoords="data", color=KLIPPAIN_COLORS['red_pink'], rotation=90, fontsize=10, - verticalalignment='bottom', horizontalalignment='right') - unpaired_peak_count +=1 - + ax.annotate( + f'Peak {unpaired_peak_count + 1}', + (signal2.freqs[peak], t[-1] * 0.05), + textcoords='data', + color=KLIPPAIN_COLORS['red_pink'], + rotation=90, + fontsize=10, + verticalalignment='bottom', + horizontalalignment='right', + ) + unpaired_peak_count += 1 + # Plot vertical lines and zones for paired peaks for idx, (peak1, peak2) in enumerate(signal1.paired_peaks): label = ALPHABET[idx] @@ -350,20 +415,24 @@ def plot_difference_spectrogram(ax, signal1, signal2, t, bins, combined_divergen ax.axvline(x_min, color=KLIPPAIN_COLORS['dark_purple'], linestyle='dotted', linewidth=1.5) ax.axvline(x_max, color=KLIPPAIN_COLORS['dark_purple'], linestyle='dotted', linewidth=1.5) ax.fill_between([x_min, x_max], 0, np.max(combined_divergent), color=KLIPPAIN_COLORS['dark_purple'], alpha=0.3) - ax.annotate(f"Peaks {label}", (x_min, t[-1]*0.05), - textcoords="data", color=KLIPPAIN_COLORS['dark_purple'], rotation=90, fontsize=10, - verticalalignment='bottom', horizontalalignment='right') + ax.annotate( + f'Peaks {label}', + (x_min, t[-1] * 0.05), + textcoords='data', + color=KLIPPAIN_COLORS['dark_purple'], + rotation=90, + fontsize=10, + verticalalignment='bottom', + horizontalalignment='right', + ) return ###################################################################### -# Custom tools +# Custom tools ###################################################################### -# Simple helper to compute a sigmoid scalling (from 0 to 100%) -def sigmoid_scale(x, k=1): - return 1 / (1 + np.exp(-k * x)) * 100 # Original Klipper function to get the PSD data of a raw accelerometer signal def compute_signal_data(data, max_freq): @@ -376,13 +445,14 @@ def compute_signal_data(data, max_freq): _, peaks, _ = detect_peaks(psd, freqs, PEAKS_DETECTION_THRESHOLD * psd.max()) return SignalData(freqs=freqs, psd=psd, peaks=peaks, paired_peaks=None, unpaired_peaks=None) - - + + ###################################################################### # Startup and main routines ###################################################################### -def belts_calibration(lognames, klipperdir="~/klipper", max_freq=200.): + +def belts_calibration(lognames, klipperdir='~/klipper', max_freq=200.0, st_version=None): set_locale() global shaper_calibrate shaper_calibrate = setup_klipper_import(klipperdir) @@ -390,7 +460,7 @@ def belts_calibration(lognames, klipperdir="~/klipper", max_freq=200.): # Parse data datas = [parse_log(fn) for fn in lognames] if len(datas) > 2: - raise ValueError("Incorrect number of .csv files used (this function needs exactly two files to compare them)!") + raise ValueError('Incorrect number of .csv files used (this function needs exactly two files to compare them)!') # Compute calibration data for the two datasets with automatic peaks detection signal1 = compute_signal_data(datas[0], max_freq) @@ -399,55 +469,67 @@ def belts_calibration(lognames, klipperdir="~/klipper", max_freq=200.): del datas # Pair the peaks across the two datasets - paired_peaks, unpaired_peaks1, unpaired_peaks2 = pair_peaks(signal1.peaks, signal1.freqs, signal1.psd, - signal2.peaks, signal2.freqs, signal2.psd) - signal1 = signal1._replace(paired_peaks = paired_peaks, unpaired_peaks = unpaired_peaks1) - signal2 = signal2._replace(paired_peaks = paired_peaks, unpaired_peaks = unpaired_peaks2) + paired_peaks, unpaired_peaks1, unpaired_peaks2 = pair_peaks( + signal1.peaks, signal1.freqs, signal1.psd, signal2.peaks, signal2.freqs, signal2.psd + ) + signal1 = signal1._replace(paired_peaks=paired_peaks, unpaired_peaks=unpaired_peaks1) + signal2 = signal2._replace(paired_peaks=paired_peaks, unpaired_peaks=unpaired_peaks2) # Compute the similarity (using cross-correlation of the PSD signals) - similarity_factor = compute_curve_similarity_factor(signal1, signal2) - print_with_c_locale(f"Belts estimated similarity: {similarity_factor:.1f}%") + similarity_factor = compute_curve_similarity_factor( + signal1.freqs, signal1.psd, signal2.freqs, signal2.psd, CURVE_SIMILARITY_SIGMOID_K + ) + print_with_c_locale(f'Belts estimated similarity: {similarity_factor:.1f}%') # Compute the MHI value from the differential spectrogram sum of gradient, salted with the similarity factor and the number of # unpaired peaks from the belts frequency profile. Be careful, this value is highly opinionated and is pretty experimental! - mhi, textual_mhi = compute_mhi(combined_sum, similarity_factor, len(signal1.unpaired_peaks) + len(signal2.unpaired_peaks)) - print_with_c_locale(f"[experimental] Mechanical Health Indicator: {textual_mhi.lower()} ({mhi:.1f}%)") + mhi, textual_mhi = compute_mhi( + combined_sum, similarity_factor, len(signal1.unpaired_peaks) + len(signal2.unpaired_peaks) + ) + print_with_c_locale(f'[experimental] Mechanical Health Indicator: {textual_mhi.lower()} ({mhi:.1f}%)') # Create graph layout - fig, (ax1, ax2) = plt.subplots(2, 1, gridspec_kw={ - 'height_ratios':[4, 3], - 'bottom':0.050, - 'top':0.890, - 'left':0.085, - 'right':0.966, - 'hspace':0.169, - 'wspace':0.200 - }) + fig, (ax1, ax2) = plt.subplots( + 2, + 1, + gridspec_kw={ + 'height_ratios': [4, 3], + 'bottom': 0.050, + 'top': 0.890, + 'left': 0.085, + 'right': 0.966, + 'hspace': 0.169, + 'wspace': 0.200, + }, + ) fig.set_size_inches(8.3, 11.6) # Add title - title_line1 = "RELATIVE BELT CALIBRATION TOOL" - fig.text(0.12, 0.965, title_line1, ha='left', va='bottom', fontsize=20, color=KLIPPAIN_COLORS['purple'], weight='bold') + title_line1 = 'RELATIVE BELTS CALIBRATION TOOL' + fig.text( + 0.12, 0.965, title_line1, ha='left', va='bottom', fontsize=20, color=KLIPPAIN_COLORS['purple'], weight='bold' + ) try: filename = lognames[0].split('/')[-1] - dt = datetime.strptime(f"{filename.split('_')[1]} {filename.split('_')[2]}", "%Y%m%d %H%M%S") + dt = datetime.strptime(f"{filename.split('_')[1]} {filename.split('_')[2]}", '%Y%m%d %H%M%S') title_line2 = dt.strftime('%x %X') - except: - print_with_c_locale("Warning: CSV filenames look to be different than expected (%s , %s)" % (lognames[0], lognames[1])) - title_line2 = lognames[0].split('/')[-1] + " / " + lognames[1].split('/')[-1] + except Exception: + print_with_c_locale( + 'Warning: CSV filenames look to be different than expected (%s , %s)' % (lognames[0], lognames[1]) + ) + title_line2 = lognames[0].split('/')[-1] + ' / ' + lognames[1].split('/')[-1] fig.text(0.12, 0.957, title_line2, ha='left', va='top', fontsize=16, color=KLIPPAIN_COLORS['dark_purple']) # Plot the graphs plot_compare_frequency(ax1, lognames, signal1, signal2, similarity_factor, max_freq) plot_difference_spectrogram(ax2, signal1, signal2, t, bins, combined_divergent, textual_mhi, max_freq) - + # Adding a small Klippain logo to the top left corner of the figure ax_logo = fig.add_axes([0.001, 0.8995, 0.1, 0.1], anchor='NW') ax_logo.imshow(plt.imread(os.path.join(os.path.dirname(os.path.abspath(__file__)), 'klippain.png'))) ax_logo.axis('off') # Adding Shake&Tune version in the top right corner - st_version = get_git_version() - if st_version is not None: + if st_version != 'unknown': fig.text(0.995, 0.985, st_version, ha='right', va='bottom', fontsize=8, color=KLIPPAIN_COLORS['purple']) return fig @@ -455,19 +537,18 @@ def belts_calibration(lognames, klipperdir="~/klipper", max_freq=200.): def main(): # Parse command-line arguments - usage = "%prog [options] " + usage = '%prog [options] ' opts = optparse.OptionParser(usage) - opts.add_option("-o", "--output", type="string", dest="output", - default=None, help="filename of output graph") - opts.add_option("-f", "--max_freq", type="float", default=200., - help="maximum frequency to graph") - opts.add_option("-k", "--klipper_dir", type="string", dest="klipperdir", - default="~/klipper", help="main klipper directory") + opts.add_option('-o', '--output', type='string', dest='output', default=None, help='filename of output graph') + opts.add_option('-f', '--max_freq', type='float', default=200.0, help='maximum frequency to graph') + opts.add_option( + '-k', '--klipper_dir', type='string', dest='klipperdir', default='~/klipper', help='main klipper directory' + ) options, args = opts.parse_args() if len(args) < 1: - opts.error("Incorrect number of arguments") + opts.error('Incorrect number of arguments') if options.output is None: - opts.error("You must specify an output file.png to use the script (option -o)") + opts.error('You must specify an output file.png to use the script (option -o)') fig = belts_calibration(args, options.klipperdir, options.max_freq) fig.savefig(options.output, dpi=150) diff --git a/K-ShakeTune/scripts/graph_shaper.py b/src/graph_creators/graph_shaper.py old mode 100755 new mode 100644 similarity index 54% rename from K-ShakeTune/scripts/graph_shaper.py rename to src/graph_creators/graph_shaper.py index 274ec34..aec74db --- a/K-ShakeTune/scripts/graph_shaper.py +++ b/src/graph_creators/graph_shaper.py @@ -6,25 +6,28 @@ # Derived from the calibrate_shaper.py official Klipper script # Copyright (C) 2020 Dmitry Butyugin # Copyright (C) 2020 Kevin O'Connor -# Written by Frix_x#0161 # +# Highly modified and improved by Frix_x#0161 # -# Be sure to make this script executable using SSH: type 'chmod +x ./graph_shaper.py' when in the folder! - -##################################################################### -################ !!! DO NOT EDIT BELOW THIS LINE !!! ################ -##################################################################### - -import optparse, matplotlib, os +import optparse +import os from datetime import datetime -import numpy as np + +import matplotlib +import matplotlib.font_manager import matplotlib.pyplot as plt -import matplotlib.font_manager, matplotlib.ticker +import matplotlib.ticker +import numpy as np matplotlib.use('Agg') -from locale_utils import set_locale, print_with_c_locale -from common_func import compute_mechanical_parameters, compute_spectrogram, detect_peaks, get_git_version, parse_log, setup_klipper_import - +from ..helpers.common_func import ( + compute_mechanical_parameters, + compute_spectrogram, + detect_peaks, + parse_log, + setup_klipper_import, +) +from ..helpers.locale_utils import print_with_c_locale, set_locale PEAKS_DETECTION_THRESHOLD = 0.05 PEAKS_EFFECT_THRESHOLD = 0.12 @@ -32,9 +35,11 @@ MAX_SMOOTHING = 0.1 KLIPPAIN_COLORS = { - "purple": "#70088C", - "dark_purple": "#150140", - "dark_orange": "#F24130" + 'purple': '#70088C', + 'orange': '#FF8D32', + 'dark_purple': '#150140', + 'dark_orange': '#F24130', + 'red_pink': '#F2055C', } @@ -42,6 +47,7 @@ # Computation ###################################################################### + # Find the best shaper parameters using Klipper's official algorithm selection with # a proper precomputed damping ratio (zeta) and using the configured printer SQV value def calibrate_shaper(datas, max_smoothing, scv, max_freq): @@ -49,30 +55,57 @@ def calibrate_shaper(datas, max_smoothing, scv, max_freq): calibration_data = helper.process_accelerometer_data(datas) calibration_data.normalize_to_frequencies() - fr, zeta, _ = compute_mechanical_parameters(calibration_data.psd_sum, calibration_data.freq_bins) - - shaper, all_shapers = helper.find_best_shaper( - calibration_data, shapers=None, damping_ratio=zeta, - scv=scv, shaper_freqs=None, max_smoothing=max_smoothing, - test_damping_ratios=None, max_freq=max_freq, - logger=print_with_c_locale) + fr, zeta, _, _ = compute_mechanical_parameters(calibration_data.psd_sum, calibration_data.freq_bins) - print_with_c_locale("\n-> Recommended shaper is %s @ %.1f Hz (when using a square corner velocity of %.1f and a computed damping ratio of %.3f)" % (shaper.name.upper(), shaper.freq, scv, zeta)) + # If the damping ratio computation fail, we use Klipper default value instead + if zeta is None: + zeta = 0.1 - return shaper.name, all_shapers, calibration_data, fr, zeta + compat = False + try: + shaper, all_shapers = helper.find_best_shaper( + calibration_data, + shapers=None, + damping_ratio=zeta, + scv=scv, + shaper_freqs=None, + max_smoothing=max_smoothing, + test_damping_ratios=None, + max_freq=max_freq, + logger=print_with_c_locale, + ) + except TypeError: + print_with_c_locale( + '[WARNING] You seem to be using an older version of Klipper that is not compatible with all the latest Shake&Tune features!' + ) + print_with_c_locale( + 'Shake&Tune now runs in compatibility mode: be aware that the results may be slightly off, since the real damping ratio cannot be used to create the filter recommendations' + ) + compat = True + shaper, all_shapers = helper.find_best_shaper(calibration_data, max_smoothing, print_with_c_locale) + + print_with_c_locale( + '\n-> Recommended shaper is %s @ %.1f Hz (when using a square corner velocity of %.1f and a damping ratio of %.3f)' + % (shaper.name.upper(), shaper.freq, scv, zeta) + ) + + return shaper.name, all_shapers, calibration_data, fr, zeta, compat ###################################################################### # Graphing ###################################################################### -def plot_freq_response(ax, calibration_data, shapers, performance_shaper, peaks, peaks_freqs, peaks_threshold, fr, zeta, max_freq): + +def plot_freq_response( + ax, calibration_data, shapers, performance_shaper, peaks, peaks_freqs, peaks_threshold, fr, zeta, max_freq +): freqs = calibration_data.freqs psd = calibration_data.psd_sum px = calibration_data.psd_x py = calibration_data.psd_y pz = calibration_data.psd_z - + fontP = matplotlib.font_manager.FontProperties() fontP.set_size('x-small') @@ -88,36 +121,42 @@ def plot_freq_response(ax, calibration_data, shapers, performance_shaper, peaks, ax.xaxis.set_minor_locator(matplotlib.ticker.MultipleLocator(5)) ax.yaxis.set_minor_locator(matplotlib.ticker.AutoMinorLocator()) - ax.ticklabel_format(axis='y', style='scientific', scilimits=(0,0)) + ax.ticklabel_format(axis='y', style='scientific', scilimits=(0, 0)) ax.grid(which='major', color='grey') ax.grid(which='minor', color='lightgrey') ax2 = ax.twinx() ax2.yaxis.set_visible(False) - + lowvib_shaper_vibrs = float('inf') lowvib_shaper = None lowvib_shaper_freq = None lowvib_shaper_accel = 0 - + # Draw the shappers curves and add their specific parameters in the legend # This adds also a way to find the best shaper with a low level of vibrations (with a resonable level of smoothing) for shaper in shapers: - shaper_max_accel = round(shaper.max_accel / 100.) * 100. - label = "%s (%.1f Hz, vibr=%.1f%%, sm~=%.2f, accel<=%.f)" % ( - shaper.name.upper(), shaper.freq, - shaper.vibrs * 100., shaper.smoothing, - shaper_max_accel) + shaper_max_accel = round(shaper.max_accel / 100.0) * 100.0 + label = '%s (%.1f Hz, vibr=%.1f%%, sm~=%.2f, accel<=%.f)' % ( + shaper.name.upper(), + shaper.freq, + shaper.vibrs * 100.0, + shaper.smoothing, + shaper_max_accel, + ) ax2.plot(freqs, shaper.vals, label=label, linestyle='dotted') # Get the performance shaper if shaper.name == performance_shaper: performance_shaper_freq = shaper.freq - performance_shaper_vibr = shaper.vibrs * 100. + performance_shaper_vibr = shaper.vibrs * 100.0 performance_shaper_vals = shaper.vals # Get the low vibration shaper - if (shaper.vibrs * 100 < lowvib_shaper_vibrs or (shaper.vibrs * 100 == lowvib_shaper_vibrs and shaper_max_accel > lowvib_shaper_accel)) and shaper.smoothing < MAX_SMOOTHING: + if ( + shaper.vibrs * 100 < lowvib_shaper_vibrs + or (shaper.vibrs * 100 == lowvib_shaper_vibrs and shaper_max_accel > lowvib_shaper_accel) + ) and shaper.smoothing < MAX_SMOOTHING: lowvib_shaper_accel = shaper_max_accel lowvib_shaper = shaper.name lowvib_shaper_freq = shaper.freq @@ -128,21 +167,45 @@ def plot_freq_response(ax, calibration_data, shapers, performance_shaper, peaks, # and the other one is the custom "low vibration" recommendation that looks for a suitable shaper that doesn't have excessive # smoothing and that have a lower vibration level. If both recommendation are the same shaper, or if no suitable "low # vibration" shaper is found, then only a single line as the "best shaper" recommendation is added to the legend - if lowvib_shaper != None and lowvib_shaper != performance_shaper and lowvib_shaper_vibrs <= performance_shaper_vibr: - ax2.plot([], [], ' ', label="Recommended performance shaper: %s @ %.1f Hz" % (performance_shaper.upper(), performance_shaper_freq)) - ax.plot(freqs, psd * performance_shaper_vals, label='With %s applied' % (performance_shaper.upper()), color='cyan') - ax2.plot([], [], ' ', label="Recommended low vibrations shaper: %s @ %.1f Hz" % (lowvib_shaper.upper(), lowvib_shaper_freq)) + if ( + lowvib_shaper is not None + and lowvib_shaper != performance_shaper + and lowvib_shaper_vibrs <= performance_shaper_vibr + ): + ax2.plot( + [], + [], + ' ', + label='Recommended performance shaper: %s @ %.1f Hz' + % (performance_shaper.upper(), performance_shaper_freq), + ) + ax.plot( + freqs, psd * performance_shaper_vals, label='With %s applied' % (performance_shaper.upper()), color='cyan' + ) + ax2.plot( + [], + [], + ' ', + label='Recommended low vibrations shaper: %s @ %.1f Hz' % (lowvib_shaper.upper(), lowvib_shaper_freq), + ) ax.plot(freqs, psd * lowvib_shaper_vals, label='With %s applied' % (lowvib_shaper.upper()), color='lime') else: - ax2.plot([], [], ' ', label="Recommended best shaper: %s @ %.1f Hz" % (performance_shaper.upper(), performance_shaper_freq)) - ax.plot(freqs, psd * performance_shaper_vals, label='With %s applied' % (performance_shaper.upper()), color='cyan') + ax2.plot( + [], + [], + ' ', + label='Recommended best shaper: %s @ %.1f Hz' % (performance_shaper.upper(), performance_shaper_freq), + ) + ax.plot( + freqs, psd * performance_shaper_vals, label='With %s applied' % (performance_shaper.upper()), color='cyan' + ) # And the estimated damping ratio is finally added at the end of the legend - ax2.plot([], [], ' ', label="Estimated damping ratio (ζ): %.3f" % (zeta)) + ax2.plot([], [], ' ', label='Estimated damping ratio (ζ): %.3f' % (zeta)) # Draw the detected peaks and name them # This also draw the detection threshold and warning threshold (aka "effect zone") - ax.plot(peaks_freqs, psd[peaks], "x", color='black', markersize=8) + ax.plot(peaks_freqs, psd[peaks], 'x', color='black', markersize=8) for idx, peak in enumerate(peaks): if psd[peak] > peaks_threshold[1]: fontcolor = 'red' @@ -150,16 +213,28 @@ def plot_freq_response(ax, calibration_data, shapers, performance_shaper, peaks, else: fontcolor = 'black' fontweight = 'normal' - ax.annotate(f"{idx+1}", (freqs[peak], psd[peak]), - textcoords="offset points", xytext=(8, 5), - ha='left', fontsize=13, color=fontcolor, weight=fontweight) + ax.annotate( + f'{idx+1}', + (freqs[peak], psd[peak]), + textcoords='offset points', + xytext=(8, 5), + ha='left', + fontsize=13, + color=fontcolor, + weight=fontweight, + ) ax.axhline(y=peaks_threshold[0], color='black', linestyle='--', linewidth=0.5) ax.axhline(y=peaks_threshold[1], color='black', linestyle='--', linewidth=0.5) ax.fill_between(freqs, 0, peaks_threshold[0], color='green', alpha=0.15, label='Relax Region') ax.fill_between(freqs, peaks_threshold[0], peaks_threshold[1], color='orange', alpha=0.2, label='Warning Region') # Add the main resonant frequency and damping ratio of the axis to the graph title - ax.set_title("Axis Frequency Profile (ω0=%.1fHz, ζ=%.3f)" % (fr, zeta), fontsize=14, color=KLIPPAIN_COLORS['dark_orange'], weight='bold') + ax.set_title( + 'Axis Frequency Profile (ω0=%.1fHz, ζ=%.3f)' % (fr, zeta), + fontsize=14, + color=KLIPPAIN_COLORS['dark_orange'], + weight='bold', + ) ax.legend(loc='upper left', prop=fontP) ax2.legend(loc='upper right', prop=fontP) @@ -169,8 +244,8 @@ def plot_freq_response(ax, calibration_data, shapers, performance_shaper, peaks, # Plot a time-frequency spectrogram to see how the system respond over time during the # resonnance test. This can highlight hidden spots from the standard PSD graph from other harmonics def plot_spectrogram(ax, t, bins, pdata, peaks, max_freq): - ax.set_title("Time-Frequency Spectrogram", fontsize=14, color=KLIPPAIN_COLORS['dark_orange'], weight='bold') - + ax.set_title('Time-Frequency Spectrogram', fontsize=14, color=KLIPPAIN_COLORS['dark_orange'], weight='bold') + # We need to normalize the data to get a proper signal on the spectrogram # However, while using "LogNorm" provide too much background noise, using # "Normalize" make only the resonnance appearing and hide interesting elements @@ -182,19 +257,34 @@ def plot_spectrogram(ax, t, bins, pdata, peaks, max_freq): # save ~150-200MB of RAM during the "fig.savefig" operation. cm = 'inferno' norm = matplotlib.colors.LogNorm(vmin=vmin_value) - ax.imshow(pdata.T, norm=norm, cmap=cm, aspect='auto', extent=[t[0], t[-1], bins[0], bins[-1]], origin='lower', interpolation='antialiased') - - ax.set_xlim([0., max_freq]) + ax.imshow( + pdata.T, + norm=norm, + cmap=cm, + aspect='auto', + extent=[t[0], t[-1], bins[0], bins[-1]], + origin='lower', + interpolation='antialiased', + ) + + ax.set_xlim([0.0, max_freq]) ax.set_ylabel('Time (s)') ax.set_xlabel('Frequency (Hz)') - + # Add peaks lines in the spectrogram to get hint from peaks found in the first graph if peaks is not None: for idx, peak in enumerate(peaks): ax.axvline(peak, color='cyan', linestyle='dotted', linewidth=1) - ax.annotate(f"Peak {idx+1}", (peak, bins[-1]*0.9), - textcoords="data", color='cyan', rotation=90, fontsize=10, - verticalalignment='top', horizontalalignment='right') + ax.annotate( + f'Peak {idx+1}', + (peak, bins[-1] * 0.9), + textcoords='data', + color='cyan', + rotation=90, + fontsize=10, + verticalalignment='top', + horizontalalignment='right', + ) return @@ -203,7 +293,8 @@ def plot_spectrogram(ax, t, bins, pdata, peaks, max_freq): # Startup and main routines ###################################################################### -def shaper_calibration(lognames, klipperdir="~/klipper", max_smoothing=None, scv=5. , max_freq=200.): + +def shaper_calibration(lognames, klipperdir='~/klipper', max_smoothing=None, scv=5.0, max_freq=200.0, st_version=None): set_locale() global shaper_calibrate shaper_calibrate = setup_klipper_import(klipperdir) @@ -211,10 +302,12 @@ def shaper_calibration(lognames, klipperdir="~/klipper", max_smoothing=None, scv # Parse data datas = [parse_log(fn) for fn in lognames] if len(datas) > 1: - print_with_c_locale("Warning: incorrect number of .csv files detected. Only the first one will be used!") + print_with_c_locale('Warning: incorrect number of .csv files detected. Only the first one will be used!') # Compute shapers, PSD outputs and spectrogram - performance_shaper, shapers, calibration_data, fr, zeta = calibrate_shaper(datas[0], max_smoothing, scv, max_freq) + performance_shaper, shapers, calibration_data, fr, zeta, compat = calibrate_shaper( + datas[0], max_smoothing, scv, max_freq + ) pdata, bins, t = compute_spectrogram(datas[0]) del datas @@ -229,38 +322,51 @@ def shaper_calibration(lognames, klipperdir="~/klipper", max_smoothing=None, scv # Peak detection algorithm peaks_threshold = [ PEAKS_DETECTION_THRESHOLD * calibration_data.psd_sum.max(), - PEAKS_EFFECT_THRESHOLD * calibration_data.psd_sum.max() + PEAKS_EFFECT_THRESHOLD * calibration_data.psd_sum.max(), ] num_peaks, peaks, peaks_freqs = detect_peaks(calibration_data.psd_sum, calibration_data.freqs, peaks_threshold[0]) - + # Print the peaks info in the console - peak_freqs_formated = ["{:.1f}".format(f) for f in peaks_freqs] + peak_freqs_formated = ['{:.1f}'.format(f) for f in peaks_freqs] num_peaks_above_effect_threshold = np.sum(calibration_data.psd_sum[peaks] > peaks_threshold[1]) - print_with_c_locale("\nPeaks detected on the graph: %d @ %s Hz (%d above effect threshold)" % (num_peaks, ", ".join(map(str, peak_freqs_formated)), num_peaks_above_effect_threshold)) + print_with_c_locale( + '\nPeaks detected on the graph: %d @ %s Hz (%d above effect threshold)' + % (num_peaks, ', '.join(map(str, peak_freqs_formated)), num_peaks_above_effect_threshold) + ) # Create graph layout - fig, (ax1, ax2) = plt.subplots(2, 1, gridspec_kw={ - 'height_ratios':[4, 3], - 'bottom':0.050, - 'top':0.890, - 'left':0.085, - 'right':0.966, - 'hspace':0.169, - 'wspace':0.200 - }) + fig, (ax1, ax2) = plt.subplots( + 2, + 1, + gridspec_kw={ + 'height_ratios': [4, 3], + 'bottom': 0.050, + 'top': 0.890, + 'left': 0.085, + 'right': 0.966, + 'hspace': 0.169, + 'wspace': 0.200, + }, + ) fig.set_size_inches(8.3, 11.6) - + # Add a title with some test info - title_line1 = "INPUT SHAPER CALIBRATION TOOL" - fig.text(0.12, 0.965, title_line1, ha='left', va='bottom', fontsize=20, color=KLIPPAIN_COLORS['purple'], weight='bold') + title_line1 = 'INPUT SHAPER CALIBRATION TOOL' + fig.text( + 0.12, 0.965, title_line1, ha='left', va='bottom', fontsize=20, color=KLIPPAIN_COLORS['purple'], weight='bold' + ) try: filename_parts = (lognames[0].split('/')[-1]).split('_') - dt = datetime.strptime(f"{filename_parts[1]} {filename_parts[2]}", "%Y%m%d %H%M%S") + dt = datetime.strptime(f'{filename_parts[1]} {filename_parts[2]}', '%Y%m%d %H%M%S') title_line2 = dt.strftime('%x %X') + ' -- ' + filename_parts[3].upper().split('.')[0] + ' axis' - title_line3 = '| Square corner velocity: ' + str(scv) + 'mm/s' - title_line4 = '| Max allowed smoothing: ' + str(max_smoothing) - except: - print_with_c_locale("Warning: CSV filename look to be different than expected (%s)" % (lognames[0])) + if compat: + title_line3 = '| Compatibility mode with older Klipper,' + title_line4 = '| and no custom S&T parameters are used!' + else: + title_line3 = '| Square corner velocity: ' + str(scv) + 'mm/s' + title_line4 = '| Max allowed smoothing: ' + str(max_smoothing) + except Exception: + print_with_c_locale('Warning: CSV filename look to be different than expected (%s)' % (lognames[0])) title_line2 = lognames[0].split('/')[-1] title_line3 = '' title_line4 = '' @@ -269,7 +375,9 @@ def shaper_calibration(lognames, klipperdir="~/klipper", max_smoothing=None, scv fig.text(0.58, 0.946, title_line4, ha='left', va='top', fontsize=10, color=KLIPPAIN_COLORS['dark_purple']) # Plot the graphs - plot_freq_response(ax1, calibration_data, shapers, performance_shaper, peaks, peaks_freqs, peaks_threshold, fr, zeta, max_freq) + plot_freq_response( + ax1, calibration_data, shapers, performance_shaper, peaks, peaks_freqs, peaks_threshold, fr, zeta, max_freq + ) plot_spectrogram(ax2, t, bins, pdata, peaks_freqs, max_freq) # Adding a small Klippain logo to the top left corner of the figure @@ -278,8 +386,7 @@ def shaper_calibration(lognames, klipperdir="~/klipper", max_smoothing=None, scv ax_logo.axis('off') # Adding Shake&Tune version in the top right corner - st_version = get_git_version() - if st_version is not None: + if st_version != 'unknown': fig.text(0.995, 0.985, st_version, ha='right', va='bottom', fontsize=8, color=KLIPPAIN_COLORS['purple']) return fig @@ -287,25 +394,24 @@ def shaper_calibration(lognames, klipperdir="~/klipper", max_smoothing=None, scv def main(): # Parse command-line arguments - usage = "%prog [options] " + usage = '%prog [options] ' opts = optparse.OptionParser(usage) - opts.add_option("-o", "--output", type="string", dest="output", - default=None, help="filename of output graph") - opts.add_option("-f", "--max_freq", type="float", default=200., - help="maximum frequency to graph") - opts.add_option("-s", "--max_smoothing", type="float", default=None, - help="maximum shaper smoothing to allow") - opts.add_option("--scv", "--square_corner_velocity", type="float", - dest="scv", default=5., help="square corner velocity") - opts.add_option("-k", "--klipper_dir", type="string", dest="klipperdir", - default="~/klipper", help="main klipper directory") + opts.add_option('-o', '--output', type='string', dest='output', default=None, help='filename of output graph') + opts.add_option('-f', '--max_freq', type='float', default=200.0, help='maximum frequency to graph') + opts.add_option('-s', '--max_smoothing', type='float', default=None, help='maximum shaper smoothing to allow') + opts.add_option( + '--scv', '--square_corner_velocity', type='float', dest='scv', default=5.0, help='square corner velocity' + ) + opts.add_option( + '-k', '--klipper_dir', type='string', dest='klipperdir', default='~/klipper', help='main klipper directory' + ) options, args = opts.parse_args() if len(args) < 1: - opts.error("Incorrect number of arguments") + opts.error('Incorrect number of arguments') if options.output is None: - opts.error("You must specify an output file.png to use the script (option -o)") + opts.error('You must specify an output file.png to use the script (option -o)') if options.max_smoothing is not None and options.max_smoothing < 0.05: - opts.error("Too small max_smoothing specified (must be at least 0.05)") + opts.error('Too small max_smoothing specified (must be at least 0.05)') fig = shaper_calibration(args, options.klipperdir, options.max_smoothing, options.scv, options.max_freq) fig.savefig(options.output, dpi=150) diff --git a/src/graph_creators/graph_vibrations.py b/src/graph_creators/graph_vibrations.py new file mode 100644 index 0000000..ff551ee --- /dev/null +++ b/src/graph_creators/graph_vibrations.py @@ -0,0 +1,840 @@ +#!/usr/bin/env python3 + +################################################## +#### DIRECTIONAL VIBRATIONS PLOTTING SCRIPT ###### +################################################## +# Written by Frix_x#0161 # + +import math +import optparse +import os +import re +from collections import defaultdict +from datetime import datetime + +import matplotlib +import matplotlib.font_manager +import matplotlib.gridspec +import matplotlib.pyplot as plt +import matplotlib.ticker +import numpy as np + +matplotlib.use('Agg') + +from ..helpers.common_func import ( + compute_mechanical_parameters, + detect_peaks, + identify_low_energy_zones, + parse_log, + setup_klipper_import, +) +from ..helpers.locale_utils import print_with_c_locale, set_locale + +PEAKS_DETECTION_THRESHOLD = 0.05 +PEAKS_RELATIVE_HEIGHT_THRESHOLD = 0.04 +CURVE_SIMILARITY_SIGMOID_K = 0.5 +SPEEDS_VALLEY_DETECTION_THRESHOLD = 0.7 # Lower is more sensitive +SPEEDS_AROUND_PEAK_DELETION = 3 # to delete +-3mm/s around a peak +ANGLES_VALLEY_DETECTION_THRESHOLD = 1.1 # Lower is more sensitive + +KLIPPAIN_COLORS = { + 'purple': '#70088C', + 'orange': '#FF8D32', + 'dark_purple': '#150140', + 'dark_orange': '#F24130', + 'red_pink': '#F2055C', +} + + +###################################################################### +# Computation +###################################################################### + + +# Call to the official Klipper input shaper object to do the PSD computation +def calc_freq_response(data): + helper = shaper_calibrate.ShaperCalibrate(printer=None) + return helper.process_accelerometer_data(data) + + +# Calculate motor frequency profiles based on the measured Power Spectral Density (PSD) measurements for the machine kinematics +# main angles and then create a global motor profile as a weighted average (from their own vibrations) of all calculated profiles +def compute_motor_profiles(freqs, psds, all_angles_energy, measured_angles=None, energy_amplification_factor=2): + if measured_angles is None: + measured_angles = [0, 90] + + motor_profiles = {} + weighted_sum_profiles = np.zeros_like(freqs) + total_weight = 0 + conv_filter = np.ones(20) / 20 + + # Creating the PSD motor profiles for each angles + for angle in measured_angles: + # Calculate the sum of PSDs for the current angle and then convolve + sum_curve = np.sum(np.array([psds[angle][speed] for speed in psds[angle]]), axis=0) + motor_profiles[angle] = np.convolve(sum_curve / len(psds[angle]), conv_filter, mode='same') + + # Calculate weights + angle_energy = ( + all_angles_energy[angle] ** energy_amplification_factor + ) # First weighting factor is based on the total vibrations of the machine at the specified angle + curve_area = ( + np.trapz(motor_profiles[angle], freqs) ** energy_amplification_factor + ) # Additional weighting factor is based on the area under the current motor profile at this specified angle + total_angle_weight = angle_energy * curve_area + + # Update weighted sum profiles to get the global motor profile + weighted_sum_profiles += motor_profiles[angle] * total_angle_weight + total_weight += total_angle_weight + + # Creating a global average motor profile that is the weighted average of all the PSD motor profiles + global_motor_profile = weighted_sum_profiles / total_weight if total_weight != 0 else weighted_sum_profiles + + return motor_profiles, global_motor_profile + + +# Since it was discovered that there is no non-linear mixing in the stepper "steps" vibrations, instead of measuring +# the effects of each speeds at each angles, this function simplify it by using only the main motors axes (X/Y for Cartesian +# printers and A/B for CoreXY) measurements and project each points on the [0,360] degrees range using trigonometry +# to "sum" the vibration impact of each axis at every points of the generated spectrogram. The result is very similar at the end. +def compute_dir_speed_spectrogram(measured_speeds, data, kinematics='cartesian', measured_angles=None): + if measured_angles is None: + measured_angles = [0, 90] + + # We want to project the motor vibrations measured on their own axes on the [0, 360] range + spectrum_angles = np.linspace(0, 360, 720) # One point every 0.5 degrees + spectrum_speeds = np.linspace(min(measured_speeds), max(measured_speeds), len(measured_speeds) * 6) + spectrum_vibrations = np.zeros((len(spectrum_angles), len(spectrum_speeds))) + + def get_interpolated_vibrations(data, speed, speeds): + idx = np.clip(np.searchsorted(speeds, speed, side='left'), 1, len(speeds) - 1) + lower_speed = speeds[idx - 1] + upper_speed = speeds[idx] + lower_vibrations = data.get(lower_speed, 0) + upper_vibrations = data.get(upper_speed, 0) + return lower_vibrations + (speed - lower_speed) * (upper_vibrations - lower_vibrations) / ( + upper_speed - lower_speed + ) + + # Precompute trigonometric values and constant before the loop + angle_radians = np.deg2rad(spectrum_angles) + cos_vals = np.cos(angle_radians) + sin_vals = np.sin(angle_radians) + sqrt_2_inv = 1 / math.sqrt(2) + + # Compute the spectrum vibrations for each angle and speed combination + for target_angle_idx, (cos_val, sin_val) in enumerate(zip(cos_vals, sin_vals)): + for target_speed_idx, target_speed in enumerate(spectrum_speeds): + if kinematics == 'cartesian': + speed_1 = np.abs(target_speed * cos_val) + speed_2 = np.abs(target_speed * sin_val) + elif kinematics == 'corexy': + speed_1 = np.abs(target_speed * (cos_val + sin_val) * sqrt_2_inv) + speed_2 = np.abs(target_speed * (cos_val - sin_val) * sqrt_2_inv) + + vibrations_1 = get_interpolated_vibrations(data[measured_angles[0]], speed_1, measured_speeds) + vibrations_2 = get_interpolated_vibrations(data[measured_angles[1]], speed_2, measured_speeds) + spectrum_vibrations[target_angle_idx, target_speed_idx] = vibrations_1 + vibrations_2 + + return spectrum_angles, spectrum_speeds, spectrum_vibrations + + +def compute_angle_powers(spectrogram_data): + angles_powers = np.trapz(spectrogram_data, axis=1) + + # Since we want to plot it on a continuous polar plot later on, we need to append parts of + # the array to start and end of it to smooth transitions when doing the convolution + # and get the same value at modulo 360. Then we return the array without the extras + extended_angles_powers = np.concatenate([angles_powers[-9:], angles_powers, angles_powers[:9]]) + convolved_extended = np.convolve(extended_angles_powers, np.ones(15) / 15, mode='same') + + return convolved_extended[9:-9] + + +def compute_speed_powers(spectrogram_data, smoothing_window=15): + min_values = np.amin(spectrogram_data, axis=0) + max_values = np.amax(spectrogram_data, axis=0) + var_values = np.var(spectrogram_data, axis=0) + + # rescale the variance to the same range as max_values to plot it on the same graph + var_values = var_values / var_values.max() * max_values.max() + + # Create a vibration metric that is the product of the max values and the variance to quantify the best + # speeds that have at the same time a low global energy level that is also consistent at every angles + vibration_metric = max_values * var_values + + # utility function to pad and smooth the data avoiding edge effects + conv_filter = np.ones(smoothing_window) / smoothing_window + window = int(smoothing_window / 2) + + def pad_and_smooth(data): + data_padded = np.pad(data, (window,), mode='edge') + smoothed_data = np.convolve(data_padded, conv_filter, mode='valid') + return smoothed_data + + # Stack the arrays and apply padding and smoothing in batch + data_arrays = np.stack([min_values, max_values, var_values, vibration_metric]) + smoothed_arrays = np.array([pad_and_smooth(data) for data in data_arrays]) + + return smoothed_arrays + + +# Function that filter and split the good_speed ranges. The goal is to remove some zones around +# additional detected small peaks in order to suppress them if there is a peak, even if it's low, +# that's probably due to a crossing in the motor resonance pattern that still need to be removed +def filter_and_split_ranges(all_speeds, good_speeds, peak_speed_indices, deletion_range): + # Process each range to filter out and split based on peak indices + filtered_good_speeds = [] + for start, end, energy in good_speeds: + start_speed, end_speed = all_speeds[start], all_speeds[end] + # Identify peaks that intersect with the current speed range + intersecting_peaks_indices = [ + idx for speed, idx in peak_speed_indices.items() if start_speed <= speed <= end_speed + ] + + if not intersecting_peaks_indices: + filtered_good_speeds.append((start, end, energy)) + else: + intersecting_peaks_indices.sort() + current_start = start + + for peak_index in intersecting_peaks_indices: + before_peak_end = max(current_start, peak_index - deletion_range) + if current_start < before_peak_end: + filtered_good_speeds.append((current_start, before_peak_end, energy)) + current_start = peak_index + deletion_range + 1 + + if current_start < end: + filtered_good_speeds.append((current_start, end, energy)) + + # Sorting by start point once and then merge overlapping ranges + sorted_ranges = sorted(filtered_good_speeds, key=lambda x: x[0]) + merged_ranges = [sorted_ranges[0]] + + for current in sorted_ranges[1:]: + last_merged_start, last_merged_end, last_merged_energy = merged_ranges[-1] + if current[0] <= last_merged_end: + new_end = max(last_merged_end, current[1]) + new_energy = min(last_merged_energy, current[2]) + merged_ranges[-1] = (last_merged_start, new_end, new_energy) + else: + merged_ranges.append(current) + + return merged_ranges + + +# This function allow the computation of a symmetry score that reflect the spectrogram apparent symmetry between +# measured axes on both the shape of the signal and the energy level consistency across both side of the signal +def compute_symmetry_analysis(all_angles, spectrogram_data, measured_angles=None): + if measured_angles is None: + measured_angles = [0, 90] + + total_spectrogram_angles = len(all_angles) + half_spectrogram_angles = total_spectrogram_angles // 2 + + # Extend the spectrogram by adding half to the beginning (in order to not get an out of bounds error later) + extended_spectrogram = np.concatenate((spectrogram_data[-half_spectrogram_angles:], spectrogram_data), axis=0) + + # Calculate the split index directly within the slicing + midpoint_angle = np.mean(measured_angles) + split_index = int(midpoint_angle * (total_spectrogram_angles / 360) + half_spectrogram_angles) + half_segment_length = half_spectrogram_angles // 2 + + # Slice out the two segments of the spectrogram and flatten them for comparison + segment_1_flattened = extended_spectrogram[split_index - half_segment_length : split_index].flatten() + segment_2_flattened = extended_spectrogram[split_index : split_index + half_segment_length].flatten() + + # Compute the correlation coefficient between the two segments of spectrogram + correlation = np.corrcoef(segment_1_flattened, segment_2_flattened)[0, 1] + percentage_correlation_biased = (100 * np.power(correlation, 0.75)) + 10 + + return np.clip(0, 100, percentage_correlation_biased) + + +###################################################################### +# Graphing +###################################################################### + + +def plot_angle_profile_polar(ax, angles, angles_powers, low_energy_zones, symmetry_factor): + angles_radians = np.deg2rad(angles) + + ax.set_title('Polar angle energy profile', fontsize=14, color=KLIPPAIN_COLORS['dark_orange'], weight='bold') + ax.set_theta_zero_location('E') + ax.set_theta_direction(1) + + ax.plot(angles_radians, angles_powers, color=KLIPPAIN_COLORS['purple'], zorder=5) + ax.fill(angles_radians, angles_powers, color=KLIPPAIN_COLORS['purple'], alpha=0.3) + ax.set_xlim([0, np.deg2rad(360)]) + ymax = angles_powers.max() * 1.05 + ax.set_ylim([0, ymax]) + ax.set_thetagrids([theta * 15 for theta in range(360 // 15)]) + + ax.text( + 0, + 0, + f'Symmetry: {symmetry_factor:.1f}%', + ha='center', + va='center', + color=KLIPPAIN_COLORS['red_pink'], + fontsize=12, + fontweight='bold', + zorder=6, + ) + + for _, (start, end, _) in enumerate(low_energy_zones): + ax.axvline( + angles_radians[start], + angles_powers[start] / ymax, + color=KLIPPAIN_COLORS['red_pink'], + linestyle='dotted', + linewidth=1.5, + ) + ax.axvline( + angles_radians[end], + angles_powers[end] / ymax, + color=KLIPPAIN_COLORS['red_pink'], + linestyle='dotted', + linewidth=1.5, + ) + ax.fill_between( + angles_radians[start:end], angles_powers[start:end], angles_powers.max() * 1.05, color='green', alpha=0.2 + ) + + ax.xaxis.set_minor_locator(matplotlib.ticker.AutoMinorLocator()) + ax.yaxis.set_minor_locator(matplotlib.ticker.AutoMinorLocator()) + ax.grid(which='major', color='grey') + ax.grid(which='minor', color='lightgrey') + + # Polar plot doesn't follow the gridspec margin, so we adjust it manually here + pos = ax.get_position() + new_pos = [pos.x0 - 0.01, pos.y0 - 0.01, pos.width, pos.height] + ax.set_position(new_pos) + + return + + +def plot_global_speed_profile( + ax, + all_speeds, + sp_min_energy, + sp_max_energy, + sp_variance_energy, + vibration_metric, + num_peaks, + peaks, + low_energy_zones, +): + ax.set_title('Global speed energy profile', fontsize=14, color=KLIPPAIN_COLORS['dark_orange'], weight='bold') + ax.set_xlabel('Speed (mm/s)') + ax.set_ylabel('Energy') + ax2 = ax.twinx() + ax2.yaxis.set_visible(False) + + ax.plot(all_speeds, sp_min_energy, label='Minimum', color=KLIPPAIN_COLORS['dark_purple'], zorder=5) + ax.plot(all_speeds, sp_max_energy, label='Maximum', color=KLIPPAIN_COLORS['purple'], zorder=5) + ax.plot(all_speeds, sp_variance_energy, label='Variance', color=KLIPPAIN_COLORS['orange'], zorder=5, linestyle='--') + ax2.plot( + all_speeds, + vibration_metric, + label=f'Vibration metric ({num_peaks} bad peaks)', + color=KLIPPAIN_COLORS['red_pink'], + zorder=5, + ) + + ax.set_xlim([all_speeds.min(), all_speeds.max()]) + ax.set_ylim([0, sp_max_energy.max() * 1.15]) + + y2min = -(vibration_metric.max() * 0.025) + y2max = vibration_metric.max() * 1.07 + ax2.set_ylim([y2min, y2max]) + + if peaks is not None and len(peaks) > 0: + ax2.plot(all_speeds[peaks], vibration_metric[peaks], 'x', color='black', markersize=8, zorder=10) + for idx, peak in enumerate(peaks): + ax2.annotate( + f'{idx+1}', + (all_speeds[peak], vibration_metric[peak]), + textcoords='offset points', + xytext=(5, 5), + fontweight='bold', + ha='left', + fontsize=13, + color=KLIPPAIN_COLORS['red_pink'], + zorder=10, + ) + + for idx, (start, end, _) in enumerate(low_energy_zones): + # ax2.axvline(all_speeds[start], color=KLIPPAIN_COLORS['red_pink'], linestyle='dotted', linewidth=1.5, zorder=8) + # ax2.axvline(all_speeds[end], color=KLIPPAIN_COLORS['red_pink'], linestyle='dotted', linewidth=1.5, zorder=8) + ax2.fill_between( + all_speeds[start:end], + y2min, + vibration_metric[start:end], + color='green', + alpha=0.2, + label=f'Zone {idx+1}: {all_speeds[start]:.1f} to {all_speeds[end]:.1f} mm/s', + ) + + ax.xaxis.set_minor_locator(matplotlib.ticker.AutoMinorLocator()) + ax.yaxis.set_minor_locator(matplotlib.ticker.AutoMinorLocator()) + ax.grid(which='major', color='grey') + ax.grid(which='minor', color='lightgrey') + + fontP = matplotlib.font_manager.FontProperties() + fontP.set_size('small') + ax.legend(loc='upper left', prop=fontP) + ax2.legend(loc='upper right', prop=fontP) + + return + + +def plot_angular_speed_profiles(ax, speeds, angles, spectrogram_data, kinematics='cartesian'): + ax.set_title('Angular speed energy profiles', fontsize=14, color=KLIPPAIN_COLORS['dark_orange'], weight='bold') + ax.set_xlabel('Speed (mm/s)') + ax.set_ylabel('Energy') + + # Define mappings for labels and colors to simplify plotting commands + angle_settings = { + 0: ('X (0 deg)', 'purple', 10), + 90: ('Y (90 deg)', 'dark_purple', 5), + 45: ('A (45 deg)' if kinematics == 'corexy' else '45 deg', 'orange', 10), + 135: ('B (135 deg)' if kinematics == 'corexy' else '135 deg', 'dark_orange', 5), + } + + # Plot each angle using settings from the dictionary + for angle, (label, color, zorder) in angle_settings.items(): + idx = np.searchsorted(angles, angle, side='left') + ax.plot(speeds, spectrogram_data[idx], label=label, color=KLIPPAIN_COLORS[color], zorder=zorder) + + ax.set_xlim([speeds.min(), speeds.max()]) + max_value = max(spectrogram_data[angle].max() for angle in [0, 45, 90, 135]) + ax.set_ylim([0, max_value * 1.1]) + + ax.xaxis.set_minor_locator(matplotlib.ticker.AutoMinorLocator()) + ax.yaxis.set_minor_locator(matplotlib.ticker.AutoMinorLocator()) + ax.grid(which='major', color='grey') + ax.grid(which='minor', color='lightgrey') + + fontP = matplotlib.font_manager.FontProperties() + fontP.set_size('small') + ax.legend(loc='upper right', prop=fontP) + + return + + +def plot_motor_profiles(ax, freqs, main_angles, motor_profiles, global_motor_profile, max_freq): + ax.set_title('Motor frequency profile', fontsize=14, color=KLIPPAIN_COLORS['dark_orange'], weight='bold') + ax.set_ylabel('Energy') + ax.set_xlabel('Frequency (Hz)') + + ax2 = ax.twinx() + ax2.yaxis.set_visible(False) + + # Global weighted average motor profile + ax.plot(freqs, global_motor_profile, label='Combined', color=KLIPPAIN_COLORS['purple'], zorder=5) + max_value = global_motor_profile.max() + + # Mapping of angles to axis names + angle_settings = {0: 'X', 90: 'Y', 45: 'A', 135: 'B'} + + # And then plot the motor profiles at each measured angles + for angle in main_angles: + profile_max = motor_profiles[angle].max() + if profile_max > max_value: + max_value = profile_max + label = f'{angle_settings[angle]} ({angle} deg)' if angle in angle_settings else f'{angle} deg' + ax.plot(freqs, motor_profiles[angle], linestyle='--', label=label, zorder=2) + + ax.set_xlim([0, max_freq]) + ax.set_ylim([0, max_value * 1.1]) + ax.ticklabel_format(axis='y', style='scientific', scilimits=(0, 0)) + + # Then add the motor resonance peak to the graph and print some infos about it + motor_fr, motor_zeta, motor_res_idx, lowfreq_max = compute_mechanical_parameters(global_motor_profile, freqs, 30) + if lowfreq_max: + print_with_c_locale( + '[WARNING] There are a lot of low frequency vibrations that can alter the readings. This is probably due to the test being performed at too high an acceleration!' + ) + print_with_c_locale( + 'Try lowering the ACCEL value and/or increasing the SIZE value before restarting the macro to ensure that only constant speeds are being recorded and that the dynamic behavior of the machine is not affecting the measurements' + ) + if motor_zeta is not None: + print_with_c_locale( + 'Motors have a main resonant frequency at %.1fHz with an estimated damping ratio of %.3f' + % (motor_fr, motor_zeta) + ) + else: + print_with_c_locale( + 'Motors have a main resonant frequency at %.1fHz but it was impossible to estimate a damping ratio.' + % (motor_fr) + ) + + ax.plot(freqs[motor_res_idx], global_motor_profile[motor_res_idx], 'x', color='black', markersize=10) + ax.annotate( + 'R', + (freqs[motor_res_idx], global_motor_profile[motor_res_idx]), + textcoords='offset points', + xytext=(15, 5), + ha='right', + fontsize=14, + color=KLIPPAIN_COLORS['red_pink'], + weight='bold', + ) + + ax2.plot([], [], ' ', label='Motor resonant frequency (ω0): %.1fHz' % (motor_fr)) + if motor_zeta is not None: + ax2.plot([], [], ' ', label='Motor damping ratio (ζ): %.3f' % (motor_zeta)) + else: + ax2.plot([], [], ' ', label='No damping ratio computed') + + ax.xaxis.set_minor_locator(matplotlib.ticker.AutoMinorLocator()) + ax.yaxis.set_minor_locator(matplotlib.ticker.AutoMinorLocator()) + ax.grid(which='major', color='grey') + ax.grid(which='minor', color='lightgrey') + + fontP = matplotlib.font_manager.FontProperties() + fontP.set_size('small') + ax.legend(loc='upper left', prop=fontP) + ax2.legend(loc='upper right', prop=fontP) + + return + + +def plot_vibration_spectrogram_polar(ax, angles, speeds, spectrogram_data): + angles_radians = np.radians(angles) + + # Assuming speeds defines the radial distance from the center, we need to create a meshgrid + # for both angles and speeds to map the spectrogram data onto a polar plot correctly + radius, theta = np.meshgrid(speeds, angles_radians) + + ax.set_title( + 'Polar vibrations heatmap', fontsize=14, color=KLIPPAIN_COLORS['dark_orange'], weight='bold', va='bottom' + ) + ax.set_theta_zero_location('E') + ax.set_theta_direction(1) + + ax.pcolormesh(theta, radius, spectrogram_data, norm=matplotlib.colors.LogNorm(), cmap='inferno', shading='auto') + ax.set_thetagrids([theta * 15 for theta in range(360 // 15)]) + ax.tick_params(axis='y', which='both', colors='white', labelsize='medium') + ax.set_ylim([0, max(speeds)]) + + # Polar plot doesn't follow the gridspec margin, so we adjust it manually here + pos = ax.get_position() + new_pos = [pos.x0 - 0.01, pos.y0 - 0.01, pos.width, pos.height] + ax.set_position(new_pos) + + return + + +def plot_vibration_spectrogram(ax, angles, speeds, spectrogram_data, peaks): + ax.set_title('Vibrations heatmap', fontsize=14, color=KLIPPAIN_COLORS['dark_orange'], weight='bold') + ax.set_xlabel('Speed (mm/s)') + ax.set_ylabel('Angle (deg)') + + ax.imshow( + spectrogram_data, + norm=matplotlib.colors.LogNorm(), + cmap='inferno', + aspect='auto', + extent=[speeds[0], speeds[-1], angles[0], angles[-1]], + origin='lower', + interpolation='antialiased', + ) + + # Add peaks lines in the spectrogram to get hint from peaks found in the first graph + if peaks is not None and len(peaks) > 0: + for idx, peak in enumerate(peaks): + ax.axvline(speeds[peak], color='cyan', linewidth=0.75) + ax.annotate( + f'Peak {idx+1}', + (speeds[peak], angles[-1] * 0.9), + textcoords='data', + color='cyan', + rotation=90, + fontsize=10, + verticalalignment='top', + horizontalalignment='right', + ) + + return + + +def plot_motor_config_txt(fig, motors, differences): + motor_details = [(motors[0], 'X motor'), (motors[1], 'Y motor')] + + distance = 0.12 + if motors[0].get_property('autotune_enabled'): + distance = 0.24 + config_blocks = [ + f"| {lbl}: {mot.get_property('motor').upper()} on {mot.get_property('tmc').upper()} @ {mot.get_property('voltage')}V {mot.get_property('run_current')}A" + for mot, lbl in motor_details + ] + config_blocks.append('| TMC Autotune enabled') + else: + config_blocks = [ + f"| {lbl}: {mot.get_property('tmc').upper()} @ {mot.get_property('run_current')}A" + for mot, lbl in motor_details + ] + config_blocks.append('| TMC Autotune not detected') + + for idx, block in enumerate(config_blocks): + fig.text( + 0.40, 0.990 - 0.015 * idx, block, ha='left', va='top', fontsize=10, color=KLIPPAIN_COLORS['dark_purple'] + ) + + tmc_registers = motors[0].get_registers() + idx = -1 + for idx, (register, settings) in enumerate(tmc_registers.items()): + settings_str = ' '.join(f'{k}={v}' for k, v in settings.items()) + tmc_block = f'| {register.upper()}: {settings_str}' + fig.text( + 0.40 + distance, + 0.990 - 0.015 * idx, + tmc_block, + ha='left', + va='top', + fontsize=10, + color=KLIPPAIN_COLORS['dark_purple'], + ) + + if differences is not None: + differences_text = f'| Y motor diff: {differences}' + fig.text( + 0.40 + distance, + 0.990 - 0.015 * (idx + 1), + differences_text, + ha='left', + va='top', + fontsize=10, + color=KLIPPAIN_COLORS['dark_purple'], + ) + + +###################################################################### +# Startup and main routines +###################################################################### + + +def extract_angle_and_speed(logname): + try: + match = re.search(r'an(\d+)_\d+sp(\d+)_\d+', os.path.basename(logname)) + if match: + angle = match.group(1) + speed = match.group(2) + else: + raise ValueError(f'File {logname} does not match expected format. Clean your /tmp folder and start again!') + except AttributeError as err: + raise ValueError( + f'File {logname} does not match expected format. Clean your /tmp folder and start again!' + ) from err + return float(angle), float(speed) + + +def vibrations_profile( + lognames, klipperdir='~/klipper', kinematics='cartesian', accel=None, max_freq=1000.0, st_version=None, motors=None +): + set_locale() + global shaper_calibrate + shaper_calibrate = setup_klipper_import(klipperdir) + + if kinematics == 'cartesian': + main_angles = [0, 90] + elif kinematics == 'corexy': + main_angles = [45, 135] + else: + raise ValueError('Only Cartesian and CoreXY kinematics are supported by this tool at the moment!') + + psds = defaultdict(lambda: defaultdict(list)) + psds_sum = defaultdict(lambda: defaultdict(list)) + target_freqs_initialized = False + + for logname in lognames: + data = parse_log(logname) + angle, speed = extract_angle_and_speed(logname) + freq_response = calc_freq_response(data) + first_freqs = freq_response.freq_bins + psd_sum = freq_response.psd_sum + + if not target_freqs_initialized: + target_freqs = first_freqs[first_freqs <= max_freq] + target_freqs_initialized = True + + psd_sum = psd_sum[first_freqs <= max_freq] + first_freqs = first_freqs[first_freqs <= max_freq] + + # Store the interpolated PSD and integral values + psds[angle][speed] = np.interp(target_freqs, first_freqs, psd_sum) + psds_sum[angle][speed] = np.trapz(psd_sum, first_freqs) + + measured_angles = sorted(psds_sum.keys()) + measured_speeds = sorted({speed for angle_speeds in psds_sum.values() for speed in angle_speeds.keys()}) + + for main_angle in main_angles: + if main_angle not in measured_angles: + raise ValueError('Measurements not taken at the correct angles for the specified kinematics!') + + # Precompute the variables used in plot functions + all_angles, all_speeds, spectrogram_data = compute_dir_speed_spectrogram( + measured_speeds, psds_sum, kinematics, main_angles + ) + all_angles_energy = compute_angle_powers(spectrogram_data) + sp_min_energy, sp_max_energy, sp_variance_energy, vibration_metric = compute_speed_powers(spectrogram_data) + motor_profiles, global_motor_profile = compute_motor_profiles(target_freqs, psds, all_angles_energy, main_angles) + + # symmetry_factor = compute_symmetry_analysis(all_angles, all_angles_energy) + symmetry_factor = compute_symmetry_analysis(all_angles, spectrogram_data, main_angles) + print_with_c_locale(f'Machine estimated vibration symmetry: {symmetry_factor:.1f}%') + + # Analyze low variance ranges of vibration energy across all angles for each speed to identify clean speeds + # and highlight them. Also find the peaks to identify speeds to avoid due to high resonances + num_peaks, vibration_peaks, peaks_speeds = detect_peaks( + vibration_metric, + all_speeds, + PEAKS_DETECTION_THRESHOLD * vibration_metric.max(), + PEAKS_RELATIVE_HEIGHT_THRESHOLD, + 10, + 10, + ) + formated_peaks_speeds = ['{:.1f}'.format(pspeed) for pspeed in peaks_speeds] + print_with_c_locale( + 'Vibrations peaks detected: %d @ %s mm/s (avoid setting a speed near these values in your slicer print profile)' + % (num_peaks, ', '.join(map(str, formated_peaks_speeds))) + ) + + good_speeds = identify_low_energy_zones(vibration_metric, SPEEDS_VALLEY_DETECTION_THRESHOLD) + if good_speeds is not None: + deletion_range = int(SPEEDS_AROUND_PEAK_DELETION / (all_speeds[1] - all_speeds[0])) + peak_speed_indices = {pspeed: np.where(all_speeds == pspeed)[0][0] for pspeed in set(peaks_speeds)} + + # Filter and split ranges based on peak indices, avoiding overlaps + good_speeds = filter_and_split_ranges(all_speeds, good_speeds, peak_speed_indices, deletion_range) + + # Add some logging about the good speeds found + print_with_c_locale(f'Lowest vibrations speeds ({len(good_speeds)} ranges sorted from best to worse):') + for idx, (start, end, _) in enumerate(good_speeds): + print_with_c_locale(f'{idx+1}: {all_speeds[start]:.1f} to {all_speeds[end]:.1f} mm/s') + + # Angle low energy valleys identification (good angles ranges) and print them to the console + good_angles = identify_low_energy_zones(all_angles_energy, ANGLES_VALLEY_DETECTION_THRESHOLD) + if good_angles is not None: + print_with_c_locale(f'Lowest vibrations angles ({len(good_angles)} ranges sorted from best to worse):') + for idx, (start, end, energy) in enumerate(good_angles): + print_with_c_locale( + f'{idx+1}: {all_angles[start]:.1f}° to {all_angles[end]:.1f}° (mean vibrations energy: {energy:.2f}% of max)' + ) + + # Create graph layout + fig, ((ax1, ax2, ax3), (ax4, ax5, ax6)) = plt.subplots( + 2, + 3, + gridspec_kw={ + 'height_ratios': [1, 1], + 'width_ratios': [4, 8, 6], + 'bottom': 0.050, + 'top': 0.890, + 'left': 0.040, + 'right': 0.985, + 'hspace': 0.166, + 'wspace': 0.138, + }, + ) + + # Transform ax3 and ax4 to polar plots + ax1.remove() + ax1 = fig.add_subplot(2, 3, 1, projection='polar') + ax4.remove() + ax4 = fig.add_subplot(2, 3, 4, projection='polar') + + # Set the global .png figure size + fig.set_size_inches(20, 11.5) + + # Add title + title_line1 = 'MACHINE VIBRATIONS ANALYSIS TOOL' + fig.text( + 0.060, 0.965, title_line1, ha='left', va='bottom', fontsize=20, color=KLIPPAIN_COLORS['purple'], weight='bold' + ) + try: + filename_parts = (lognames[0].split('/')[-1]).split('_') + dt = datetime.strptime(f"{filename_parts[1]} {filename_parts[2].split('-')[0]}", '%Y%m%d %H%M%S') + title_line2 = dt.strftime('%x %X') + if accel is not None: + title_line2 += ' at ' + str(accel) + ' mm/s² -- ' + kinematics.upper() + ' kinematics' + except Exception: + print_with_c_locale('Warning: CSV filenames appear to be different than expected (%s)' % (lognames[0])) + title_line2 = lognames[0].split('/')[-1] + fig.text(0.060, 0.957, title_line2, ha='left', va='top', fontsize=16, color=KLIPPAIN_COLORS['dark_purple']) + + # Add the motors infos to the top of the graph + if motors is not None and len(motors) == 2: + differences = motors[0].compare_to(motors[1]) + plot_motor_config_txt(fig, motors, differences) + if differences is not None and kinematics == 'corexy': + print_with_c_locale(f'Warning: motors have different TMC configurations!\n{differences}') + + # Plot the graphs + plot_angle_profile_polar(ax1, all_angles, all_angles_energy, good_angles, symmetry_factor) + plot_vibration_spectrogram_polar(ax4, all_angles, all_speeds, spectrogram_data) + + plot_global_speed_profile( + ax2, + all_speeds, + sp_min_energy, + sp_max_energy, + sp_variance_energy, + vibration_metric, + num_peaks, + vibration_peaks, + good_speeds, + ) + plot_angular_speed_profiles(ax3, all_speeds, all_angles, spectrogram_data, kinematics) + plot_vibration_spectrogram(ax5, all_angles, all_speeds, spectrogram_data, vibration_peaks) + + plot_motor_profiles(ax6, target_freqs, main_angles, motor_profiles, global_motor_profile, max_freq) + + # Adding a small Klippain logo to the top left corner of the figure + ax_logo = fig.add_axes([0.001, 0.924, 0.075, 0.075], anchor='NW') + ax_logo.imshow(plt.imread(os.path.join(os.path.dirname(os.path.abspath(__file__)), 'klippain.png'))) + ax_logo.axis('off') + + # Adding Shake&Tune version in the top right corner + if st_version != 'unknown': + fig.text(0.995, 0.985, st_version, ha='right', va='bottom', fontsize=8, color=KLIPPAIN_COLORS['purple']) + + return fig + + +def main(): + # Parse command-line arguments + usage = '%prog [options] ' + opts = optparse.OptionParser(usage) + opts.add_option('-o', '--output', type='string', dest='output', default=None, help='filename of output graph') + opts.add_option( + '-c', '--accel', type='int', dest='accel', default=None, help='accel value to be printed on the graph' + ) + opts.add_option('-f', '--max_freq', type='float', default=1000.0, help='maximum frequency to graph') + opts.add_option( + '-k', '--klipper_dir', type='string', dest='klipperdir', default='~/klipper', help='main klipper directory' + ) + opts.add_option( + '-m', + '--kinematics', + type='string', + dest='kinematics', + default='cartesian', + help='machine kinematics configuration', + ) + options, args = opts.parse_args() + if len(args) < 1: + opts.error('No CSV file(s) to analyse') + if options.output is None: + opts.error('You must specify an output file.png to use the script (option -o)') + if options.kinematics not in ['cartesian', 'corexy']: + opts.error('Only cartesian and corexy kinematics are supported by this tool at the moment!') + + fig = vibrations_profile(args, options.klipperdir, options.kinematics, options.accel, options.max_freq) + fig.savefig(options.output, dpi=150) + + +if __name__ == '__main__': + main() diff --git a/K-ShakeTune/scripts/klippain.png b/src/graph_creators/klippain.png similarity index 100% rename from K-ShakeTune/scripts/klippain.png rename to src/graph_creators/klippain.png diff --git a/src/helpers/__init__.py b/src/helpers/__init__.py new file mode 100644 index 0000000..e69de29 diff --git a/src/helpers/common_func.py b/src/helpers/common_func.py new file mode 100644 index 0000000..14b9a24 --- /dev/null +++ b/src/helpers/common_func.py @@ -0,0 +1,228 @@ +#!/usr/bin/env python3 + +# Common functions for the Shake&Tune package +# Written by Frix_x#0161 # + +import math +import os +import sys +from importlib import import_module +from pathlib import Path + +import numpy as np +from git import GitCommandError, Repo +from scipy.signal import spectrogram + + +def parse_log(logname): + with open(logname) as f: + for header in f: + if not header.startswith('#'): + break + if not header.startswith('freq,psd_x,psd_y,psd_z,psd_xyz'): + # Raw accelerometer data + return np.loadtxt(logname, comments='#', delimiter=',') + # Power spectral density data or shaper calibration data + raise ValueError( + 'File %s does not contain raw accelerometer data and therefore ' + 'is not supported by Shake&Tune. Please use the official Klipper ' + 'script to process it instead.' % (logname,) + ) + + +def setup_klipper_import(kdir): + kdir = os.path.expanduser(kdir) + sys.path.append(os.path.join(kdir, 'klippy')) + return import_module('.shaper_calibrate', 'extras') + + +# This is used to print the current S&T version on top of the png graph file +def get_git_version(): + try: + # Get the absolute path of the script, resolving any symlinks + # Then get 2 times to parent dir to be at the git root folder + script_path = Path(__file__).resolve() + repo_path = script_path.parents[1] + repo = Repo(repo_path) + + try: + version = repo.git.describe('--tags') + except GitCommandError: + # If no tag is found, use the simplified commit SHA instead + version = repo.head.commit.hexsha[:7] + return version + + except Exception: + return None + + +# This is Klipper's spectrogram generation function adapted to use Scipy +def compute_spectrogram(data): + N = data.shape[0] + Fs = N / (data[-1, 0] - data[0, 0]) + # Round up to a power of 2 for faster FFT + M = 1 << int(0.5 * Fs - 1).bit_length() + window = np.kaiser(M, 6.0) + + def _specgram(x): + return spectrogram( + x, fs=Fs, window=window, nperseg=M, noverlap=M // 2, detrend='constant', scaling='density', mode='psd' + ) + + d = {'x': data[:, 1], 'y': data[:, 2], 'z': data[:, 3]} + f, t, pdata = _specgram(d['x']) + for axis in 'yz': + pdata += _specgram(d[axis])[2] + return pdata, t, f + + +# Compute natural resonant frequency and damping ratio by using the half power bandwidth method with interpolated frequencies +def compute_mechanical_parameters(psd, freqs, min_freq=None): + max_under_min_freq = False + + if min_freq is not None: + min_freq_index = np.searchsorted(freqs, min_freq, side='left') + if min_freq_index >= len(freqs): + return None, None, None, max_under_min_freq + if np.argmax(psd) < min_freq_index: + max_under_min_freq = True + else: + min_freq_index = 0 + + # Consider only the part of the signal above min_freq + psd_above_min_freq = psd[min_freq_index:] + if len(psd_above_min_freq) == 0: + return None, None, None, max_under_min_freq + + max_power_index_above_min_freq = np.argmax(psd_above_min_freq) + max_power_index = max_power_index_above_min_freq + min_freq_index + fr = freqs[max_power_index] + max_power = psd[max_power_index] + + half_power = max_power / math.sqrt(2) + indices_below = np.where(psd[:max_power_index] <= half_power)[0] + indices_above = np.where(psd[max_power_index:] <= half_power)[0] + + # If we are not able to find points around the half power, we can't compute the damping ratio and return None instead + if len(indices_below) == 0 or len(indices_above) == 0: + return fr, None, max_power_index, max_under_min_freq + + idx_below = indices_below[-1] + idx_above = indices_above[0] + max_power_index + freq_below_half_power = freqs[idx_below] + (half_power - psd[idx_below]) * ( + freqs[idx_below + 1] - freqs[idx_below] + ) / (psd[idx_below + 1] - psd[idx_below]) + freq_above_half_power = freqs[idx_above - 1] + (half_power - psd[idx_above - 1]) * ( + freqs[idx_above] - freqs[idx_above - 1] + ) / (psd[idx_above] - psd[idx_above - 1]) + + bandwidth = freq_above_half_power - freq_below_half_power + bw1 = math.pow(bandwidth / fr, 2) + bw2 = math.pow(bandwidth / fr, 4) + + try: + zeta = math.sqrt(0.5 - math.sqrt(1 / (4 + 4 * bw1 - bw2))) + except ValueError: + # If a math problem arise such as a negative sqrt term, we also return None instead for damping ratio + return fr, None, max_power_index, max_under_min_freq + + return fr, zeta, max_power_index, max_under_min_freq + + +# This find all the peaks in a curve by looking at when the derivative term goes from positive to negative +# Then only the peaks found above a threshold are kept to avoid capturing peaks in the low amplitude noise of a signal +def detect_peaks(data, indices, detection_threshold, relative_height_threshold=None, window_size=5, vicinity=3): + # Smooth the curve using a moving average to avoid catching peaks everywhere in noisy signals + kernel = np.ones(window_size) / window_size + smoothed_data = np.convolve(data, kernel, mode='valid') + mean_pad = [np.mean(data[:window_size])] * (window_size // 2) + smoothed_data = np.concatenate((mean_pad, smoothed_data)) + + # Find peaks on the smoothed curve + smoothed_peaks = ( + np.where((smoothed_data[:-2] < smoothed_data[1:-1]) & (smoothed_data[1:-1] > smoothed_data[2:]))[0] + 1 + ) + smoothed_peaks = smoothed_peaks[smoothed_data[smoothed_peaks] > detection_threshold] + + # Additional validation for peaks based on relative height + valid_peaks = smoothed_peaks + if relative_height_threshold is not None: + valid_peaks = [] + for peak in smoothed_peaks: + peak_height = smoothed_data[peak] - np.min( + smoothed_data[max(0, peak - vicinity) : min(len(smoothed_data), peak + vicinity + 1)] + ) + if peak_height > relative_height_threshold * smoothed_data[peak]: + valid_peaks.append(peak) + + # Refine peak positions on the original curve + refined_peaks = [] + for peak in valid_peaks: + local_max = peak + np.argmax(data[max(0, peak - vicinity) : min(len(data), peak + vicinity + 1)]) - vicinity + refined_peaks.append(local_max) + + num_peaks = len(refined_peaks) + + return num_peaks, np.array(refined_peaks), indices[refined_peaks] + + +# The goal is to find zone outside of peaks (flat low energy zones) in a signal +def identify_low_energy_zones(power_total, detection_threshold=0.1): + valleys = [] + + # Calculate the a "mean + 1/4" and standard deviation of the entire power_total + mean_energy = np.mean(power_total) + (np.max(power_total) - np.min(power_total)) / 4 + std_energy = np.std(power_total) + + # Define a threshold value as "mean + 1/4" minus a certain number of standard deviations + threshold_value = mean_energy - detection_threshold * std_energy + + # Find valleys in power_total based on the threshold + in_valley = False + start_idx = 0 + for i, value in enumerate(power_total): + if not in_valley and value < threshold_value: + in_valley = True + start_idx = i + elif in_valley and value >= threshold_value: + in_valley = False + valleys.append((start_idx, i)) + + # If the last point is still in a valley, close the valley + if in_valley: + valleys.append((start_idx, len(power_total) - 1)) + + max_signal = np.max(power_total) + + # Calculate mean energy for each valley as a percentage of the maximum of the signal + valley_means_percentage = [] + for start, end in valleys: + if not np.isnan(np.mean(power_total[start:end])): + valley_means_percentage.append((start, end, (np.mean(power_total[start:end]) / max_signal) * 100)) + + # Sort valleys based on mean percentage values + sorted_valleys = sorted(valley_means_percentage, key=lambda x: x[2]) + + return sorted_valleys + + +# Calculate or estimate a "similarity" factor between two PSD curves and scale it to a percentage. This is +# used here to quantify how close the two belts path behavior and responses are close together. +def compute_curve_similarity_factor(x1, y1, x2, y2, sim_sigmoid_k=0.6): + # Interpolate PSDs to match the same frequency bins and do a cross-correlation + y2_interp = np.interp(x1, x2, y2) + cross_corr = np.correlate(y1, y2_interp, mode='full') + + # Find the peak of the cross-correlation and compute a similarity normalized by the energy of the signals + peak_value = np.max(cross_corr) + similarity = peak_value / (np.sqrt(np.sum(y1**2) * np.sum(y2_interp**2))) + + # Apply sigmoid scaling to get better numbers and get a final percentage value + scaled_similarity = sigmoid_scale(-np.log(1 - similarity), sim_sigmoid_k) + + return scaled_similarity + + +# Simple helper to compute a sigmoid scalling (from 0 to 100%) +def sigmoid_scale(x, k=1): + return 1 / (1 + np.exp(-k * x)) * 100 diff --git a/src/helpers/filemanager.py b/src/helpers/filemanager.py new file mode 100644 index 0000000..9ac2d75 --- /dev/null +++ b/src/helpers/filemanager.py @@ -0,0 +1,38 @@ +#!/usr/bin/env python3 + +# Common file management functions for the Shake&Tune package +# Written by Frix_x#0161 # + +import os +import time +from pathlib import Path + + +def wait_file_ready(filepath: Path, timeout: int = 60) -> None: + file_busy = True + loop_count = 0 + + while file_busy: + if loop_count >= timeout: + raise TimeoutError(f'Klipper is taking too long to release the CSV file ({filepath})!') + + # Try to open the file in write-only mode to check if it is in use + # If we successfully open and close the file, it is not in use + try: + fd = os.open(filepath, os.O_WRONLY) + os.close(fd) + file_busy = False + except OSError: + # If OSError is caught, it indicates the file is still being used + pass + except Exception: + # If another exception is raised, it's not a problem, we just loop again + pass + + loop_count += 1 + time.sleep(1) + + +def ensure_folders_exist(folders: list[Path]) -> None: + for folder in folders: + folder.mkdir(parents=True, exist_ok=True) diff --git a/K-ShakeTune/scripts/locale_utils.py b/src/helpers/locale_utils.py similarity index 73% rename from K-ShakeTune/scripts/locale_utils.py rename to src/helpers/locale_utils.py index ef4018c..611ecbd 100644 --- a/K-ShakeTune/scripts/locale_utils.py +++ b/src/helpers/locale_utils.py @@ -6,6 +6,7 @@ import locale + # Set the best locale for time and date formating (generation of the titles) def set_locale(): try: @@ -15,16 +16,19 @@ def set_locale(): except locale.Error: locale.setlocale(locale.LC_TIME, 'C') + # Print function to avoid problem in Klipper console (that doesn't support special characters) due to locale settings def print_with_c_locale(*args, **kwargs): try: original_locale = locale.getlocale() locale.setlocale(locale.LC_ALL, 'C') except locale.Error as e: - print("Warning: Failed to set a basic locale. Special characters may not display correctly in Klipper console:", e) + print( + 'Warning: Failed to set a basic locale. Special characters may not display correctly in Klipper console:', e + ) finally: - print(*args, **kwargs) # Proceed with printing regardless of locale setting success + print(*args, **kwargs) # Proceed with printing regardless of locale setting success try: locale.setlocale(locale.LC_ALL, original_locale) except locale.Error as e: - print("Warning: Failed to restore the original locale setting:", e) + print('Warning: Failed to restore the original locale setting:', e) diff --git a/src/helpers/motorlogparser.py b/src/helpers/motorlogparser.py new file mode 100644 index 0000000..4e6e743 --- /dev/null +++ b/src/helpers/motorlogparser.py @@ -0,0 +1,205 @@ +#!/usr/bin/env python3 + +# Classes to parse the Klipper log and parse the TMC dump to extract the relevant information +# Written by Frix_x#0161 # + +import re +from pathlib import Path +from typing import Any, Dict, List, Optional, Union + + +class Motor: + def __init__(self, name: str): + self._name: str = name + self._registers: Dict[str, Dict[str, Any]] = {} + self._properties: Dict[str, Any] = {} + + def set_register(self, register: str, value: Any) -> None: + # Special parsing for CHOPCONF to extract meaningful values + if register == 'CHOPCONF': + # Add intpol=0 if missing from the register dump + if 'intpol=' not in value: + value += ' intpol=0' + # Simplify the microstep resolution format + mres_match = re.search(r'mres=\d+\((\d+)usteps\)', value) + if mres_match: + value = re.sub(r'mres=\d+\(\d+usteps\)', f'mres={mres_match.group(1)}', value) + + # Special parsing for CHOPCONF to avoid pwm_ before each values + if register == 'PWMCONF': + parts = value.split() + new_parts = [] + for part in parts: + key, val = part.split('=', 1) + if key.startswith('pwm_'): + key = key[4:] + new_parts.append(f'{key}={val}') + value = ' '.join(new_parts) + + # General cleaning to remove extraneous labels and colons and parse the whole into Motor _registers + cleaned_values = re.sub(r'\b\w+:\s+\S+\s+', '', value) + + # Then fill the registers while merging all the thresholds into the same THRS virtual register + if register in ['TPWMTHRS', 'TCOOLTHRS']: + existing_thrs = self._registers.get('THRS', {}) + new_values = self._parse_register_values(cleaned_values) + merged_values = {**existing_thrs, **new_values} + self._registers['THRS'] = merged_values + else: + self._registers[register] = self._parse_register_values(cleaned_values) + + def _parse_register_values(self, register_string: str) -> Dict[str, Any]: + parsed = {} + parts = register_string.split() + for part in parts: + if '=' in part: + k, v = part.split('=', 1) + parsed[k] = v + return parsed + + def get_register(self, register: str) -> Optional[Dict[str, Any]]: + return self._registers.get(register) + + def get_registers(self) -> Dict[str, Dict[str, Any]]: + return self._registers + + def set_property(self, property: str, value: Any) -> None: + self._properties[property] = value + + def get_property(self, property: str) -> Optional[Any]: + return self._properties.get(property) + + def __str__(self): + return f'Stepper: {self._name}\nKlipper config: {self._properties}\nTMC Registers: {self._registers}' + + # Return the other motor properties and registers that are different from the current motor + def compare_to(self, other: 'Motor') -> Optional[Dict[str, Dict[str, Any]]]: + differences = {'properties': {}, 'registers': {}} + + # Compare properties + all_keys = self._properties.keys() | other._properties.keys() + for key in all_keys: + val1 = self._properties.get(key) + val2 = other._properties.get(key) + if val1 != val2: + differences['properties'][key] = val2 + + # Compare registers + all_keys = self._registers.keys() | other._registers.keys() + for key in all_keys: + reg1 = self._registers.get(key, {}) + reg2 = other._registers.get(key, {}) + if reg1 != reg2: + reg_diffs = {} + sub_keys = reg1.keys() | reg2.keys() + for sub_key in sub_keys: + reg_val1 = reg1.get(sub_key) + reg_val2 = reg2.get(sub_key) + if reg_val1 != reg_val2: + reg_diffs[sub_key] = reg_val2 + if reg_diffs: + differences['registers'][key] = reg_diffs + + # Clean up: remove empty sections if there are no differences + if not differences['properties']: + del differences['properties'] + if not differences['registers']: + del differences['registers'] + + if not differences: + return None + + return differences + + +class MotorLogParser: + _section_pattern: str = r'DUMP_TMC stepper_(x|y)' + _register_patterns: Dict[str, str] = { + 'CHOPCONF': r'CHOPCONF:\s+\S+\s+(.*)', + 'PWMCONF': r'PWMCONF:\s+\S+\s+(.*)', + 'COOLCONF': r'COOLCONF:\s+(.*)', + 'TPWMTHRS': r'TPWMTHRS:\s+\S+\s+(.*)', + 'TCOOLTHRS': r'TCOOLTHRS:\s+\S+\s+(.*)', + } + + def __init__(self, filepath: Path, config_string: Optional[str] = None): + self._filepath = filepath + + self._motors: List[Motor] = [] + self._config = self._parse_config(config_string) if config_string else {} + + self._parse_registers() + + def _parse_config(self, config_string: str) -> Dict[str, Any]: + config = {} + entries = config_string.split('|') + for entry in entries: + if entry: + key, value = entry.split(':') + config[key.strip()] = self._convert_value(value.strip()) + return config + + def _convert_value(self, value: str) -> Union[int, float, bool, str]: + if value.isdigit(): + return int(value) + try: + return float(value) + except ValueError: + if value.lower() in ['true', 'false']: + return value.lower() == 'true' + return value + + def _parse_registers(self) -> None: + with open(self._filepath, 'r') as file: + log_content = file.read() + + sections = re.split(self._section_pattern, log_content) + + # Detect only the latest dumps from the log (to ignore potential previous and outdated dumps) + last_sections: Dict[str, int] = {} + for i in range(1, len(sections), 2): + stepper_name = 'stepper_' + sections[i].strip() + last_sections[stepper_name] = i + + for stepper_name, index in last_sections.items(): + content = sections[index + 1] + motor = Motor(stepper_name) + + # Apply general properties from config string + for key, value in self._config.items(): + if stepper_name in key: + prop_key = key.replace(stepper_name + '_', '') + motor.set_property(prop_key, value) + elif 'autotune' in key: + motor.set_property(key, value) + + # Parse TMC registers + for key, pattern in self._register_patterns.items(): + match = re.search(pattern, content) + if match: + values = match.group(1).strip() + motor.set_register(key, values) + + self._motors.append(motor) + + # Find and return the motor by its name + def get_motor(self, motor_name: str) -> Optional[Motor]: + for motor in self._motors: + if motor._name == motor_name: + return motor + return None + + # Get all the motor list at once + def get_motors(self) -> List[Motor]: + return self._motors + + +# # Usage example: +# config_string = "stepper_x_tmc:tmc2240|stepper_x_run_current:0.9|stepper_x_hold_current:0.9|stepper_y_tmc:tmc2240|stepper_y_run_current:0.9|stepper_y_hold_current:0.9|autotune_enabled:True|stepper_x_motor:ldo-35sth48-1684ah|stepper_x_voltage:|stepper_y_motor:ldo-35sth48-1684ah|stepper_y_voltage:|" +# parser = MotorLogParser('/path/to/your/logfile.log', config_string) + +# stepper_x = parser.get_motor('stepper_x') +# stepper_y = parser.get_motor('stepper_y') + +# print(stepper_x) +# print(stepper_y) diff --git a/src/is_workflow.py b/src/is_workflow.py new file mode 100755 index 0000000..5847fda --- /dev/null +++ b/src/is_workflow.py @@ -0,0 +1,424 @@ +#!/usr/bin/env python3 + +############################################ +###### INPUT SHAPER KLIPPAIN WORKFLOW ###### +############################################ +# Written by Frix_x#0161 # + +# This script is designed to be used with gcode_shell_commands directly from Klipper +# Use the provided Shake&Tune macros instead! + + +import abc +import argparse +import tarfile +import traceback +from datetime import datetime +from pathlib import Path +from typing import Callable, Optional + +from git import GitCommandError, Repo +from matplotlib.figure import Figure + +import src.helpers.filemanager as fm +from src.graph_creators.analyze_axesmap import axesmap_calibration +from src.graph_creators.graph_belts import belts_calibration +from src.graph_creators.graph_shaper import shaper_calibration +from src.graph_creators.graph_vibrations import vibrations_profile +from src.helpers.locale_utils import print_with_c_locale +from src.helpers.motorlogparser import MotorLogParser + + +class Config: + KLIPPER_FOLDER = Path.home() / 'klipper' + KLIPPER_LOG_FOLDER = Path.home() / 'printer_data/logs' + RESULTS_BASE_FOLDER = Path.home() / 'printer_data/config/K-ShakeTune_results' + RESULTS_SUBFOLDERS = {'belts': 'belts', 'shaper': 'inputshaper', 'vibrations': 'vibrations'} + + @staticmethod + def get_results_folder(type: str) -> Path: + return Config.RESULTS_BASE_FOLDER / Config.RESULTS_SUBFOLDERS[type] + + @staticmethod + def get_git_version() -> str: + try: + # Get the absolute path of the script, resolving any symlinks + # Then get 1 times to parent dir to be at the git root folder + script_path = Path(__file__).resolve() + repo_path = script_path.parents[1] + repo = Repo(repo_path) + try: + version = repo.git.describe('--tags') + except GitCommandError: + version = repo.head.commit.hexsha[:7] # If no tag is found, use the simplified commit SHA instead + return version + except Exception as e: + print_with_c_locale(f'Warning: unable to retrieve Shake&Tune version number: {e}') + return 'unknown' + + @staticmethod + def parse_arguments() -> argparse.Namespace: + parser = argparse.ArgumentParser(description='Shake&Tune graphs generation script') + parser.add_argument( + '-t', + '--type', + dest='type', + choices=['belts', 'shaper', 'vibrations', 'axesmap'], + required=True, + help='Type of output graph to produce', + ) + parser.add_argument( + '--accel', + type=int, + default=None, + dest='accel_used', + help='Accelerometion used for vibrations profile creation or axes map calibration', + ) + parser.add_argument( + '--chip_name', + type=str, + default='adxl345', + dest='chip_name', + help='Accelerometer chip name used for vibrations profile creation or axes map calibration', + ) + parser.add_argument( + '--max_smoothing', + type=float, + default=None, + dest='max_smoothing', + help='Maximum smoothing to allow for input shaper filter recommendations', + ) + parser.add_argument( + '--scv', + '--square_corner_velocity', + type=float, + default=5.0, + dest='scv', + help='Square corner velocity used to compute max accel for input shapers filter recommendations', + ) + parser.add_argument( + '-m', + '--kinematics', + dest='kinematics', + default='cartesian', + choices=['cartesian', 'corexy'], + help='Machine kinematics configuration used for the vibrations profile creation', + ) + parser.add_argument( + '--metadata', + type=str, + default=None, + dest='metadata', + help='Motor configuration metadata printed on the vibrations profiles', + ) + parser.add_argument( + '-c', + '--keep_csv', + action='store_true', + default=False, + dest='keep_csv', + help='Whether to keep the raw CSV files after processing in addition to the PNG graphs', + ) + parser.add_argument( + '-n', + '--keep_results', + type=int, + default=3, + dest='keep_results', + help='Number of results to keep in the result folder after each run of the script', + ) + parser.add_argument('--dpi', type=int, default=150, dest='dpi', help='DPI of the output PNG files') + parser.add_argument('-v', '--version', action='version', version=f'Shake&Tune {Config.get_git_version()}') + + return parser.parse_args() + + +class GraphCreator(abc.ABC): + def __init__(self, keep_csv: bool, dpi: int): + self._keep_csv = keep_csv + self._dpi = dpi + + self._graph_date = datetime.now().strftime('%Y%m%d_%H%M%S') + self._version = Config.get_git_version() + + self._type = None + self._folder = None + + def _setup_folder(self, graph_type: str) -> None: + self._type = graph_type + self._folder = Config.get_results_folder(graph_type) + + def _move_and_prepare_files( + self, + glob_pattern: str, + min_files_required: Optional[int] = None, + custom_name_func: Optional[Callable[[Path], str]] = None, + ) -> list[Path]: + tmp_path = Path('/tmp') + globbed_files = list(tmp_path.glob(glob_pattern)) + + # If min_files_required is not set, use the number of globbed files as the minimum + min_files_required = min_files_required or len(globbed_files) + + if not globbed_files: + raise FileNotFoundError(f'no CSV files found in the /tmp folder to create the {self._type} graphs!') + if len(globbed_files) < min_files_required: + raise FileNotFoundError(f'{min_files_required} CSV files are needed to create the {self._type} graphs!') + + lognames = [] + for filename in sorted(globbed_files, key=lambda f: f.stat().st_mtime, reverse=True)[:min_files_required]: + fm.wait_file_ready(filename) + custom_name = custom_name_func(filename) if custom_name_func else filename.name + new_file = self._folder / f'{self._type}_{self._graph_date}_{custom_name}.csv' + filename.rename(new_file) + fm.wait_file_ready(new_file) + lognames.append(new_file) + return lognames + + def _save_figure_and_cleanup(self, fig: Figure, lognames: list[Path], axis_label: Optional[str] = None) -> None: + axis_suffix = f'_{axis_label}' if axis_label else '' + png_filename = self._folder / f'{self._type}_{self._graph_date}{axis_suffix}.png' + fig.savefig(png_filename, dpi=self._dpi) + + if self._keep_csv: + self._archive_files(lognames) + else: + self._remove_files(lognames) + + def _archive_files(self, _: list[Path]) -> None: + return + + def _remove_files(self, lognames: list[Path]) -> None: + for csv in lognames: + csv.unlink(missing_ok=True) + + @abc.abstractmethod + def create_graph(self) -> None: + pass + + @abc.abstractmethod + def clean_old_files(self, keep_results: int) -> None: + pass + + +class BeltsGraphCreator(GraphCreator): + def __init__(self, keep_csv: bool = False, dpi: int = 150): + super().__init__(keep_csv, dpi) + + self._setup_folder('belts') + + def create_graph(self) -> None: + lognames = self._move_and_prepare_files( + glob_pattern='raw_data_axis*.csv', + min_files_required=2, + custom_name_func=lambda f: f.stem.split('_')[3].upper(), + ) + fig = belts_calibration( + lognames=[str(path) for path in lognames], + klipperdir=str(Config.KLIPPER_FOLDER), + st_version=self._version, + ) + self._save_figure_and_cleanup(fig, lognames) + + def clean_old_files(self, keep_results: int = 3) -> None: + # Get all PNG files in the directory as a list of Path objects + files = sorted(self._folder.glob('*.png'), key=lambda f: f.stat().st_mtime, reverse=True) + + if len(files) <= keep_results: + return # No need to delete any files + + # Delete the older files + for old_file in files[keep_results:]: + file_date = '_'.join(old_file.stem.split('_')[1:3]) + for suffix in ['A', 'B']: + csv_file = self._folder / f'belts_{file_date}_{suffix}.csv' + csv_file.unlink(missing_ok=True) + old_file.unlink() + + +class ShaperGraphCreator(GraphCreator): + def __init__(self, keep_csv: bool = False, dpi: int = 150): + super().__init__(keep_csv, dpi) + + self._max_smoothing = None + self._scv = None + + self._setup_folder('shaper') + + def configure(self, scv: float, max_smoothing: float = None) -> None: + self._scv = scv + self._max_smoothing = max_smoothing + + def create_graph(self) -> None: + if not self._scv: + raise ValueError('scv must be set to create the input shaper graph!') + + lognames = self._move_and_prepare_files( + glob_pattern='raw_data*.csv', + min_files_required=1, + custom_name_func=lambda f: f.stem.split('_')[3].upper(), + ) + fig = shaper_calibration( + lognames=[str(path) for path in lognames], + klipperdir=str(Config.KLIPPER_FOLDER), + max_smoothing=self._max_smoothing, + scv=self._scv, + st_version=self._version, + ) + self._save_figure_and_cleanup(fig, lognames, lognames[0].stem.split('_')[-1]) + + def clean_old_files(self, keep_results: int = 3) -> None: + # Get all PNG files in the directory as a list of Path objects + files = sorted(self._folder.glob('*.png'), key=lambda f: f.stat().st_mtime, reverse=True) + + if len(files) <= 2 * keep_results: + return # No need to delete any files + + # Delete the older files + for old_file in files[2 * keep_results :]: + csv_file = old_file.with_suffix('.csv') + csv_file.unlink(missing_ok=True) + old_file.unlink() + + +class VibrationsGraphCreator(GraphCreator): + def __init__(self, keep_csv: bool = False, dpi: int = 150): + super().__init__(keep_csv, dpi) + + self._kinematics = None + self._accel = None + self._chip_name = None + self._motors = None + + self._setup_folder('vibrations') + + def configure(self, kinematics: str, accel: float, chip_name: str, metadata: str) -> None: + self._kinematics = kinematics + self._accel = accel + self._chip_name = chip_name + + parser = MotorLogParser(Config.KLIPPER_LOG_FOLDER / 'klippy.log', metadata) + self._motors = parser.get_motors() + + def _archive_files(self, lognames: list[Path]) -> None: + tar_path = self._folder / f'{self._type}_{self._graph_date}.tar.gz' + with tarfile.open(tar_path, 'w:gz') as tar: + for csv_file in lognames: + tar.add(csv_file, arcname=csv_file.name, recursive=False) + + def create_graph(self) -> None: + if not self._accel or not self._chip_name or not self._kinematics: + raise ValueError('accel, chip_name and kinematics must be set to create the vibrations profile graph!') + + lognames = self._move_and_prepare_files( + glob_pattern=f'{self._chip_name}-*.csv', + min_files_required=None, + custom_name_func=lambda f: f.name.replace(self._chip_name, self._type), + ) + fig = vibrations_profile( + lognames=[str(path) for path in lognames], + klipperdir=str(Config.KLIPPER_FOLDER), + kinematics=self._kinematics, + accel=self._accel, + st_version=self._version, + motors=self._motors, + ) + self._save_figure_and_cleanup(fig, lognames) + + def clean_old_files(self, keep_results: int = 3) -> None: + # Get all PNG files in the directory as a list of Path objects + files = sorted(self._folder.glob('*.png'), key=lambda f: f.stat().st_mtime, reverse=True) + + if len(files) <= keep_results: + return # No need to delete any files + + # Delete the older files + for old_file in files[keep_results:]: + old_file.unlink() + tar_file = old_file.with_suffix('.tar.gz') + tar_file.unlink(missing_ok=True) + + +class AxesMapFinder: + def __init__(self, accel: float, chip_name: str): + self._accel = accel + self._chip_name = chip_name + + self._graph_date = datetime.now().strftime('%Y%m%d_%H%M%S') + + self._type = 'axesmap' + self._folder = Config.RESULTS_BASE_FOLDER + + def find_axesmap(self) -> None: + tmp_folder = Path('/tmp') + globbed_files = list(tmp_folder.glob(f'{self._chip_name}-*.csv')) + + if not globbed_files: + raise FileNotFoundError('no CSV files found in the /tmp folder to find the axes map!') + + # Find the CSV files with the latest timestamp and wait for it to be released by Klipper + logname = sorted(globbed_files, key=lambda f: f.stat().st_mtime, reverse=True)[0] + fm.wait_file_ready(logname) + + results = axesmap_calibration( + lognames=[str(logname)], + accel=self._accel, + ) + + result_filename = self._folder / f'{self._type}_{self._graph_date}.txt' + with result_filename.open('w') as f: + f.write(results) + + +def main(): + options = Config.parse_arguments() + fm.ensure_folders_exist( + folders=[Config.RESULTS_BASE_FOLDER / subfolder for subfolder in Config.RESULTS_SUBFOLDERS.values()] + ) + + print_with_c_locale(f'Shake&Tune version: {Config.get_git_version()}') + + graph_creators = { + 'belts': (BeltsGraphCreator, None), + 'shaper': (ShaperGraphCreator, lambda gc: gc.configure(options.scv, options.max_smoothing)), + 'vibrations': ( + VibrationsGraphCreator, + lambda gc: gc.configure(options.kinematics, options.accel_used, options.chip_name, options.metadata), + ), + 'axesmap': (AxesMapFinder, None), + } + + creator_info = graph_creators.get(options.type) + if not creator_info: + print_with_c_locale('Error: invalid graph type specified!') + return + + # Instantiate the graph creator + graph_creator_class, configure_func = creator_info + graph_creator = graph_creator_class(options.keep_csv, options.dpi) + + # Configure it if needed + if configure_func: + configure_func(graph_creator) + + # And then run it + try: + graph_creator.create_graph() + except FileNotFoundError as e: + print_with_c_locale(f'FileNotFound error: {e}') + return + except TimeoutError as e: + print_with_c_locale(f'Timeout error: {e}') + return + except Exception as e: + print_with_c_locale(f'Error while generating the graphs: {e}') + traceback.print_exc() + return + + print_with_c_locale(f'{options.type} graphs created successfully!') + graph_creator.clean_old_files(options.keep_results) + print_with_c_locale(f'Cleaned output folder to keep only the last {options.keep_results} results!') + + +if __name__ == '__main__': + main()