Flexible exoskeleton control method and system with multi-source intent mapping and comfort constraints

Abstract

The application discloses a flexible exoskeleton control method and system of multi-source intention mapping and comfort constraint, relates to the technical field of automatic control, and aims at solving the problem of difficult compensation of existing Bowden cable exoskeleton hysteresis. The method comprises the following steps: collecting joint motion and electromyography and myoelectricity dual-mode physiological signals and preprocessing and extracting features; obtaining human body intention torque through physiological feature mapping, time alignment and safety check; combining an extended state observer to correct model parameters and weighting and fusing multi-source torques; adaptively updating impedance parameters with comfort constraint based on motion working conditions, synthesizing torque instructions through a PI closed loop and converting the torque instructions into joint target displacement; relying on a rate-related asymmetric RA-PI hysteresis model to identify parameters on-line and construct an inverse compensator, realizing displacement feedforward compensation and closed loop correction, and completing high-precision cascade closed loop control of the exoskeleton. The application effectively suppresses cable hysteresis disturbance, adapts to multiple motion working conditions, improves wearing comfort and motion following precision, and has strong robustness.

Description

Technical Field

[0001] This invention relates to the field of automatic control technology, and discloses a flexible exoskeleton control method and system based on multi-source intent mapping and comfort constraints. Background Technology

[0002] Soft exosuits are wearable assistive devices primarily constructed from fabrics, elastic materials, and cables. Compared to traditional rigid exosuits, they offer advantages such as lighter weight, greater freedom of movement, and natural compatibility with human kinesiology, showing broad application prospects in rehabilitation training, walking assistance, and reducing the burden of industrial labor. Soft exosuits typically use Bowden cables as the force transmission medium, transmitting the driving force generated by a distal motor through steel cables within a flexible sheath to the joint ends. However, Bowden cable transmission mechanisms inherently exhibit nonlinear hysteresis: due to Coulomb and viscous friction between the cable and the inner wall of the sheath, the displacement at the input and output ends of the cable exhibits a nonlinear mapping relationship related to the motion history path, manifesting as a hysteresis loop on the displacement-displacement characteristic curve.

[0003] Existing methods for modeling Bowden line hysteresis mainly include the Bouc-Wen model, the Preisach model, and the Prandtl-Ishlinskii (PI) model. The Bouc-Wen model describes hysteresis in the form of differential equations, but parameter identification relies on numerical optimization and lacks analytical invertibility. The Preisach model models hysteresis using operator superposition, offering high accuracy but high computational complexity, and its inverse model requires numerical iteration. The classic PI model, based on weighted superposition of the Play operator, has significant advantages in terms of structural simplicity and analytical invertibility. However, its original form has two inherent limitations: the hysteresis characteristics of the loading and unloading paths are constrained to a symmetrical form, while the actual Bowden line exhibits significant directional asymmetry in curved sheaths due to differences in curvature and preload; the model output depends only on the historical path of the input displacement and is independent of the motion rate, while the viscous friction component increases in the actual Bowden line at high speeds, and the hysteresis loop width expands significantly with speed. These two limitations result in insufficient modeling accuracy of the classic PI model under dynamic conditions of flexible exoskeletons, leading to large hysteresis compensation residuals.

[0004] Accurate estimation of human-machine interaction torque is a necessary feedback for closed-loop force / impedance control. Rigid exoskeletons usually install force / torque sensors at the joints to directly measure the interaction force, but the compliant structure and distributed force transmission characteristics of flexible exoskeletons make the installation of rigid force sensors face the following difficulties: (1) The rigid structure of the sensor is incompatible with the flexible fabric substrate, which destroys the wearing comfort and freedom of movement of the flexible exoskeleton; (2) The force transmission path on the flexible structure is uncertain, and the local force at the sensor installation position cannot represent the overall interaction torque; (3) Creep and relaxation of flexible materials cause severe zero-point drift of the force sensor. Summary of the Invention

[0005] This invention provides a flexible exoskeleton control method and system based on multi-source intent mapping and comfort constraints to solve the problems of difficult hysteresis compensation and poor torque tracking performance of existing Bowden line exoskeletons.

[0006] This invention employs the following technical solution: a flexible exoskeleton control method based on multi-source intent mapping and comfort constraints, comprising the following steps: S1. Simultaneously collect joint motion signals, surface electromyography signals, and muscle sound signals from flexible exoskeleton wearers, perform preprocessing for noise reduction and feature extraction respectively, and obtain dual-modal physiological feature signals and joint angle motion information; S2. Based on the mapping of initial motion intention torque using dual-modal physiological characteristic signals, combined with electromechanical delay completion time alignment, the output intention torque with safety constraints is verified. S3. Establish an exoskeleton joint impedance dynamics model, observe the joint motion state and human-machine interaction torque through a nonlinear expansion state observer, obtain the estimated value of the interaction torque, correct the model parameters online, and obtain the fused torque feedback signal by weighted fusion of multi-source torque information. S4. Construct a target impedance model with wearable comfort quantitative constraints, and adaptively and smoothly update the four impedance parameters, namely virtual inertia, virtual damping, virtual stiffness and impedance mixing weight, according to real-time motion conditions. Perform comfort constraint projection and boundary limiting correction on the impedance parameters. Based on the intention torque and fused feedback torque, synthesize the torque reference command through PI closed-loop correction, and convert it into joint target displacement command after gain stability verification. S5. Based on the rate-dependent asymmetric RA-PI hysteresis model, online identification of hysteresis parameters is performed. Based on the identified hysteresis parameters, an inverse hysteresis compensator is constructed to perform feedforward hysteresis compensation on the target displacement command. After displacement closed-loop correction and amplitude limiting, the motor control signal is output to drive the exoskeleton to complete impedance-torque cascade closed-loop control.

[0007] Further, S1 specifically includes collecting surface electromyography (sEMG) signals and myosalpingography (MMG) signals from the target muscle group in the human body; The electromyography (EMG) signal was processed by baseline removal, fourth-order Butterworth bandpass filtering, and power frequency notch filtering to obtain a noise-reduced effective EMG signal. The muscle sound signal is first synthesized into triaxial component magnitudes, and then denoised using a fourth-order Butterworth bandpass filter to obtain the denoised effective muscle sound signal. Extracting bimodal physiological feature signals from high signal-to-noise ratio electromyography (EMG) and muscle sound (SMS) effective signals; Based on the preprocessed dual-modal effective signal, a fixed window length is used. Fixed window step The sliding window mechanism extracts signal data within the window and extracts multi-dimensional signal features. The sEMG time-domain features include root mean square value, waveform length, and zero-crossing rate, calculated as follows: ; In the formula, For the first The root mean square eigenvalue of the sEMG signal under a single muscle channel. For window length, For the time index of discrete sampling points within the window, The starting index of the current window; ; In the formula, For the first Waveform length of sEMG signal under a single muscle channel. For the first Channel, First The noise-reduced effective electromyographic signal at each sampling point; ; In the formula, For the first Zero-crossing rate of sEMG signal under a single muscle channel. For indicator functions, The threshold for determining zero crossing and ; The time-frequency domain characteristics of MMG include the root mean square value and the peak frequency, in order: ; In the formula, For the first The root mean square value of the MMG signal under a single muscle channel; ; In the formula, For the first The peak frequency of MMG signal within a single muscle channel. For discrete frequency points in the frequency domain, and These are the lower and upper frequency limits of the effective frequency band for the MMG signal, respectively. For the first Channel MMG signal at frequency Short-time power spectral density at; Two types of cross-modal coupling core features are constructed, namely sEMG-MMG electromechanical coupling coefficient and electromechanical delay; The sEMG-MMG electromechanical coupling coefficient for: ; In the formula, The first one in the current sliding window Effective electromyographic signals with noise reduction under muscle channels The mean, The first one in the current sliding window Effective signal for noise reduction muscle sounds under a single muscle channel The mean; The electromechanical delay estimation for: ; In the formula, For the normalized cross-correlation function, For the number of delayed samples, This represents the maximum upper limit for the number of delayed samples; All features from all channels are concatenated to construct a fixed-dimensional bimodal physiological feature vector. : ; The joint angle information refers to the joint angular displacement. Joint angular velocity and angular acceleration .

