At the normal hill (hill size = 86 m) in Akita Japan, take-offs of eight male professional jumpers (23 ± 7.7 years old, 176 ± 4.3 cm, 62 ± 6.5 kg) from Chinese national ski jumping team were filmed after being informed. The trial protocol was reviewed and approved by the ethics committee of Capital University of Physical Education and Sports (approval number: 2022A41). Twenty-seven jumps (three or four jumps per participant) were recorded at 120 fps using two cameras positioning in front of the table edge. A 3D radial calibration frame (QFS-28 DLT Calibration Frame, HuiAnMing Sciences Co., Ltd., Beijing, China) was set up on the table and filmed for camera calibration. The equipment layout is shown in Figure 1A.

Videos from two cameras were synchronized based on the frame in which jumpers’ heels left the table edge. Video digitization was implemented through an automatic algorithm assisted with manual corrections. A 2D pose estimation algorithm (Motion-3D V1.2.5, Fastmove, Dalian, China) was used to identify 20 joint points of human body automatically (Figure 1B). Joint points with weak estimation accuracy were then corrected manually. Based on images of the calibration frame, 3D trajectories of 20 joint points were calculated by triangulation and will be used as the raw kinematics data for the following musculoskeletal simulation.

It is critical for our simulation to use a suitable musculoskeletal model. The model developed by Lai et al. (Lai et al., 2017) was used as base model because of its’ refined lower extremity musculature. This model modified the range of knee flexion and properties of several muscles on the basis of the model published by Rajagopal et al., making it suitable for movements involving substantial knee and knee flexion, such as ski jumping take-off. In addition, since the external forces in this simulation would be calculated based on full-body kinematics, the model we used should also be capable of full-body movements tracking. Thus, the lumbar spine of the FBLS model (Raabe and Chaudhari, 2016), the head and neck of the MASI model (Cazzola et al., 2017), and the shoulder construction from the model developed by Blana et al. (Blana et al., 2008) were combined into the original Lai model to enable it to track full-body movements. While combining these models, inertial properties of each body were also adjusted to meet the mass distribution of a human body. Additionally, since the raw kinematic data was not enough to determine the movements of subtalar and metatarsal-phalangeal joints, these two joints were locked in the simulation.

The general model constructed above was then scaled based on athletes’ height and weight through scaling module of OpenSim (Scott L et al., 2007). Since these professional athletes have higher muscle strength than ordinary people, it is necessary to increase the muscle strength. The maximum isometric forces of all muscles were first doubled and then multiplied by a scaling factor based on the height and weight of each athlete (Handsfield et al., 2014).

Due to equipment limitations, we were unable to measure the GRF directly during the take-off. However, we developed a general method based on inverse dynamics iteration to estimate the GRF from kinematic data. Given kinematics and kinetics, inverse dynamics solves the equations of motion and calculates generalized forces (Kuo, 1998), which include the compensating forces and torques directly acting on the pelvis to satisfy Newton’s second law, called residuals. These residuals represent dynamic inconsistencies between kinematics and kinetics, and consist of six components, namely, the components of force and torque along three orthogonal axes, denoted by F x r e s , F y r e s , F z r e s , M x r e s , M y r e s , and M z r e s . Assuming equal external loads on both legs, the GRF also consists of six components, denoted by F x l / r , F y l / r , F z l / r , M x l / r , M y l / r , and M z l / r respectively, where the superscripts l / r represent the left or right leg. By iteratively performing inverse dynamics and adding the residuals to the GRF, we were able to obtain a GRF that matches the kinematic data and minimizes residuals. Initially, all GRF components are set to 0 and applied at proximal second metatarsal, and inverse dynamic analysis is performed to calculate the residuals. Then, the following formula is performed to add the residual components to the GRF components of the left and right legs evenly and update the GRF: F x l / r , s = F x l / r , s − 1 + 1 2 F x r e s , s F y l / r , s = F y l / r , s − 1 + 1 2 F y r e s , s F z l / r , s = F z l / r , s − 1 + 1 2 F z r e s , s M x l / r , s = M x l / r , s − 1 + 1 2 M x r e s , s M y l / r , s = M y l / r , s − 1 + 1 2 M y r e s , s M z l / r , s = M z l / r , s − 1 + 1 2 M z r e s , s (1) where superscript s is the index of iterative steps. In each iteration, the mean root mean square (RMS) of residuals through the time period of simulation is calculated. Finally, when the mean RMS is less than the threshold 0.001 (N for forces and Nm for torques), the residual is considered close enough to 0 and the GRF is output as result.

