Before converting time domain signals into frequency domain signals, some data preprocessing steps are usually required to ensure accurate and meaningful frequency domain representations, which mainly include data normalization and Blackman-Harris window.

1) Data normalization: The data normalization of the time domain signal can ensure that the amplitude of the signal is within a reasonable range and avoid the numerical stability problems caused by too large or too small signal amplitude. In this paper, the Min-Max normalization method as Formula (1) is adopted, which scales the data to a specified range by linear transformation.

After data normalization, the Blackman-Harris window function is applied to the time-domain signal, which is equivalent to introducing a window in the time-domain signal to reduce the discontinuity of the signal at the window boundary, thereby reducing spectrum leakage and improving the accuracy of spectrum analysis. The steps usually involve point-by-point multiplication, that is, multiplying each sample of the signal with the corresponding sample of the window function.

The general form of a Blackman-Harris window is as Formula (2) is adopted: w n = a 0 − a 1 ⁡ cos 2 π n N − 1 + a 2 ⁡ cos 4 π n N − 1 − a 3 ⁡ cos 6 π n N − 1 (2) where N is the window length, n is the index of the sample in the window, n = 1 , 2 , … , N , and a 0 , a 1 , a 2 , and a 3 are the coefficients.

Suppose there is a time domain signal x n , and the signal y n after applying the window function can be expressed as Formula (3): y n = x n ⋅ w n (3) where ⋅ stands for point-by-point multiplication.

The conversion of the time-domain signal to the frequency-domain representation helps to better understand the frequency component of the signal. Frequency domain analysis has several advantages: First, it provides a more compact representation of the main signal components. Second, it allows the detection of attacks that change the frequency of the usual operating pattern. Finally, attacks that typically evade existing time-domain detection methods are more easily detected in frequency-domain analysis.

Compared with other methods that convert time domain signals into frequency domain signals, Mixed-Radix FFT is a flexible and efficient discrete Fourier transform (DFT) calculation method. By decomgenerating the length of DFT into the product of different prime factors, Mixed-Radix FFT can make more effective use of computing resources, reduce computing complexity and improve computing efficiency, this makes it potentially more advantageous in real-time applications, especially for cases where fast calculations are required, and in addition, the Mixed-Radix FFT can handle transformations of multiple lengths and is not limited to dealing with DFT whose length is a power of 2, so it is more flexible.

The basic formula of DFT is to transform the time domain signal into the frequency domain signal by Fourier transform of the signal sequence. The DFT formula is as Formula (4): X k = ∑ m = 0 M − 1 x m ⋅ e − 2 π i M ⋅ k ⋅ m (4) where x m represents the time domain sequence, X k represents the frequency domain sequence, m represents the index in the time domain, k represents the index in the frequency domain, M represents the length of the input sequence, which is also the length of the output sequence of the DFT, i is an imaginary unit, full i 2 = − 1 , e − 2 π i M ⋅ k ⋅ m is the rotation factor.

The idea of Mixed-Radix FFT is to divide the calculation task of DFT into smaller subtasks to reduce the complexity of calculation, and it realizes the DFT calculation of the whole signal by recursively applying decomposition and Butterfly Operation. The key steps of the algorithm include:

Decomposition: Decomposition of a DFT of length R into multiple smaller DFTS, such as Formula (5): R = R 1 ⋅ R 2 ⋅ … ⋅ R l (5)

Groups: The calculation of the DFT is broken down into multiple small-scale DFT calculations. Each small DFT length is one of the prime factors obtained by decomposition. A set of radix is selected, usually 2, 3, and 5 are selected as radix, so that R i can be an integer power of 2, 3, and 5 respectively.

Butterfly Operation: Each small DFT is calculated using the butterfly operation. The butterfly operation involves multiplying and adding complex numbers, and its specific form depends on the length R i of a small DFT. Each Butterfly Operation involves two inputs and two outputs, one of which is the real part and the other is the imaginary part. The mathematical representation is as Formula (6): Y 0 = X 0 + W R l ⋅ X 1 Y 1 = X 0 − W R l ⋅ X 1 (6) where W R l is the rotation factor, which can be calculated from W R l = e − 2 π i R ⋅ l . Y 0 and Y 1 are the outputs of the butterfly operation, X 0 and X 1 are the inputs of the butterfly operation, and l is the frequency index of the current calculation.

Combination: The results of all small-scale DFT calculations are combined to get the final DFT result. This usually involves an appropriate weighted sum.