[0008] Furthermore, the initial motion intention torque mapped based on dual-modal physiological feature signals includes: mapping the dual-modal physiological feature vector... The Dimensional components Iteratively update the sliding window and calculate the mean and variance: ; In the formula, For the first The moving average of the dimensional feature components within the current window. Normalized forgetting factor and , For the first The moving average of the dimensional feature components over the previous window. For the first The original feature values ​​of the dimensional feature components in the current window; ; In the formula, For the first The sliding variance of the 3D feature components in the current window. For the first The sliding variance of the feature components in the previous window; Then, adaptive Z-score normalization is performed to obtain the normalized feature vector. ; Calculate the first based on the electromechanical coupling coefficient. sEMG trust level of the channel : ; In the formula, It is the hyperbolic tangent function. For the first Muscle channels, current sliding window The electromechanical coupling coefficient, The electromechanical coupling degradation threshold, The sigmoid sharpness factor and , ; Calculate the first based on the electromechanical coupling coefficient. MMG Trust Level of the Channel : ; Construct the first Channel weighted eigenvectors : ; In the formula, For the first sEMG confidence level of individual muscle channels , , For the first Three sEMG normalized features of each muscle channel For the first MMG trust level for each muscle channel , For the first Two MMG-normalized features of each muscle channel, For the first The sixth normalized feature of the muscle channel For the first The 7th normalized feature of the muscle channel; Fusion feature vectors : ; In the formula, For the first The single-channel weighted feature vectors corresponding to each muscle channel, and , Fixing a single muscle channel includes 7 features. for Total feature dimensions after aggregating all channels; A Gaussian kernel ridge regression nonlinear mapping model was constructed, and offline training was performed using wearer-calibrated action data. The optimal model dual coefficients were obtained by solving for regularization constraints. Calculate the kernel vector Output the initial intended torque .

[0009] Furthermore, the intentional torque with safety constraints output after verification, combined with electromechanical delay completion time alignment, includes: According to the obtained Calculate the multi-channel weighted average electromechanical delay : ; right Apply first-order exponential smoothing: ; In the formula, This is the weighted average electromechanical delay value of the current sliding window after smoothing. The smoothing coefficient and , The interval between the previous and current windows is [value]. Weighted average electromechanical delay value after smoothing For window stepping; For joint torque Perform a time right shift operation to obtain the aligned joint torque values. : ; In the formula, The sampling period; right By applying safety saturation limits, the desired torque value of the safety constraint is obtained. : ; In the formula, The upper limit of the safe torque for the joint, This is the saturation limiting function.

[0010] Furthermore, the establishment of the exoskeleton joint impedance dynamics model, and the observation of joint motion state and human-machine interaction torque through a nonlinear expansion state observer to obtain estimated interaction torque values, includes: Input displacement via data acquisition cable With output displacement Differential encoder displacement acquisition and second-order Butterworth low-pass filtering are performed to obtain the output of the first... At each sampling time, the filtered cable input displacement , No. At each sampling time, the output displacement of the cable after filtering is... and the instantaneous motion rate of the cable input end ; A lumped-parameter dynamic model of a single joint of a flexible exoskeleton containing Bowden cable drive was established. The dynamic model was then converted into standard state-space equations to complete the joint dynamics modeling of the exoskeleton. A nonlinear extended state observer was used to monitor the system state and unknown disturbances in real time, and the observed joint angular velocities were obtained by solving the equations. Observation of human-computer interaction torque ; The online correction model parameters include: Calculate the root mean square value of the observed output error. : ; In the formula, For fixed parameters, For the first The observation output error at the sampling time, when Time-triggered parameter correction Adjust the threshold for the parameters; After trigger parameter correction, a joint angle acceleration regression equation is constructed. The Tikhonov regularized least squares method is used to solve the regression equation to obtain the parameter vector estimate. ; Based on parameter vector estimates The estimated equivalent moment of inertia was calculated. Equivalent elastic stiffness estimate Equivalent damping estimate ; By imposing physical constraints on the estimated equivalent moment of inertia, a new equivalent moment of inertia is obtained. : ; In the formula, For the new equivalent moment of inertia, and These are the minimum and maximum constraint values ​​for the equivalent moment of inertia, respectively. A fixed-weight soft update iterative method is used to smoothly update the equivalent rotational inertia, with a soft weighting factor as the weighting factor. Integrating historical equivalent rotational inertia before the update With the new equivalent moment of inertia after constraint treatment The final stable equivalent moment of inertia is obtained. The dynamic model parameters are adaptively corrected online, including soft update of the weight factors. The value is set to 0.8 to complete the iterative update of the equivalent moment of inertia.

[0011] Furthermore, the weighted fusion of multi-source torque information to obtain the fused torque feedback signal includes: Two torque information paths are set: Path A is the torque for active motion intention. Path B is the estimated value of the interaction torque. This enables independent output through dual channels. Time-scale alignment of the torque information of the two paths: ; In the formula, For the first At the sampling time, after time-scale alignment, the final active motion intention torque is obtained; Calculate the consistency residuals of the two torque paths and dynamic trust : ; ; In the formula, It is the minimum value. The adjustment coefficient of the sigmoid mapping and , For the first The consistency residual obtained after smoothing at the sampling time. The consistency residual threshold and ; Construct a multi-source torque complementary fusion formula: ; In the formula, This is the estimated torque value after fusion; When consistency residual Exceeding the residual safety threshold continuously And continuous sampling period threshold When this occurs, a safety degradation is triggered, locking the online parameter correction update operation and outputting only the fused torque estimate. .

[0012] Furthermore, S4 specifically includes: Construct a comfort constraint impedance model: ; In the formula, For the desired virtual inertia and , This is the second derivative of the joint angle deviation. For the desired virtual damping, The first derivative of the joint angle deviation. For desired virtual stiffness and , For joint angle difference and , This refers to the actual angle of the joint. For joint reference angle, For desired virtual stiffness and , The desired interaction torque; Combined with joint movement speed angular acceleration Temporal characteristics of electromyography Based on a preset angular velocity threshold Angular acceleration threshold Low threshold of electromyography High threshold of electromyography It distinguishes between four human movement conditions: support, extension, swing, and impact; and presets the priority of each condition and the most... Short dwell time constraints suppress frequent changes in operating conditions and output stable motion conditions. ; Based on steady motion conditions The virtual inertia is updated in real time using an exponential smoothing formula: ; In the formula, This is the real-time virtual inertia updated at the current moment. The smoothing coefficient is the parameter. , This is the virtual inertia from the previous moment. Preset values ​​for the baseline virtual inertia corresponding to each motion state; Virtual damping Virtual stiffness Impedance torque mixing ratio Synchronous updates using a smooth formula; Virtual damping Virtual stiffness Perform comfort constraint projection correction on virtual inertia Apply a lower bound constraint, Amplitude is limited to minimum comfortable virtual inertia Above, regarding virtual stiffness Apply upper bound constraints, The amplitude is limited to the maximum upper limit of virtual stiffness. Below; when At that time, check and correct the damping ratio parameter; Based on the actual joint position, reference trajectory, and observed angular velocity, the angular deviation, deviation angular velocity, and deviation angular acceleration are discretized and solved, and then substituted into the target impedance model to calculate the desired interaction torque. Calculate torque tracking error ; ; Resultant torque reference command: ; In the formula, This is the original unlimited torque reference command. For proportional gain, For integral gain, For the impedance torque mixing ratio and , For the first time in history Torque tracking error at sampling points; The final saturation limit is: ; In the formula, This serves as the final effective torque reference command. To achieve maximum safe output torque; by Actual angular displacement of exoskeleton joints As input, the torque command is converted into a cable target displacement command. : ; In the formula, This is the lever arm coefficient. This is the equivalent stiffness coefficient of the cable.

