2.1.2. Distribution Functions

In Figure 7 in the column of graphs on the left, the temporal variations of the forces that load the supports of the twelve working bodies of the version with a working width of 1 m of the MCLS are shown. Additionally, in Figure 7, the histogram of the frequencies of the force values is represented in the right column. The polygon of frequencies is built using the histograms. When the number of classes in the histograms becomes very large (theoretically to infinity), the frequency polygon tends to reflect the real distribution density of each force. The frequency histograms are constructed with 64 classes using the formula of Mosteller and Tukey [45]. We took into account the length of the data strings obtained experimentally in calculating the number of classes in the histograms.

**Figure 7.** The time dependence of the load in the supports of the active bodies (**left column**) and the frequency distribution of the load values (**right column**). The four rows of figures correspond to the four lines of active organs, starting from the top with the ones closest to the tractor.

The elementary Formula (2) was used to calculate frequency histograms.

$$f\_i = \sum\_{j=1}^{N} H\_i(\mathbf{x}\_j), \ i = \mathbf{1}, \dots, n\_{\text{lb}}, j = \mathbf{1}, \dots, N \tag{2}$$

where:

$$H\_{\mathbf{i}}(\mathbf{x}\_{\mathbf{j}}) = \begin{cases} \mathbf{1}, \; h\_{\mathbf{i}} \le \mathbf{x}\_{\mathbf{j}} < h\_{\mathbf{i}+1} \\ \mathbf{0}, \mathbf{x}\_{\mathbf{j}} < h\_{\mathbf{i}} \text{ or } \mathbf{x}\_{\mathbf{j}} \ge h\_{\mathbf{i}+1} \end{cases}, \; \mathbf{i} = \mathbf{1}, \ldots, n\_{\mathbf{h}\prime} \quad \mathbf{j} = \mathbf{1}, \ldots, \mathbf{N}, \tag{3}$$

$$H\_{\mathfrak{n}\_{\mathbf{h}}}(\mathbf{x}\_{\mathbf{j}}) = \begin{cases} \mathbf{1}, \, \mathfrak{h}\_{\mathfrak{n}\_{\mathbf{h}}-\mathbf{1}} \le \mathbf{x}\_{\mathbf{j}} \le \mathfrak{h}\_{\mathfrak{n}\_{\mathbf{h}}}\\ \mathbf{0}, \mathbf{x}\_{\mathbf{j}} \langle \mathfrak{h}\_{\mathfrak{n}\_{\mathbf{h}}-\mathbf{1}} \text{ or } \mathbf{x}\_{\mathbf{j}} \rangle \mathbf{h}\_{\mathfrak{n}\_{\mathbf{h}}} \, \, \, \mathbf{j} = \mathbf{1}, \ldots, N \end{cases} \tag{4}$$

*N* is the number of observations from the examined random sequence, and *nh* is the number of classes in the frequency histogram. Then the probability density histogram was calculated according to Formula (5):