Before the data is input into the machine learning model, PCA is used to reduce the dimension of the data set, which can remove redundant information in the data, retain the most important information in the data set, and reduce the training time of the model, high-dimensional data sets are usually accompanied by more computing overhead, and through dimensionality reduction, the training process of the model can be accelerated. Overall, PCA can help simplify data, improve model performance, reduce computational costs, and provide better interpretability. In addition, PCA assumes that the main information of the data is concentrated in the direction of large variance, and can better play its advantages of reducing the data dimension and extracting the main information when processing the data with strong linear correlation. The process of converting the time domain signal to the frequency domain signal is usually achieved through the FFT, which is linear, so the obtained frequency domain signal can be considered as the data of linear structure. The basic steps of PCA are as follows:

Data standardization: The data is standardized to ensure that each feature contributes equally to the principal component.

Calculate the covariance matrix: Calculate the covariance matrix of the normalized data. The covariance matrix reflects the correlation between different features as Formula (7). C o v X , Y = ∑ i = 1 d x i − X ¯ y i − Y ¯ d − 1 (7) where, X and Y are two features, x i and y i are the eigenvalues of the i samples of the two characteristics respectively, X ¯ and Y ¯ are their average values respectively, and d is the number of samples.

Calculation of eigenvalues and eigenvectors: The eigenvalue decomposition of covariance matrix is carried out to obtain the eigenvalues and corresponding eigenvectors. The eigenvalue of the covariance matrix represents the variance of the data in the direction of the corresponding eigenvector. The eigenvectors represent these directions. The idea of PCA is to select the eigenvector corresponding to the maximum eigenvalue, that is, to select the direction that can maintain the variance of the original data to the greatest extent.

Select principal components: The eigenvalues are arranged in order from largest to smallest, and the eigenvector with the largest first h eigenvalues is selected as the principal component.

Projection: The original data is projected onto the selected principal component to obtain a reduced data set.

In order to capture the potential relationship between different network traffic data and the frequency dependency in the data, and realize efficient detection of network attacks and accurate classification of network attack types, we proposed a Dual-Head Output model based on parallel CNN-BiLSTM, as shown in Figure 2. The parallel CNN-BiLSTM network, as a shared backbone network of multi-head output structures, can simultaneously process different parts of the input data and make full use of the capability of parallel computation, which significantly improves the computational efficiency and speeds up the model training and reasoning process. In addition, the parallel CNN and BiLSTM layers help to extract and integrate data features simultaneously, which allows us to capture information related to various aspects of the data and provide efficient model representations that improve prediction accuracy. Specifically, the CNN component is used to extract inherent features between different data types within a certain number of frequencies. At the same time, BiLSTM captures deeper frequency features by considering information in both “forward” and “backward” directions. The parallel architecture of CNN and BiLSTM allows the independent extraction of inherent features from various data types and then concatenation of these features into a final feature vector. The specially designed Dual-Head Output structure shares the underlying feature representation extracted from the input data by the parallel CNN-BiLSTM network, which helps to improve the efficiency and generalization performance of the model, and can accurately classify the types of network attacks while realizing the detection of network attacks. In order to improve the accuracy of the model for attack detection and type classification at the same time, it is necessary to customize the degree of attention for different tasks. We introduce a Self-Attention Mechanism in the middle layer of the two output headers respectively. The Self-Attention Mechanism is suitable for the task interested in the relationship between the elements in the input sequence, and can capture the global dependency. The middle layer usually contains more information and abstract features, and the Self-Attention Mechanism can flexibly and dynamically focus on the task-related part of the feature vector output from the backbone network, which can effectively improve the performance of the model on different tasks.

Remark 1

Multi-head output networks refer to neural network architectures that incorporate multiple output heads, each dedicated to a specific task or objective. Multi-task learning involves training a single model to perform multiple tasks concurrently. The key distinction lies in their architectural design and training strategies. Multi-head Output Networks adopts a modularized structure where each task corresponds to its own output head. During training, the model optimizes the losses from all output heads jointly. This modular design is advantageous when tasks are independent or exhibit significant differences, as it allows for task-specific fine-tuning without compromising the shared foundation. In contrast, multi-task learning takes a holistic approach by jointly training the model on all tasks. The model optimizes a joint objective that encompasses all tasks, promoting the discovery of shared representations. This holistic approach is particularly effective when tasks are interrelated or share underlying structures, fostering a collaborative learning process that benefits multiple objectives simultaneously.