[0013] Furthermore, the online identification of hysteresis parameters based on the rate-dependent asymmetric RA-PI hysteresis model, the construction of a hysteresis inverse compensator based on the identified hysteresis parameters, the feedforward hysteresis compensation for the target displacement command, and the generation of motor control commands by combining displacement closed-loop correction, including: Configure multiple sets of fixed symmetric thresholds Play operator After recursive calculation, combined with linear gain With weight Weighted combination is used to construct a classic PI hysteresis model and output the predicted hysteresis displacement. ; Construct a rate-dependent asymmetric RA-PI hysteresis model adapted to the drive characteristics of flexible Bowden cables, abandoning the traditional symmetric single threshold model. The format uses independent forward loading thresholds. With reverse unloading threshold Composition of asymmetric threshold pairs , This corresponds to the differentiated hysteresis characteristics of the bidirectional motion of cable loading and tensioning, and unloading and release; by introducing a value range as... Asymmetric factors Achieve symmetric threshold To asymmetric threshold pairs Adaptive conversion based on asymmetric factors The different values ​​of the cable's bidirectional hysteresis dead zone width are dynamically adjusted to vary. Define cable input displacement increment Determine the real-time movement direction of the cable based on displacement increments. Positive and negative state matching asymmetric Play operator Differentiated output rules, where bidirectional thresholds are satisfied It degenerates into a traditional symmetrical hysteresis form; Introducing instantaneous motion rate at the cable input end Related rate envelope function The hysteresis loop expansion width is dynamically modulated at different motion speeds to characterize the differential characteristics of static hysteresis and high-speed dynamic hysteresis. The rate envelope function... Overall output of the RA-PI hysteresis model They are respectively: ; In the formula, The upper bound of the rate gain and , The rate sensitivity coefficient and ; ; In the formula, For linear gain, For the first The weights of the asymmetric Play operator Let be the rate envelope function. and These are the hysteresis thresholds for the forward loading and reverse unloading of the asymmetric Play operator, respectively. Model output This means outputting displacement at the far end of the Bowden line. The predicted value; The set of parameters to be identified in the RA-PI model Expanded to: ; The RA-PI hysteresis model is linearized, and a regression vector is constructed. With linear parameter vector Combined with the actual output displacement of the cable Constructing identification error Using a forgetting factor The recursive least squares algorithm updates the gain vector periodically. Linear parameter estimation With covariance matrix This enables high-frequency updates of linear parameters; For asymmetric factors Rate correlation coefficient , Every Data from a sliding window is extracted for each sampling period. The optimal nonlinear parameters are solved iteratively using the LM nonlinear least squares algorithm. By integrating the linear and nonlinear parameters, a complete set of RA-PI model identification parameters is obtained. ; A RA-PI inverse hysteresis feedforward compensator is constructed using a parameter analytical inversion method, based on the identified parameters at the current time. , , , , Define the parameters of the inverse model: ; ; ; ; In the formula, This represents the linear gain corresponding to the RA-PI inverse model. In the RA-PI inverse model, the first... The weights of the inverse asymmetric Play operator For the RA-PI model, the first The weights of the asymmetric Play operator For the RA-PI inverse model, the first The loading direction threshold of the asymmetric Play operator, i.e., the forward threshold of the RA-PI inverse model. For the RA-PI inverse model, the first The unloading direction threshold of an asymmetric Play operator, i.e., the negative threshold of the RA-PI inverse model; Based on the solved inverse model parameters, construct a system with inverse threshold pairs. Constrained inverse asymmetric Play hysteresis operator Based on the expected displacement increment The positive and negative states are distinguished by bidirectional operation rules, enabling differentiated inverse hysteresis solutions for cable loading and unloading directions; Solving for the hysteresis feedforward compensation input displacement by combining the results of the inverse asymmetric Play hysteresis operator operation. for: ; In the formula, To preset the desired cable output displacement, It is the inverse asymmetric Play hysteresis operator; By superimposing a feedback correction term on the feedforward compensation, the basic motor displacement command is obtained. : ; In the formula, For basic motor displacement commands For proportional gain, For integral gain, The residual tracking error between the expected output and the actual output; Will Replace with the system torque closed-loop output The final motor control command is obtained. : .

[0014] Furthermore, the output motor control signal after displacement closed-loop correction and amplitude limiting drives the exoskeleton to complete impedance-torque cascade closed-loop control, including: Based on the Routh stability criterion, the torque closed-loop gain constraint condition is set: ; ; In the formula, For torque proportional gain, For torque integral gain, For the equivalent moment of inertia, For equivalent damping, The cable damping coefficient; The safe torque proportional gain is obtained by correcting the torque proportional gain using the gain adaptive safety constraint formula. : ; In the formula, For gain stability margin and ; Based on the preset maximum cable movement speed With sampling period Calculate the displacement increment limiting threshold : ; Initial control commands for the motor are based on the displacement increment limiting threshold. Incremental clamping and limiting processing is performed to obtain the final system output motor control command. : ; The final control signal is output to drive the exoskeleton.

[0015] This invention also provides a flexible exoskeleton control system with multi-source intent mapping and comfort constraints, including: a signal acquisition module, an intent torque analysis module, a dynamic observation and torque fusion module, an adaptive impedance closed-loop module, and a hysteresis identification and compensation module; The signal acquisition module is used to simultaneously acquire joint motion signals, surface electromyography signals, and muscle sound signals. After preprocessing and feature extraction, it outputs dual-modal physiological feature signals and joint motion information. The intention torque analysis module is used to obtain the motion intention torque based on the dual-modal physiological characteristic signal, complete the electromechanical delay time alignment and safety verification, and output the intention torque with safety constraints. The dynamic observation and torque fusion module uses a nonlinear extended state observer to obtain the estimated value of the interactive torque, corrects the dynamic model parameters online, and weights and fuses the multi-source torques to obtain the fused torque feedback signal. The adaptive impedance closed-loop module is used to construct an impedance model with comfort constraints, adaptively update impedance parameters according to motion conditions, generate torque reference commands by combining the intended torque and the fused torque, and convert them into cable target displacement commands. The hysteresis identification and compensation module is used to identify hysteresis parameters based on the rate-dependent asymmetric RA-PI hysteresis model, construct an inverse hysteresis compensator to feed forward compensation for the target displacement command, and combine displacement closed-loop correction to output motor control signals to drive the flexible exoskeleton.

[0016] Compared with the prior art, the present invention has the following beneficial effects: This invention proposes a Bowden line hysteresis modeling and analytical inverse compensation method based on the rate-dependent asymmetric improved Prandtl-Ishlinskii model (RA-PI). By using direction-dependent asymmetric threshold pairs and rate-dependent exponential envelope functions, it overcomes the dual limitations of symmetry and rate independence of the classical PI model, while maintaining analytical invertibility to achieve real-time feedforward compensation.

[0017] This invention proposes a motion intent continuous decoding method based on the fusion of sEMG and MMG dual-modal signals. It introduces the electromechanical coupling coefficient (EMCC) as a cross-modal signal quality index to drive the dynamic adjustment of modal confidence, and uses cross-correlation analysis to estimate electromechanical delay in real time to achieve adaptive time alignment.

[0018] This invention proposes a torque estimation method for human-computer interaction using a nonlinear extended state observer (NLESO), establishes an online correction mechanism for nominal parameters triggered by observation errors, and constructs a dual-path complementary fusion and anomaly detection architecture of physical measurement chain and physiological intention chain.

[0019] This invention proposes a flexible exoskeleton control method based on multi-source intent mapping and comfort constraints. Wearing comfort is quantified into a constraint domain of impedance parameters. Based on motion state recognition, segmented impedance parameters are adaptively adjusted. A three-layer in-depth stability defense architecture is constructed through Routh criterion gain constraints, adaptive safety margin, and hardware-level displacement limiting. Attached Figure Description

[0020] Figure 1 This is a flowchart of the technology of the present invention. Detailed Implementation

[0021] The present invention will be further illustrated below with reference to embodiments. These embodiments are for illustrative purposes only and are not intended to limit the invention in any way. It should be understood that the described embodiments are merely some, not all, of the embodiments described in this application. All other embodiments obtained by those skilled in the art based on the embodiments in this application without inventive effort are within the scope of protection of this application.