$$\rho\_{\mathbf{i}} = \frac{\mathbf{f}\_{\mathbf{i}}}{\mathbf{N}'} \text{ i } = 1, \dots, r \text{ } \mathbf{n}\_{\mathbf{h}} \tag{5}$$

where {*ρi*}*i*=1,...,*nh* is the series of probability densities. The series of probabilities (or cumulative probabilities) is calculated using Formula (6).

$$p\_i = \sum\_{k=1}^{l} \rho\_{k\prime} \quad i = 1, \ldots, n\_h \tag{6}$$

In Figure 8, the probability density and the probability that the random variable (the series of numerical values over time recorded at each of the twelve measurement locations) will take different values in the experimental range are shown.

**Figure 8.** The probability density and the probability of taking certain values from the experimental range for the forces recorded in the twelve measurement locations.

The probability densities calculated as shown above do not resemble those found in statistics books or in articles in which the experimental results are approximated, by hypothesis, with one of the classical statistical distributions (normal, exponential, Student, Fischer, etc.). Probability distributions differ from the statistical densities with which experimental curves are usually approximated. Precision is lost by approximating with the classical statistical distributions of the probability densities, and we cannot precisely estimate the local error in relation to the experimental data. For this reason, we preferred, considering our goals, to interpolate the empirical probability densities through cubic spline

functions since we wanted maximum precision and were not interested in using the results in immediate theoretic models. Regarding the probability curves, they have shapes very similar to the classical ones, but for the same reasons, we proceeded with the probabilities by interpolating with spline functions. In Figure 8, the experimental curves interpolated by spline functions are given in each graph for the twelve series of 4000 numerical data points.

The probabilities interpolated with spline functions are used to calculate the probability of the occurrence of a force greater than a fixed force (related to the limit resistance characteristics of the material of the supports of the working body).

The distribution functions of the twelve random variables described by the numerical sequences in Figure 7 are calculated in the same way as the probability densities and probabilities whose graphic representations are given in Figure 8. Therefore, the distribution functions are expressed as interpolations by cubic spline functions, for which the graphic representations are given in Figure 9. The graphic representations of the twelve distribution functions of the sequences of numerical values coming from the measurement points located on the supports of the twelve working bodies are grouped three by one in the graph, corresponding to the three working bodies in each line, with the counting starting from the back of the tractor.

**Figure 9.** The graphic representation of the distribution functions of the twelve numerical data sequences was obtained by measuring at the locations indicated in Figure 5.

For complete agreement with the theoretical definition, the interpolated distribution functions are extended with the value zero to the left of the minimum value and with the value one to the right of the maximum value of the corresponding random variable.

#### 2.1.3. Characteristic Functions of Random Vibrations

From a theoretical point of view, the first tested properties of random vibrations are *stationarity* and *ergodicity*. The definitions of these notions can be found in all the literature dedicated to random vibrations, for example, in [2,4–9,42,43,46]. In general, the response of structures to random vibrations is studied through experiments that have strictly elaborated procedural standards [47].

A synthetic characterisation of random signals is given in [3]. The non-deterministic signal has specific characteristics: mean, dispersion, global average, global dispersion, histogram, power spectral density, etc. The signal can have a certain degree of predictability in its evolution over time. Depending on specific characteristics, the non-deterministic signal can be:

*Stationary*, if the mean and dispersion do not depend on time but are constant;

*Ergodic,* if the mean per portion does not differ from the global mean;

*White noise*, with a constant power spectral density throughout the frequency band.

Testing the stationarity and ergodicity of the experimental signals (sequences of finite length) cannot be conducted using the definitions of the mean value and the autocorrelation function of the random variables, for example, in [42]. The definitions in [42] assume processes for crossing the limit after the number of samples. We do not have an infinite number of samples, and if we did, we would be constrained by costs to limit the duration of the experiences as much as possible. For this reason, the calculations are made for a finite number of samples, limited by the number of samples in the entire sequence.

In Figure 10, we see the behaviour of the average value of the random variable given by the registration with the ch4 code (coming from the working body in the first line after the tractor, on the extreme left). The tendency to decrease the amplitude is observed with an increase in the number of samples in the sequence considered for calculating the average value. Therefore, one can only suspect asymptotic behaviour, but such behaviour cannot be stated strictly theoretically.

**Figure 10.** The behaviour of average values increases as the number of samples used in their calculation increases.

In Figure 11, the terms below the limit of the definition of the autocorrelation function are calculated for the subsequence of the ch4 signal, with the graphic representation limited to only five values of the gap. There is a tendency to decrease the amplitude of the autocorrelation function with the gap value and a weak asymptotic tendency, which does not allow us to make statements about the stationarity and/or ergodicity of the signal.

**Figure 11.** Autocorrelations of the ch4 signal for five offset (*τ*) values.

In Figure 12, the autocorrelation curves are given for each of the twelve signals collected from the supports of the working bodies of the 1 m working width version of the MCLS. The curves are calculated using the *lcorr* function of the programme [48]. The calculation can also be performed directly by programming the three simple formulas given in [49]. These curves can be considered substitutes in the field of real signals of finite length for the assessment of the stationarity or ergodicity of signals. In Figure 12, it is first observed that all twelve signals tend to asymptote in time towards the value 0. First, this means that the signals are weakly autocorrelated (another argument for their randomness). The fact that they are not constant in time shows that the signals are neither stationary nor ergodic.

**Figure 12.** The autocorrelations of the twelve signals recorded in the experiment were calculated using [50].