In order to improve the convergence speed of the neural network, reduce the sensitivity to the input feature scale, and improve the generalization ability of the model. When using a deep learning model, the data whose dimensionality has been reduced by PCA needs to be normalized again to ensure that the model can learn and adapt better. The detailed operations were shown in Formula 1.

The normalized data X n o r is modeled as tensor X i n , and two corresponding data labels are modeled, including attack detection label Y d and attack type classification label Y c , which correspond to the attack detection output Y ^ d and type classification output Y ^ c of Dual-Head Output respectively, as shown in the Formula (8)–Formula (13). X n o r = X i n 1 , X i n 2 , … , X i n R T (8) X i n = X i n 1 , X i n 2 , … , X i n S T (9) Y d = Y d 1 , Y d 2 , … , Y d S T (10) Y c = Y c 1 , Y c 2 , … , Y c S T (11) Y ^ d = Y ^ d 1 , Y ^ d 2 , … , Y ^ d S T (12) Y ^ c = Y ^ c 1 , Y ^ c 2 , … , Y ^ c S T (13) where X n o r ∈ R R × h , R represents the first dimension of the time domain data, that is, the total number of samples, and the first dimension of the frequency domain data transformed by DFT, that is, the total number of frequencies; h is the number of data features retained after PCA reduces the data dimension. X i n s ∈ R L × h is the s -th element of the data set X i n ∈ R S × L × h that will be fed to the neural network, Data set X i n is divided into S samples, and each sample contains L frequency numbers, It should be noted that when each sample is input into CNN network, it needs to be transposed and then input into BiLSTM network at the same time. Y d ∈ R S × 1 is a set of numbers, each set of numbers has S numbers, and their values are 0 or 1, 0 indicates that the network is not attacked, and 1 indicates that the network is attacked. Y c ∈ R S × b contains S groups of numbers, each group of numbers has b values of 0 or 1, the e -th number in each group corresponds to the e -th type of network attack, e = 1 , 2 , … , b , there are b types of network attacks the network may suffer, 0 indicates that the network is not attacked by the e -th type, and 1 indicates that the network is attacked by the e -th type. Y ^ d s ∈ R 1 × 1 means that if you input the s -th sample into the neural network model, it will output a real number between 0 and 1, according to the size of this real number, determine whether the network is under attack. Y ^ c s ∈ R 1 × b means that feeding the s -th sample to the neural network model will output b real numbers between 0 and 1 and add up to 1, according to the size of these real numbers to determine what kind of network attacks the network may suffer, the larger the number, the more likely the network is to be subjected to this type of attack.

Suppose that samples with a batch size of f are input into the neural network each time, f ≤ S , then the attack detection label and attack classification label of the i -th sample are y d i ∈ R 1 × 1 and y c i ∈ R 1 × b respectively, and the output value of attack detection and attack classification are y ^ d i ∈ R 1 × 1 and y ^ c i ∈ R 1 × b respectively, then the loss function L d of attack detection and the loss function L c of attack type classification are set to as Formula (14), Formula (15): L d = − 1 f ∑ i = 1 f y d i ⁡ log y ^ d i + 1 − y d i log 1 − y ^ d i (14) L c = − 1 f ∑ i = 1 f ∑ j = 1 b y c i j ⁡ log y ^ c i j (15) where y c i j ∈ R 1 × 1 represents the label of the j -th attack type in the i -th sample in each batch, and y ^ c i j ∈ R 1 × 1 represents the output value of the j -th attack type corresponding to the i -th sample in each batch.

We set a hyperparameter p as the threshold, 0 ≤ p ≤ 1 , when y ^ d i > p , the intermediate parameter g is set to 1, which proves that the network is under attack, otherwise it is set to 0, which proves that the network is not under attack. Let q = ⁡ max y ^ c i j , j = 1 , 2 , … , b , q corresponds to the j -th type of network attack, let v = q ⋅ g , Only when v ≠ 0 occurs can it be proven that the network has experienced the j -th type of network attack. Figure 3 shows the above determination process.