[0022] Example 1 A flexible exoskeleton control method based on multi-source intent mapping and comfort constraints includes the following steps: S1. Simultaneously collect joint motion signals, surface electromyography signals, and muscle sound signals from flexible exoskeleton wearers, perform preprocessing for noise reduction and feature extraction respectively, and obtain dual-modal physiological feature signals and joint angle motion information; Surface electromyography (sEMG) and myophone (MMG) signals were collected from target muscle groups in the human body and simultaneously acquired and preprocessed. Differential surface electromyography (EMG) electrode pairs are placed on the surface of the target muscle belly along the muscle fiber direction. A miniature triaxial accelerometer is placed adjacent to the electrode pairs in the same muscle. All channels share the same sampling clock and sampling frequency. Sampling period ; set up Total number of muscle channels (in this embodiment) ), No. The raw sEMG signals of each channel are denoted as follows: The original MMG signal is denoted as , For discrete time step index; For sEMG channel preprocessing, the sEMG signal contains power frequency interference (50Hz), motion artifacts (<20Hz), and high-frequency instrument noise. The first step is to perform mean removal processing: ; In the formula For the first After the first sEMG channel undergoes mean-removing preprocessing, the second... Electromyography signal output values ​​at each sampling time point For the first The raw surface electromyography signals of each sEMG channel, For the first The raw signal of each sEMG channel recently The moving average of each sampling point and Used to eliminate baseline drift ( (i.e., a 1-second window). Next, a fourth-order Butterworth bandpass filter is used for bandpass filtering, and the passband... The cascaded second-order Biquad form is as follows: ; In the formula, For the first After the first sEMG channel is filtered by a fourth-order Butterworth bandpass filter, the second... The output signal at each discrete sampling time; For processing sEMG signals, the first One bandpass filter unit; A second-order IIR notch filter is used. Eliminate 50Hz power frequency interference: ; In the formula, For notch filters, The notch frequency and , The sampling frequency for signal acquisition. To control the bandwidth of the notch filter and This embodiment =0.95, corresponding to a 3 dB bandwidth of approximately 5 Hz; The preprocessed sEMG signal was obtained after eliminating 50Hz power frequency interference. ; The composite amplitude of the three-axis raw MMG signal is obtained. The raw MMG signal of each channel: ; In the formula, For the first The first muscle channel, the first The raw MMG signal at discrete time steps For the first The first muscle channel, the first A discrete time step, MEMS accelerometer The raw MMG signal acquired from the direction. For the first The first muscle channel, the first A discrete time step, MEMS accelerometer The raw MMG signal acquired from the direction. For the first The first muscle channel, the first A discrete time step, MEMS accelerometer The raw MMG signal acquired from the direction; A fourth-order Butterworth bandpass filter was used to filter the original MMG signal. Perform bandpass filtering, wherein the passband of the bandpass filter is ; The fourth-order Butterworth bandpass filter is implemented by cascading two second-order sections, and the filtering formula is as follows: ; In the formula, For the first Each muscle channel, after being denoised using a fourth-order Butterworth bandpass filter, yields an effective denoised muscle sound signal. For the first A second-order section transfer function for bandpass filtering of MMG signals; Collect all sEMG signals after preprocessing of each channel and MMG signal This forms a signal set: ; In the formula, This represents the total number of muscle channels; Thus, dual-modal physiological characteristic signals are obtained; Use fixed window length (This embodiment takes) =200, i.e., 200ms), fixed window step. (This embodiment takes) =50, which is 50ms, corresponding to an intention update frequency of 20Hz. A sliding window is used to extract timing data frame by frame, with the starting index of the current window set to . The data segment within the window is denoted as and We extract sEMG time-domain features and MMG time-frequency features respectively, and construct cross-modal coupling features and electromechanical delay features. The sEMG time-domain features include the root mean square value, waveform length, and zero-crossing rate, calculated as follows: ; In the formula, For the first The root mean square eigenvalue of the sEMG signal under a single muscle channel. For window length, For the time index of discrete sampling points within the window, The starting index of the current window; ; In the formula, For the first Waveform length of sEMG signal under a single muscle channel. For the first Channel, First The noise-reduced effective electromyographic signal at each sampling point; ; In the formula, For the first Zero-crossing rate of sEMG signal under muscle channels For indicator functions, The threshold for determining zero crossing and ; The MMG time-frequency domain characteristics include the root mean square value and the peak frequency, in order: ; In the formula, For the first The root mean square value of the MMG signal under a single muscle channel; ; In the formula, For the first The peak frequency of MMG signal within a single muscle channel. For discrete frequency points in the frequency domain, and These are the lower and upper frequency limits of the effective frequency band for the MMG signal, respectively. For the first Channel MMG signal at frequency Short-time power spectral density at; Two types of cross-modal coupling core features are constructed, namely sEMG-MMG electromechanical coupling coefficient and electromechanical delay; The electromechanical coupling coefficient of sEMG-MMG for: ; In the formula, The first one in the current sliding window Effective electromyographic signals with noise reduction under muscle channels The mean, The first one in the current sliding window Effective signal for noise reduction muscle sounds under a single muscle channel The mean; The electromechanical delay estimation for: ; In the formula, For the normalized cross-correlation function, For the number of delayed samples, This is the maximum upper limit of the number of delayed samples (in this embodiment) (i.e., 80ms); All All features from each channel are concatenated to construct a fixed-dimensional bimodal physiological feature vector. : ; In this embodiment Therefore ; The joint angle information refers to the joint angular displacement. Joint angular velocity and angular acceleration .

[0023] S2. Based on the mapping of initial motion intention torque using dual-modal physiological characteristic signals, combined with electromechanical delay completion time alignment, the output intention torque with safety constraints is verified. For bimodal physiological feature vectors The Dimensional components Iterative updates of the moving average mean and variance: ; In the formula, For the first The moving average of the dimensional feature components within the current window. Normalized forgetting factor and (This embodiment) ), For the first The moving average of the dimensional feature components over the previous window. For the first The original feature values ​​of the dimensional feature components in the current window; ; In the formula, For the first The sliding variance of the 3D feature components in the current window. For the first The sliding variance of the feature components in the previous window; For characteristic components Perform adaptive Z-score normalization to obtain the normalized feature vector. : ; In the formula, To prevent small constants from being divided by zero and Normalized features ; Calculate the first sEMG trust level of the channel : ; In the formula, It is the hyperbolic tangent function. For the first One muscle channel, current sliding window The electromechanical coupling coefficient, Electromechanical coupling degradation threshold (in this embodiment) ), The sigmoid sharpness factor and (In this embodiment) ), ; when hour, Tendency towards 1, complete trust in sEMG; when hour, The value tends to 0, so sEMG is unreliable and relies on MMG. Calculate the first MMG Trust Level of the Channel : ; Construct the first Channel weighted eigenvectors : ; In the formula, For the first sEMG confidence level of individual muscle channels , , For the first Three sEMG normalized features of each muscle channel For the first MMG trust level for each muscle channel , For the first Two MMG-normalized features of each muscle channel, For the first The sixth normalized feature of the muscle channel For the first The 7th normalized feature of the muscle channel; Fusion feature vectors : ; In the formula, For the first The single-channel weighted feature vectors corresponding to each muscle channel, and , Fixing a single muscle channel includes 7 features. for Total feature dimensions after aggregating all channels; Execution during the calibration phase A predefined action, with reference torque recorded via an external torque sensor. ; Collect training set and construct matrix sum vector : ; In the formula, Corresponding to the 1st, 2nd to 3rd in sequence The fused feature vector of each training window; ; In the formula, Corresponding to the 1st, 2nd to 3rd in sequence Reference torque values ​​for each training window; Constructing the kernel matrix Element shift Gaussian kernel function: ; In the formula, The first in the kernel matrix line, number The element values ​​of the column, For kernel bandwidth parameters and ; Determined by the median of the calibrated data: ; Calculate the dual coefficient vector of kernel ridge regression : ; In the formula, The dual coefficient vector obtained after training with kernel ridge regression and , The regularization parameter and (This embodiment) ); Calculate the fusion features of the current window kernel vectors of all training samples : ; In the formula, The fused feature vector corresponding to the current window The fused feature vector with the first training sample Gaussian kernel similarity between them The fused feature vector corresponding to the current window The fused feature vector with the first training sample Gaussian kernel similarity between them; Output an estimate of the initial intended torque. : ; In the formula, For continuous intention torque output based on sEMG-MMG dual-modal fusion; According to the obtained Calculate the multi-channel weighted average electromechanical delay : ; right Apply first-order exponential smoothing: ; In the formula, This is the weighted average electromechanical delay value of the current sliding window after smoothing. Smoothing coefficient and , The interval between the previous and current windows is [value]. Weighted average electromechanical delay value after smoothing For window stepping; right Perform a time right shift operation to obtain the aligned joint torque values. : ; In the formula, The sampling period; right By applying safety saturation limits, the desired torque value of the safety constraint is obtained. : ; ; In the formula, The upper limit of the safe torque for the joint (ankle joint in this embodiment) ), This is the saturation limiting function.