In order to verify this method under standard in-lab motion capture trial procedure, an imitation jump trial was implemented. After being informed the trial procedure, three healthy adult men (29 ± 2.2 years old, 172 ± 0.04 cm, 83 ± 2.83 kg) from local college agreed to participate. The trial protocol was reviewed and approved by the ethics committee of Capital University of Physical Education and Sports (approval number: 2022A41). Participants were asked to perform 3 jumps each. Thirty-nine reflective markers were attached to multiple anatomical landmarks all over the body. Raw marker kinematics data were collected by a 12-camera motion capture system (Optitrack Prime 22, NaturalPoint, Corvallis, United States) with 200 fps, while the GRF were collected synchronously with 1,000 Hz by two force plates (Multicomponent Force Plate 9260AA, Kistler, Winterthur, Switzerland). GRF estimated by the above method was compared with the result of force plates.

The high-speed nature of ski jumping makes the aerodynamic effects significant enough to be considered in the simulation. In our case, aerodynamic forces during take-off were estimated based on body kinematics and CFD simulations.

CFD simulations were carried out under various possible ski jumping postures. Seven representative posture angles with a significant impact on jumper aerodynamic properties were selected (Figure 2A, including the angle between legs θ 2 , the angle between skis θ 1 , the angle between arms and trunk δ 4 , the angle between trunk and thighs δ 3 , knee angle δ 2 , ankle angle δ 1 , and the angle of attack α ). Based on jumpers’ possible postures during the take-off and early flight phases, several values were selected for each angle to form a total of 278 postures (Table 1). Under each posture, geometric models of the jumper and skis were constructed using Gambit (V 2.4.6). The jumper model was 1.70 m high. The size of the skis was 2.57 × 0.22 m. Triangular unstructured surface grid was generated over the model surface. A tetrahedral unstructured grid was then generated in the computational domain with a size of 15 m (x) × 15 m (y) × 8 m (z) (Figure 2B). There were about 1.5 million grids in total. No-slip wall, velocity inlet, and pressure outlet boundary conditions (Figure 2C) were set up in Fluent (V 16.1). The inlet flow velocity was set to be 23 m/s (according to FIS Certificates of the jumping hill). The incompressible flow model was adopted. The mathematical model adopts the Reynolds-averaged N-S equation and the k-ε turbulence model. Using the first-order upwind style for spatial discretization, the SIMPLE algorithm was used to iteratively compute until convergence. Aerodynamic forces and moments were calculated by integrating the pressure field over the model surface. The model was also divided into 10 parts: left and right calves, thighs, forearms, and upper arms as well as trunk and head, to integrate separately. All aerodynamic forces and moments were imported into a ski jumping aerodynamic database.

Aerodynamic forces were estimated by matching the measured body motion sequences with the postures in the database. The above posture angles during the take-off were calculated based on model and marker set kinematics. Wherein, ski postures were calculated by markers on heels and toes. The angles in each frame were linearly interpolated in the multidimensional posture space to calculate the corresponding aerodynamic force. Flow velocity and body size were set as constant values in the above CFD simulations. To reflect these effects, the following empirical formula was used: F A D = ρ 2 c A w 2 (2) where F A D is the aerodynamic force, ρ the air density, c the aerodynamic coefficient which is related to the body postures and clothing, A the cross-sectional area, and w the speed of the body. A is assumed to be linearly related to the height of the athlete. Therefore, the aerodynamic forces calculated above were scaled linearly by height and quadratically by velocity. The resulting aerodynamics were used as external loads in musculoskeletal simulations as shown in Figure 3.