The data acquisition environment is shown in Figure 4, is composed of two parts, the physical topology and the communication topology, which ensure the two-way transmission of energy and information respectively (Ren et al., 2023; Teng et al., 2023). The physical topology consists of 17 units. It includes three energy carriers of electricity, gas and heat, four energy loads of electricity load, gas load, heat load and colding load, three energy storage devices of electricity storage, heat storage and gas storage, four energy conversion devices of electrical refrigerant, electrical boiler, CHP and refrigerating machine, but also through photovoltaic panels and solar water heaters to use the renewable energy of solar energy, in addition, the voltage is also changed by transformers.

In the multi-energy system network, the collection of network traffic data mainly involves the equipment and protocol related to network traffic. Multi-energy system communication protocols include Message Queuing Telemetry Transport (MQTT), Constrained Application Protocol (CoAP), OPC Unified Architecture (OPC UA), and Distributed Network Protocol 3 (DNP3). The main communication protocol related to network security is NetFlow, which is used for network traffic monitoring and analysis. In the data collection process of this article, we used a Netflow-enabled switch (HUAWEI-LSS7G48TX6E0), a hardware device that generates NetFlow data and provides detailed information about network traffic. Network traffic data is collected using Switched Port Analyzer (SPAN), which is a local traffic mirroring technology usually provided by a switch. With SPAN, we can select multiple source ports and then copy all traffic on those ports to a destination port dedicated to monitoring. This allows us to get all the data flows through the source port.

To collect normal and abnormal traffic, we configure mirroring on four Netflow-enabled switches and prepare four computers to receive the mirrored traffic. Configure SPAN on the logged in switch, connect the computer’s network adapter to the destination port configured with SPAN, and then use Wireshark to extract traffic from each network. Because normal traffic and attack traffic are often mixed together in real networks, in order to better simulate the real environment, it helps the model better learn to work in the real environment. For normal traffic data and attack traffic data, we do not collect them separately, but collect them successively, and label the collected data to make it clear whether each sample is normal traffic or attack traffic, so as to facilitate subsequent model training and evaluation. It is important to note that we only conduct one type of cyber attack at a time, that is, the attack traffic data we collect contains information about only one type of cyber attack, not a mixture of attacks.

In this study, we choose to generate and collect network traffic data involving common port scanning attack, DoS attack and FDI attack. When the attack data needs to be collected, the attack traffic can be generated by connecting the attack host to the data collection environment. The information on attack traffic is summarized in Table 1.

1) For port scanning attack, we chose the easy-to-use nessus with automated scanning and rules engine and the open source NMAP with multi-platform support, the attack types include SYN Scan, UDP Scan, XMAS Scan, and ACK Scan.

In order to further verify the effectiveness of the proposed method, an ablation analysis was carried out in this case. Table 6; Figure 13 compare the evaluation metrics of attack detection in these different cases, namely, proposed method, no Blackman-Harris window (Case 1), no time-frequency domain transformation based on Mixed-Radix FFT (Case 2), and no PCA dimensionality reduction (Case 3). Table 7; Figure 14 compare the evaluation metrics of attack type classification in different cases, including port scanning attack, DoS attack and FDI attack. The following conclusions can be drawn:

The results show that our method has the highest evaluation metrics compared with no Blackman-Harris window, no time-frequency domain transformation based on Mixed-Radix FFT and no PCA dimensionality reduction. This is because Blackman-Harris window can effectively reduce spectrum leakage and improve the accuracy of time-frequency domain transformation. The Mixed-Radix FFT method has both accuracy and computational efficiency, and provides more accurate time-frequency information for complex signals. PCA dimensionality reduction helps to reduce the input feature dimension, improve the efficiency of model training, and filter out redundant information, so as to optimize the performance of the neural network and make it more suitable for processing the data after time-frequency domain transformation. This integrated approach can improve the accuracy of network attack detection and attack type classification.

The evaluation metrics of Case 2 is slightly lower than that of Case 1, and Case 3 has the lowest evaluation metrics. It can be concluded that PCA dimensionality reduction has the greatest impact on the accuracy of network attack detection and attack type classification, because PCA dimensionality reduction directly reduces the number of features and can significantly affect the accuracy. The impact of the time-frequency domain transformation based on Mixed-Radix FFT and Blackman-Harris window on the accuracy of network attack detection and attack type classification is similar, the time-frequency domain transformation based on Mixed-Radix FFT is relatively important, because the time-frequency domain transformation and window function are more related to signal processing. The influence on the detection of network attacks may be more subtle, in addition, the time-frequency domain transformation based on Mixed-Radix FFT has great advantages in capturing the local characteristics of the signal, which helps to improve the performance of the model.