[0024] S3. Establish an exoskeleton joint impedance dynamics model, observe the joint motion state and human-machine interaction torque through a nonlinear expansion state observer, obtain the estimated value of the interaction torque, correct the model parameters online, and obtain the fused torque feedback signal by weighted fusion of multi-source torque information. Data acquisition of Bowden cable input end displacement Output displacement A second-order Butterworth low-pass filter is used for... and The displacement signal after smoothing and filtering and : ; ; In the formula, For unit delay operators, The coefficients are discrete transfer function coefficients, determined by the cutoff frequency. and sampling frequency Determined through bilinear transformation, The selection principle is as follows: Higher than the maximum operating frequency of the Bowden line Below the main noise frequency band, (In this embodiment) ); Calculate the first-order differential rate of the input displacement: ; In the formula, The instantaneous speed of movement at the cable input end. For the first The cable input displacement after filtering at each sampling time. For the first The cable input displacement after filtering at each sampling time. For discrete time step index; Establish the rotational motion dynamics equations from the motor side to the joint side: ; In the formula, The equivalent rotational inertia of the joint and , Joint angular acceleration, The driving torque applied by the Bowden line at the joint ends, This is the equivalent viscous damping coefficient. The joint angular velocity, The equivalent elastic stiffness coefficient and , This refers to joint angular displacement. Human-computer interaction torque, This is the gravitational torque; gravitational torque Modeling: ; In the formula, For foot equivalent quality, It is the acceleration due to gravity. This is the distance from the joint axis to the center of mass of the foot. The angle of deviation of the joint angle from the vertical direction when standing naturally; Bowden line driving torque modeling: ; In the formula, This is the equivalent axial stiffness of the cable. This is the equivalent axial damping of the cable. For joint lever arm; Rewrite the dynamic equations in state-space form and define the state variables. Input quantity for: ; In the formula, For Bowden line driving torque, This is the gravitational torque; The state-space equations are: ; In the formula, The system state coefficient matrix, To control the input allocation matrix, To synthesize control input quantity, Assign a matrix to the perturbation input. For interactive torque; Define extended state vector : ; In the formula, , which is the normalized interaction torque acceleration component; The extended state-space equations are: ; ; ; In the formula, This is a time-varying external disturbance signal. This is the upper bound of the disturbance amplitude; The system is configured to measure joint angular displacement as its output. ; Construct a nonlinear extended state observer. The state update of the nonlinear extended state observer satisfies the following formula: ; ; ; In the formula, These are the observed values ​​of joint angular displacement. These are the observed values ​​of the joint angular velocity. For the normalized observations of the interaction torque, For the observation output error and , For nonlinear functions, when hour, ,when hour, , , , It is a non-linear exponent and , The threshold of the linear segment and (This embodiment) ), , , For observer gain, and , For the observer band and ; The discretization update formula for the observer, implemented using the discretization method, is as follows: ; ; ; In the formula, , , The first Extended state vector at each sampling time , , The corresponding observed values, This represents the observation output error at the current moment; Physical interaction torque estimate : ; Online correction of nominal parameters and adaptive compensation of observer based on angular velocity-angular acceleration regression; Calculate the root mean square value of the observed output error. : ; In the formula, For fixed parameters, For the first The observation output error at the sampling time, when Time-triggered parameter correction Adjust the threshold for the parameters; The relationship between joint angular acceleration and relevant parameters is established, and the regression equation is as follows: ; In the formula, is a regression vector and , are parameter vectors and ; The regression equation is solved using the Tikhonov regularized least squares method to obtain the parameter vector estimates. : ; In the formula, For the regression matrix, This is a regularization adjustment parameter. It is a 3-order identity matrix. It is the angular acceleration vector; Perform parameter recovery and update based on the parameter vector estimates. Calculate the estimated equivalent moment of inertia. Equivalent elastic stiffness estimate Equivalent damping estimate : ; ; ; By imposing physical constraints on the equivalent moment of inertia, a new equivalent moment of inertia is obtained. : ; In the formula, For the new equivalent moment of inertia, and These are the minimum and maximum constraint values ​​for the equivalent moment of inertia, respectively. The equivalent moment of inertia is updated using a soft update method to obtain the final updated equivalent moment of inertia. : ; In the formula, For soft update weight factors and , The equivalent moment of inertia before the update. This is the new equivalent moment of inertia obtained through physical constraints; The weighted fusion of multi-source torque information to obtain the fused torque feedback signal includes: Two torque information paths are set: Path A is the torque for the intention of active motion. Path B is the estimated value of the interaction torque. This enables independent output through dual channels. Time-scale alignment of the torque information of the two paths: ; In the formula, For the first At the sampling time, after time-scale alignment, the final active motion intention torque is obtained; Construct a multi-source torque complementary fusion formula: ; In the formula, This is the estimated torque value after fusion. For NLESO trust level, This is an estimate of the physical interaction torque; Calculate the consistency residuals of the two torque paths and dynamic trust : ; ; In the formula, It is the minimum value. The adjustment coefficient of the sigmoid mapping and , For the first The consistency residual obtained after smoothing at the sampling time. The consistency residual threshold and ; Set residual safety threshold Sampling period threshold , This corresponds to a sampling duration of 0.5 seconds; When consistency residual Continuous exceedance value And continuous When this occurs, a safety degradation is triggered, locking the online parameter correction update operation and only outputting the fused torque estimate. .

[0025] S4. Construct a target impedance model with wearable comfort quantitative constraints, and adaptively and smoothly update impedance parameters such as virtual inertia, virtual damping, virtual stiffness, and impedance mixing weights according to real-time motion conditions. Perform comfort constraint projection and boundary limiting correction on the impedance parameters. Based on the intended torque and fused feedback torque, synthesize the torque reference command through PI closed-loop correction, and convert it into joint target displacement command after gain stability verification. Construct a comfort constraint impedance model: ; In the formula, For the desired virtual inertia and , This is the second derivative of the joint angle deviation. For the desired virtual damping, The first derivative of the joint angle deviation. For desired virtual stiffness and , For joint angle difference and , This refers to the actual angle of the joint. For joint reference angle, For desired virtual stiffness and , The desired interaction torque; At the discrete sampling time of the controller Discrete formulas are used to solve for the deviation and first and second-order differential components in real time: ; ; ; In the formula, for Constant joint angle deviation, for The actual angular displacement of the joint at any given moment. for Always close Section reference angular displacement, for The deviation angular velocity at any moment For the output of the nonlinear extended state observer Observed joint angular velocity values ​​at time points. for Constantly refer to the trajectory angular velocity, for Angular acceleration at any given moment For the output of the nonlinear extended state observer Observed joint angular velocity values ​​at time points. for Constantly reference trajectory angular acceleration; Combined with joint movement speed angular acceleration Temporal characteristics of electromyography Based on a preset angular velocity threshold Angular acceleration threshold Low threshold of electromyography High threshold of electromyography It distinguishes four types of human movement conditions: support, extension, swing, and impact. Determination of four motion conditions: support: ; push-up: ;swing: Impact: ; The priority order is: support > extension > swing > impact, and a minimum dwell time is introduced. (50 ms) Prevents frequent switching and outputs stable motion conditions. ; Based on steady motion conditions The virtual inertia is updated in real time using an exponential smoothing formula: ; In the formula, This is the real-time virtual inertia updated at the current moment. The smoothing coefficient is the parameter. , This is the virtual inertia from the previous moment. Preset values ​​for the baseline virtual inertia corresponding to each motion state; Virtual damping Virtual stiffness Impedance torque mixing ratio Synchronous updates using a smooth formula; Virtual damping Virtual stiffness Perform comfort constraint projection correction: ; ; In the formula, For minimum comfortable virtual inertia, For maximum true virtual stiffness; when At that time, check and correct the damping ratio parameter; Based on the actual joint position, reference trajectory, and observed angular velocity, the angular deviation, deviation angular velocity, and deviation angular acceleration are discretized and solved, and then substituted into the target impedance model to calculate the desired interaction torque. Calculate torque tracking error ; ; Resultant torque reference command: ; In the formula, This is the original unlimited torque reference command. For proportional gain, For integral gain, The impedance torque is mixed in proportion and , For the first time in history Torque tracking error at sampling points; The final saturation limit is: ; In the formula, This serves as the final effective torque reference command. To achieve maximum safe output torque; by Actual angular displacement of exoskeleton joints As input, the torque command is converted into a cable target displacement command. : ; In the formula, This is the lever arm coefficient. This is the equivalent stiffness coefficient of the cable.