All data processing procedures were implemented through OpenSim 4.1 (Stanford University, Stanford, United States) shown in Figure 4. Movements within 0.1s before release were simulated. Firstly, full-body kinematics were calculated from joint points trajectories through inverse kinematics (Lu and O'Connor, 1999). An example of model motion sequence from inverse kinematics is shown in Figure 5. Then, GRFs and aerodynamic forces were calculated through the method established above. The GRF components in the medial-lateral direction and all moment components were ignored, which allows F z l / r , M x l / r , M y l / r , and M z l / r being considered equal to zero. F y l / r was calculated iteratively by inverse dynamics. F x l / r was the friction between the skis and table with a friction coefficient of 0.018 (calculated by a friction angle of 1° according to the regulations of the International Ski Federation) (Gasser, 2018). Finally, muscle activation and force were calculated through static optimization using a criterion of minimizing the sum of squared muscle activation (Thelen and Anderson, 2006). The torque contributions of these muscles to the lower extremity joints were further calculated by multiplying muscle force by moment arm.

Measured and estimated GRFs of simulated jumps were normalized to body weight (Figure 6). They were similar in the anterior-posterior and vertical directions, while the differences were larger in the medial-lateral direction. The mean root mean square error (RMSE) between the measured and estimated GRFs was 39.92 N, 106.41 N, and 19.44 N in the anterior-posterior, vertical, and medial-lateral directions, respectively.

The mean RMSE between the coordinates of the simulated markers derived from inverse kinematics and the measured one was 3.13 cm, which was slightly larger than the reported result using a marker-based motion capture system (1.38–2.03 cm) (Glover et al., 2021).

Take-off GRF (normalized to body weight) component perpendicular to the table from our calculation was reported and compared with the GRF from 2 previous studies (Virmavirta and Koml, 1993b; Kaps et al., 1997) (Figure 7). The take-off GRFs all share a similar trend, that is, decrease as the take-off process. The GRF from our estimation reduced more rapidly than that of Virmavirta and Komi.

The calculated take-off aerodynamic resultant force (mean ± S.D.) was shown in Figure 8A, where F x and F y represent the components of the force in the x and y directions (the coordinate system is defined in Figure 1A), respectively. The drag force along the direction of the platform gradually increased from 60 N to about 90 N, while the drag force perpendicular to the table was about 10 N–20 N. The aerodynamic forces (mean ± S.D.) acting on each body part are shown in Figures 8B–F. Aerodynamic force acts on the trunk and on calves and skis primarily.

Simulated mean muscle activation was reported and compared with EMG (Figure 9) from study of Virmavirta et al. (Virmavirta et al., 2001) (EMG of jumper VML and ML on hill K-90) and that of Virmavirta and Komi (Virmavirta and Komi, 2001) (mean EMG of ten jumpers on hill K-100). The simulated activation of the tibialis anterior showed a similar trend to both the EMG data of the jumper VML and the mean EMG data of ten athletes (Figure 9A). Similarly, the simulated activation of the gluteus maximus was similar to the EMG data of jumper ML (Figure 9B), and the activation of the biceps femoris long head was similar to the mean EMG of ten athletes (Figure 9E). As for other muscles, the simulated activation was less consistent with EMG.

Simulated activation of selected lower limb muscles was shown in Figure 10. Muscles with the function of lower limb extension were selected and grouped as ankle plantar flexors, knee extensors, and hip extensors to plot. Activation level of the ankle plantar flexors was relatively low (Figure 10A). The soleus was activated until the end of the take-off, where the gastrocnemius began to activate slightly. The activation level of quadriceps decreased with take-off progress (Figure 10B). Similarly, the activation level of the hip extensors also decreased with take-off progress, among which lower bundles of gluteus maximus and biceps femoris long head showed higher activation level than the others (Figure 10C).

The torque contributions of the lower extremity muscles to the ankle, knee, and hip joints are shown in Figure 11. The relative magnitudes of the extension torques produced by the ankle plantar flexors and knee extensors on their corresponding joints (Figures 11A,B) were approximately consistent with the relative magnitudes of their activation levels. As for the hip joint (Figure 11C), the extension torque produced by the long head of the biceps femoris, gluteus maximus, and semimembranosus were similar, although the activation level of the biceps femoris was higher than the others. In addition, the biceps femoris long head and semimembranosus were found to generate large flexion torques on the knee, while the rectus femoris produces large flexion torque on the hip.