[0026] S5. Based on the rate-dependent asymmetric RA-PI hysteresis model, online identification of hysteresis parameters is performed. Based on the identified hysteresis parameters, an inverse hysteresis compensator is constructed to perform feedforward hysteresis compensation on the target displacement command. Combined with displacement closed-loop correction, motor control command is generated. After displacement closed-loop correction and amplitude limiting, the motor control signal is output to drive the exoskeleton to complete impedance-torque cascade closed-loop control. Displacement input via cable Differential encoder displacement acquisition and second-order Butterworth low-pass filtering are performed to obtain the output of the first... At each sampling time, the filtered cable input displacement Configure multiple sets of symmetric thresholds for input. Play operator : ; In the formula, For the first Symmetric thresholds of the Play operators and , For the first The sampling time of the first sampling moment The output of the Play operator For the first At the sampling time, the first The output of the Play operator, initial value ; , The total number of Play operators; Threshold sequence In the interval The distribution is uniform, and the calculation formula is: ; In the formula, Based on hysteresis error The peak value determines the prior estimate of the maximum hysteresis width. , The number of Play operators is taken in this embodiment. ; Total output of the PI hysteresis model for: ; In the formula, This represents the predicted output displacement of the Bowden line by the PI hysteresis model. For linear gain and , For the first The weights of the Play operator and ; The set of parameters to be identified in the PI hysteresis model for: ; In the formula, For the first The weights of the Play operator, For the first The weights of each Play operator; Symmetric threshold Replace with direction-dependent asymmetric threshold pairs Define the asymmetric Play operator ; Define input displacement increment : ; when The output value of the asymmetric Play operator when the cable is moving in the loading / tension direction. for: ; In the formula, For the first The sampling time of the first sampling moment The output values ​​of the asymmetric Play operator For the first The loading direction threshold of the asymmetric Play operator and , For the first The unloading direction threshold of the asymmetric Play operator and Control the hysteresis dead zone width during release. Right now For the first The sampling time of the first sampling moment Output values ​​of asymmetric Play operators; when The output value of the asymmetric Play operator when the cable is moving in the unloading / release direction. for ; In the formula, For the first The unloading direction threshold of the asymmetric Play operator and , Control the width of the hysteresis dead zone during release; when It degenerates into a symmetrical form; Construction of asymmetric threshold pairs: Introducing an asymmetric factor and ,when This indicates that the dead zone in the loading direction is narrower than that in the unloading direction; This indicates that the dead zone in the loading direction is wider than that in the unloading direction. The loading direction threshold is calculated using the following formula. and unloading direction threshold : ; ; In the formula, For the first Symmetric thresholds for each Play operator; A rate-dependent envelope function is superimposed on the output of each asymmetric Play operator. : ; In the formula, The upper bound of the rate gain and The maximum expansion ratio of the hysteresis loop relative to the static width at high speed is controlled, which is physically the upper limit of the ratio of viscous friction to Coulomb friction. The rate sensitivity coefficient and , The larger the hysteresis loop, the faster it reaches saturation and expansion; when hour It manifests as static lag; when hour This manifests as maximum expansion; The total output of the RA-PI model, combining the asymmetric Play operator and the rate-dependent envelope function. : ; In the formula, The rate envelope function; Model output This means outputting displacement at the far end of the Bowden line. The predicted value; The set of parameters to be identified in the RA-PI model Expanded to: ; Rewrite the RA-PI model output in parametric linear form and define the regression vector. : ; In the formula, , for the first The velocity envelope value at each sampling time. Right now , for the first The output value of the asymmetric Play operator at the current time. for 3D real space; Define a linear parameter vector : ; The total output of the RA-PI model It can also be expressed as: ; The output after filtering Define the identification error as the target value. : ; In the formula, For the first linear parameter vector at each sampling time The estimated value; The parameter update rule for Recursive Least Squares-FF with a forgetting factor is as follows: ; ; ; In the formula, For the gain vector, , Let covariance matrix be the variance matrix. initial value , Large positive numbers (in this embodiment) ), For unit array; Forgetting factor, (In this embodiment) ), The algorithm gradually reduces the weight of old data, thereby tracking the slow drift of hysteresis parameters caused by sheath wear or path changes; For nonlinear parameters A hierarchical identification strategy is adopted: a Levenberg-Marquardt-based nonlinear least squares optimization is performed at a lower frequency in the outer layer (in this embodiment). That is, the nonlinear least squares optimization of LM: ; In the formula, For the first Identification estimates of asymmetric factors of the asymmetric Play operator For the parameters of the rate-dependent envelope function The identification estimate, For the parameters of the rate-dependent envelope function The identification estimate, The sampling time sequence number used for local summation. The number of sampling points in the sliding window for LM nonlinear identification. For the first The actual output displacement of the Bowden line after filtering at each sampling time point. For the first The theoretical predicted output of the RA-PI model at each sampling time point linear parameter vector The estimated value; Merging identified linear parameters With all nonlinear parameters The complete set of identification parameters of the RA-PI model were obtained by assembly. : ; A RA-PI inverse hysteresis feedforward compensator is constructed using a parameter analytical inversion method, based on the identified parameters at the current time. , , , , Define the parameters of the inverse model: ; ; ; ; In the formula, This represents the linear gain corresponding to the RA-PI inverse model. In the RA-PI inverse model, the first... The weights of the inverse asymmetric Play operator For the RA-PI model, the first The weights of the asymmetric Play operator For the RA-PI inverse model, the first The loading direction threshold of the asymmetric Play operator, i.e., the forward threshold of the RA-PI inverse model. For the RA-PI inverse model, the first The unloading direction threshold of an asymmetric Play operator, i.e., the negative threshold of the RA-PI inverse model; by As input, with inverse parameter Construct the inverse asymmetric Play operator with the threshold value. ; when At that time, the output value of the inverse asymmetric Play operator for: ; In the formula, For the first The sampling time of the first sampling moment The output values ​​of the inverse asymmetric Play operator For the first At each sampling time, given the desired output displacement at the far end of the cable, Right now For the first The sampling time of the first sampling moment The output values ​​of the inverse asymmetric Play operator; when At that time, the output value of the inverse asymmetric Play operator for: ; Compensation input displacement for: ; In the formula, To preset the desired cable output displacement, It is the inverse asymmetric Play hysteresis operator; By superimposing a feedback correction term on the feedforward compensation, the basic motor displacement command is obtained. : ; In the formula, Based on the motor displacement command For proportional gain, For integral gain, The residual tracking error between the expected output and the actual output; Will Replace with the system torque closed-loop output The final motor control command is obtained. : ; Based on the Routh stability criterion, the torque closed-loop gain constraint condition is set: ; ; In the formula, For torque proportional gain, For torque integral gain, For equivalent rotational inertia, For equivalent damping, The cable damping coefficient; The torque proportional gain is obtained by correcting the gain through gain adaptive safety constraints. : ; In the formula, For gain stability margin and ; Based on the preset maximum cable movement speed With system sampling period Calculate the displacement increment limiting threshold : ; Initial control commands for the motor are based on the displacement increment limiting threshold. Incremental clamping and limiting processing is performed to obtain the final system output motor control command. : ; The final control signal is output to drive the exoskeleton.

[0027] Example 2 This embodiment provides a flexible exoskeleton control system with multi-source intent mapping and comfort constraints, including: a signal acquisition module, an intent torque analysis module, a dynamic observation and torque fusion module, an adaptive impedance closed-loop module, and a hysteresis identification and compensation module; The signal acquisition module is used to simultaneously acquire joint motion signals, surface electromyography signals, and muscle sound signals. After preprocessing and feature extraction, it outputs dual-modal physiological feature signals and joint motion information. The intention torque analysis module is used to obtain the motion intention torque based on the dual-modal physiological characteristic signal, complete the electromechanical delay time alignment and safety verification, and output the intention torque with safety constraints. The dynamic observation and torque fusion module uses a nonlinear extended state observer to obtain the estimated value of the interactive torque, corrects the dynamic model parameters online, and weights and fuses the multi-source torques to obtain the fused torque feedback signal. The adaptive impedance closed-loop module is used to construct an impedance model with comfort constraints, adaptively update impedance parameters according to motion conditions, generate torque reference commands by combining the intended torque and the fused torque, and convert them into cable target displacement commands. The hysteresis identification and compensation module is used to identify hysteresis parameters based on the rate-dependent asymmetric RA-PI hysteresis model, construct an inverse hysteresis compensator to feed forward compensation for the target displacement command, and combine displacement closed-loop correction to output motor control signals to drive the flexible exoskeleton.

[0028] This invention addresses the inherent challenges of Bowden cable drives, such as rate-dependent asymmetric hysteresis and the coupling between cable elastic deformation and external disturbances. It achieves compliant ankle-assisted control of an exoskeleton through a hierarchical RA-PI hysteresis model parameter identification, real-time human-machine interaction torque observation based on an extended state observer, and an integrated impedance-adaptive torque closed-loop control scheme. The invention employs a hierarchical optimization strategy combining RLS high-frequency linear parameter iteration and LM low-frequency nonlinear parameter window identification, significantly reducing the computational overhead of joint optimization of all parameters. This enables real-time online parameter self-tuning on an embedded hardware platform and adaptively compensates for model time-varying errors caused by Bowden cable wear and layout path offsets. By relying on a target impedance model combined with multi-source torque closed-loop correction, it optimizes human-machine interaction compliance while accurately following human movement intentions. Four-state gait partitioning allows for differentiated adjustment of impedance parameters at different movement stages, balancing wearability and tracking accuracy. The entire control architecture features rigorous modeling logic and strong engineering feasibility, effectively suppressing the negative impact of transmission hysteresis and unknown disturbances on system control performance, and significantly improving the assisted following performance and environmental adaptability of the flexible exoskeleton.

[0029] Of course, the above description is not intended to limit the present invention, and the present invention is not limited to the examples given above. Any changes, modifications, additions or substitutions made by those skilled in the art within the scope of the present invention should also fall within the protection scope of the present invention.

Claims

1. A flexible exoskeleton control method based on multi-source intent mapping and comfort constraints, characterized in that, Includes the following steps: S1. Simultaneously collect joint motion signals, surface electromyography signals, and muscle sound signals from flexible exoskeleton wearers, perform preprocessing for noise reduction and feature extraction respectively, and obtain dual-modal physiological feature signals and joint angle motion information; S2. Based on the mapping of initial motion intention torque using dual-modal physiological characteristic signals, combined with electromechanical delay completion time alignment, the output intention torque with safety constraints is verified. S3. Establish an exoskeleton joint impedance dynamics model, observe the joint motion state and human-machine interaction torque through a nonlinear expansion state observer, obtain the estimated value of the interaction torque, correct the model parameters online, and obtain the fused torque feedback signal by weighted fusion of multi-source torque information. S4. Construct a target impedance model with wearable comfort quantitative constraints, and adaptively and smoothly update the four impedance parameters, namely virtual inertia, virtual damping, virtual stiffness and impedance mixing weight, according to real-time motion conditions. Perform comfort constraint projection and boundary limiting correction on the impedance parameters. Based on the intention torque and fused feedback torque, synthesize the torque reference command through PI closed-loop correction, and convert it into joint target displacement command after gain stability verification. S5. Based on the rate-dependent asymmetric RA-PI hysteresis model, online identification of hysteresis parameters is performed. Based on the identified hysteresis parameters, an inverse hysteresis compensator is constructed to perform feedforward hysteresis compensation on the target displacement command. After displacement closed-loop correction and amplitude limiting, the motor control signal is output to drive the exoskeleton to complete impedance-torque cascade closed-loop control.

2. The flexible exoskeleton control method based on multi-source intent mapping and comfort constraints according to claim 1, characterized in that, S1 specifically includes collecting surface electromyography (sEMG) signals and myophone (MMG) signals from the target muscle group in the human body. The electromyography (EMG) signal was processed by baseline removal, fourth-order Butterworth bandpass filtering, and power frequency notch filtering to obtain a noise-reduced effective EMG signal. The muscle sound signal is first synthesized into triaxial component magnitudes, and then denoised using a fourth-order Butterworth bandpass filter to obtain the denoised effective muscle sound signal. Extracting bimodal physiological feature signals from high signal-to-noise ratio electromyography (EMG) and muscle sound (SMS) effective signals; Based on the preprocessed dual-modal effective signal, a fixed window length is used. Fixed window step The sliding window mechanism extracts signal data within the window and extracts multi-dimensional signal features. The sEMG time-domain features include root mean square value, waveform length, and zero-crossing rate, calculated as follows: ； In the formula, For the first The root mean square eigenvalue of the sEMG signal under a single muscle channel. For window length, For the time index of discrete sampling points within the window, The starting index of the current window; ； In the formula, For the first Waveform length of sEMG signal under a single muscle channel. For the first Channel, First The noise-reduced effective electromyographic signal at each sampling point; ； In the formula, For the first Zero-crossing rate of sEMG signal under muscle channels For indicator functions, The threshold for determining zero crossing and ; The time-frequency domain characteristics of MMG include the root mean square value and the peak frequency, in order: ； In the formula, For the first The root mean square value of the MMG signal under a single muscle channel; ； In the formula, For the first The peak frequency of MMG signal within a single muscle channel. For discrete frequency points in the frequency domain, and These are the lower and upper frequency limits of the effective frequency band for the MMG signal, respectively. For the first Channel MMG signal at frequency Short-time power spectral density at; Two types of cross-modal coupling core features are constructed, namely sEMG-MMG electromechanical coupling coefficient and electromechanical delay; The electromechanical coupling coefficient of sEMG-MMG for: ； In the formula, The first one in the current sliding window Effective electromyographic signals with noise reduction under muscle channels The mean, The first one in the current sliding window Effective signal for noise reduction muscle sounds under a single muscle channel The mean; The electromechanical delay estimation for: ； In the formula, For the normalized cross-correlation function, For the number of delayed samples, This represents the maximum upper limit for the number of delayed samples; All features from all channels are concatenated to construct a fixed-dimensional bimodal physiological feature vector. : ； The joint angle information refers to the joint angular displacement. Joint angular velocity and angular acceleration .

3. The flexible exoskeleton control method based on multi-source intent mapping and comfort constraints according to claim 1, characterized in that, The initial motion intention torque based on the mapping of bimodal physiological feature signals includes: mapping the bimodal physiological feature vector... The Dimensional components Iteratively update the sliding window and calculate the mean and variance: ； In the formula, For the first The moving average of the dimensional feature components within the current window. Normalized forgetting factor and , For the first The moving average of the dimensional feature components over the previous window. For the first The original feature values ​​of the dimensional feature components in the current window; ； In the formula, For the first The sliding variance of the 3D feature components in the current window. For the first The sliding variance of the feature components in the previous window; Then, adaptive Z-score normalization is performed to obtain the normalized feature vector. ; Calculate the first based on the electromechanical coupling coefficient. sEMG trust level of the channel : ； In the formula, It is the hyperbolic tangent function. For the first Muscle channels, current sliding window The electromechanical coupling coefficient, The electromechanical coupling degradation threshold, The sigmoid sharpness factor and , ; Calculate the first based on the electromechanical coupling coefficient. MMG Trust Level of the Channel : ； Construct the first Channel weighted eigenvectors : ； In the formula, For the first sEMG confidence level of individual muscle channels , , For the first Three sEMG normalized features of each muscle channel For the first MMG trust level for each muscle channel , For the first Two MMG-normalized features of each muscle channel, For the first The sixth normalized feature of the muscle channel For the first The 7th normalized feature of the muscle channel; Fusion feature vectors : ； In the formula, For the first The single-channel weighted feature vectors corresponding to each muscle channel, and , Fixing a single muscle channel includes 7 features. for Total feature dimensions after aggregating all channels; A Gaussian kernel ridge regression nonlinear mapping model was constructed, and offline training was performed using wearer-calibrated action data. The optimal model dual coefficients were obtained by solving for regularization constraints. Calculate the kernel vector Output the initial intended torque .

4. The flexible exoskeleton control method based on multi-source intent mapping and comfort constraints according to claim 1, characterized in that, The intentional torque with safety constraints, output after verification and combining electromechanical delay completion time alignment, includes: According to the obtained Calculate the multi-channel weighted average electromechanical delay : ； right Apply first-order exponential smoothing: ； In the formula, This is the weighted average electromechanical delay value of the current sliding window after smoothing. Smoothing coefficient and , The interval between the previous and current windows is [value]. Weighted average electromechanical delay value after smoothing For window stepping; For joint torque Perform a time right shift operation to obtain the aligned joint torque values. : ； In the formula, The sampling period; right By applying safety saturation limits, the desired torque value of the safety constraint is obtained. : ； In the formula, The upper limit of the safe torque for the joint, This is the saturation limiting function.

5. The flexible exoskeleton control method based on multi-source intent mapping and comfort constraints according to claim 1, characterized in that, The establishment of the exoskeleton joint impedance dynamics model, which involves observing the joint motion state and human-machine interaction torque through a nonlinear expansion state observer to obtain estimated interaction torque values, includes: Input displacement via data acquisition cable With output displacement Differential encoder displacement acquisition and second-order Butterworth low-pass filtering are performed to obtain the output of the first... At each sampling time, the filtered cable input displacement , No. At each sampling time, the output displacement of the cable after filtering is... and the instantaneous motion rate of the cable input end ; A lumped-parameter dynamic model of a single joint of a flexible exoskeleton containing Bowden cable drive was established. The dynamic model was then converted into standard state-space equations to complete the joint dynamics modeling of the exoskeleton. A nonlinear extended state observer was used to monitor the system state and unknown disturbances in real time, and the observed joint angular velocities were obtained by solving the equations. Observation of human-computer interaction torque ; The online correction model parameters include: Calculate the root mean square value of the observed output error. : ； In the formula, For fixed parameters, For the first The observation output error at the sampling time, when Time-triggered parameter correction Adjust the threshold for the parameters; After trigger parameter correction, a joint angle acceleration regression equation is constructed. The Tikhonov regularized least squares method is used to solve the regression equation to obtain the parameter vector estimate. ; Based on parameter vector estimates The estimated equivalent moment of inertia was calculated. Equivalent elastic stiffness estimate Equivalent damping estimate ; By imposing physical constraints on the estimated equivalent moment of inertia, a new equivalent moment of inertia is obtained. : ； In the formula, For the new equivalent moment of inertia, and These are the minimum and maximum constraint values ​​for the equivalent moment of inertia, respectively. A fixed-weight soft update iterative method is used to smoothly update the equivalent rotational inertia, with a soft weighting factor as the weighting factor. Integrating historical equivalent rotational inertia before the update With the new equivalent moment of inertia after constraint treatment The final stable equivalent moment of inertia is obtained. The dynamic model parameters are adaptively corrected online, including soft update of the weight factors. The value is set to 0.8 to complete the iterative update of the equivalent moment of inertia.

6. The flexible exoskeleton control method based on multi-source intent mapping and comfort constraints according to claim 1, characterized in that, The weighted fusion of multi-source torque information to obtain the fused torque feedback signal includes: Two torque information paths are set: Path A is the torque for the intention of active motion. Path B is the estimated value of the interaction torque. This enables independent output through dual channels. Time-scale alignment of the torque information of the two paths: ； In the formula, For the first At the sampling time, after time-scale alignment, the final active motion intention torque is obtained; Calculate the consistency residuals of the two torque paths and dynamic trust : ； ； In the formula, It is the minimum value. The adjustment coefficient of the sigmoid mapping and , For the first The consistency residual obtained after smoothing at the sampling time. The consistency residual threshold and ; Construct a multi-source torque complementary fusion formula: ； In the formula, This is the estimated torque value after fusion; When consistency residual Exceeding the residual safety threshold continuously And continuous sampling period threshold When this occurs, a safety degradation is triggered, locking the online parameter correction update operation and outputting only the fused torque estimate. .

7. The flexible exoskeleton control method based on multi-source intent mapping and comfort constraints according to claim 1, characterized in that, S4 specifically includes: Construct a comfort-constrained impedance model: ； In the formula, For the desired virtual inertia and , This is the second derivative of the joint angle deviation. For the desired virtual damping, The first derivative of the joint angle deviation. For desired virtual stiffness and , For joint angle difference and , This refers to the actual angle of the joint. For joint reference angle, For desired virtual stiffness and , The desired interaction torque; Combined with joint movement speed angular acceleration Temporal characteristics of electromyography Based on a preset angular velocity threshold Angular acceleration threshold Low threshold of electromyography High threshold of electromyography It distinguishes between four human motion states: support, extension, swing, and impact. By presetting the priority of each state and setting minimum dwell time constraints, it suppresses frequent changes in states and outputs stable motion states. ; Based on steady motion conditions The virtual inertia is updated in real time using an exponential smoothing formula: ； In the formula, This is the real-time virtual inertia updated at the current moment. The smoothing coefficient is the parameter and , This is the virtual inertia from the previous moment. Preset values ​​for the baseline virtual inertia corresponding to each motion state; Virtual damping Virtual stiffness Impedance torque mixing ratio Synchronous updates using a smooth formula; Virtual damping Virtual stiffness Perform comfort constraint projection correction on virtual inertia Apply a lower bound constraint, Amplitude is limited to minimum comfortable virtual inertia Above, regarding virtual stiffness Apply upper bound constraints, The amplitude is limited to the maximum upper limit of virtual stiffness. Below; when At that time, check and correct the damping ratio parameter; Based on the actual joint position, reference trajectory, and observed angular velocity, the angular deviation, deviation angular velocity, and deviation angular acceleration are discretized and solved, and then substituted into the target impedance model to calculate the desired interaction torque. Calculate torque tracking error ; ； Resultant torque reference command: ； In the formula, This is the original unlimited torque reference command. For proportional gain, For integral gain, For the impedance torque mixing ratio and , For the first time in history Torque tracking error at sampling points; The final saturation limit is: ； In the formula, This serves as the final effective torque reference command. To achieve maximum safe output torque; by Actual angular displacement of exoskeleton joints As input, the torque command is converted into a cable target displacement command. : ； In the formula, This is the lever arm coefficient. This is the equivalent stiffness coefficient of the cable.

8. The flexible exoskeleton control method based on multi-source intent mapping and comfort constraints according to claim 1, characterized in that, The online identification of hysteresis parameters based on the rate-dependent asymmetric RA-PI hysteresis model, the construction of a hysteresis inverse compensator based on the identified hysteresis parameters, the feedforward hysteresis compensation for the target displacement command, and the generation of motor control commands by combining displacement closed-loop correction include: Configure multiple sets of fixed symmetric thresholds Play operator After recursive calculation, combined with linear gain With weight Weighted combination is used to construct a classic PI hysteresis model and output the predicted hysteresis displacement. ; Construct a rate-dependent asymmetric RA-PI hysteresis model adapted to the drive characteristics of flexible Bowden cables, abandoning the traditional symmetric single threshold model. The format uses independent forward loading thresholds. With reverse unloading threshold Composition of asymmetric threshold pairs This distinguishes the differentiated hysteresis characteristics of the bidirectional motion of cable loading and tensioning, and unloading and releasing; by introducing a value range as... Asymmetric factors Achieve symmetric threshold To asymmetric threshold pairs Adaptive conversion based on asymmetric factors The different values ​​of the cable's bidirectional hysteresis dead zone width are dynamically adjusted to vary. Define cable input displacement increment Determine the real-time movement direction of the cable based on displacement increments. Positive and negative state matching asymmetric Play operator Differentiated output rules, where bidirectional thresholds are satisfied It degenerates into a traditional symmetrical hysteresis form; Introducing instantaneous motion rate at the cable input end Related rate envelope function The hysteresis loop expansion width is dynamically modulated at different motion speeds to characterize the differential characteristics of static hysteresis and high-speed dynamic hysteresis. The rate envelope function... Overall output of the RA-PI hysteresis model They are respectively: ； In the formula, The upper bound of the rate gain and , The rate sensitivity coefficient and ; ； In the formula, For linear gain, For the first The weights of the asymmetric Play operator Let be the rate envelope function. and These are the hysteresis thresholds for the forward loading and reverse unloading of the asymmetric Play operator, respectively. Model output This means outputting displacement at the far end of the Bowden line. The predicted value; The set of parameters to be identified in the RA-PI model Expanded to: ； The RA-PI hysteresis model is linearized, and a regression vector is constructed. With linear parameter vector Combined with the actual output displacement of the cable Constructing identification error Using a forgetting factor The recursive least squares algorithm updates the gain vector periodically. Linear parameter estimation With covariance matrix This enables high-frequency updates of linear parameters; For asymmetric factors Rate correlation coefficient , Every Data from a sliding window is extracted for each sampling period. The optimal nonlinear parameters are solved iteratively using the LM nonlinear least squares algorithm. By integrating the linear and nonlinear parameters, a complete set of RA-PI model identification parameters is obtained. ; A RA-PI inverse hysteresis feedforward compensator is constructed using a parameter analytical inversion method, based on the identified parameters at the current time. , , , , Define the parameters of the inverse model: ； ； ； ； In the formula, This represents the linear gain corresponding to the RA-PI inverse model. In the RA-PI inverse model, the first... The weights of the inverse asymmetric Play operator For the RA-PI model, the first The weights of the asymmetric Play operator For the RA-PI inverse model, the first The loading direction threshold of the asymmetric Play operator, i.e., the forward threshold of the RA-PI inverse model. For the RA-PI inverse model, the first The unloading direction threshold of an asymmetric Play operator, i.e., the negative threshold of the RA-PI inverse model; Based on the solved inverse model parameters, construct a system with inverse threshold pairs. Constrained inverse asymmetric Play hysteresis operator Based on the expected displacement increment The positive and negative states are distinguished by bidirectional operation rules, enabling differentiated inverse hysteresis solutions for cable loading and unloading directions; Solving for the hysteresis feedforward compensation input displacement by combining the results of the inverse asymmetric Play hysteresis operator operation. for: ； In the formula, To preset the desired cable output displacement, It is the inverse asymmetric Play hysteresis operator; By superimposing a feedback correction term on the feedforward compensation, the basic motor displacement command is obtained. : ； In the formula, Based on the motor displacement command For proportional gain, For integral gain, The residual tracking error between the expected output and the actual output; Will Replace with the system torque closed-loop output The final motor control command is obtained. : 。 9. The flexible exoskeleton control method based on multi-source intent mapping and comfort constraints according to claim 1, characterized in that, The motor control signal, after displacement closed-loop correction and amplitude limiting, drives the exoskeleton to complete impedance-torque cascade closed-loop control, including: Based on the Routh stability criterion, the torque closed-loop gain constraint condition is set: ； ； In the formula, For torque proportional gain, For torque integral gain, For the equivalent moment of inertia, For equivalent damping, The cable damping coefficient; The safe torque proportional gain is obtained by correcting the torque proportional gain using the gain adaptive safety constraint formula. : ； In the formula, For gain stability margin and ; Based on the preset maximum cable movement speed With sampling period Calculate the displacement increment limiting threshold : ； Initial control commands for the motor are based on the displacement increment limiting threshold. Incremental clamping and limiting processing is performed to obtain the final system output motor control command. : ； The final control signal is output to drive the exoskeleton.

10. A flexible exoskeleton control system with multi-source intent mapping and comfort constraints, characterized in that, The flexible exoskeleton control method using multi-source intent mapping and comfort constraints as described in any one of claims 1 to 9 includes: a signal acquisition module, an intent torque analysis module, a dynamic observation and torque fusion module, an adaptive impedance closed-loop module, and a hysteresis identification and compensation module. The signal acquisition module is used to simultaneously acquire joint motion signals, surface electromyography signals, and muscle sound signals. After preprocessing and feature extraction, it outputs dual-modal physiological feature signals and joint motion information. The intention torque analysis module is used to obtain the motion intention torque based on the dual-modal physiological characteristic signal, complete the electromechanical delay time alignment and safety verification, and output the intention torque with safety constraints. The dynamic observation and torque fusion module uses a nonlinear extended state observer to obtain the estimated value of the interactive torque, corrects the dynamic model parameters online, and weights and fuses the multi-source torques to obtain the fused torque feedback signal. The adaptive impedance closed-loop module is used to construct an impedance model with comfort constraints, adaptively update impedance parameters according to motion conditions, generate torque reference commands by combining the intended torque and the fused torque, and convert them into cable target displacement commands. The hysteresis identification and compensation module is used to identify hysteresis parameters based on the rate-dependent asymmetric RA-PI hysteresis model, construct an inverse hysteresis compensator to feed forward compensation for the target displacement command, and combine displacement closed-loop correction to output motor control signals to drive the flexible exoskeleton.

Images