The Use of Satellite Data to Determine the Changes of Hydrodynamic Parameters in the Gulf of Gdańsk via EcoFish Model

Using mathematical models alone to describe the changes in the parameters characterizing the analyzed reservoir may be insufficient due to the complexity of ocean circulation. One of the ways to improve the accuracy of models is to use data assimilation based on remote sensing methods. In this study, we tested the EcoFish numerical model that was developed for the Gulf of Gdańsk area, under the FindFish Knowledge Transfer Platform. In order to improve the model results and map local phenomena occurring in the studied water, which would be difficult to simulate using only mathematical equations, EcoFish was extended with a satellite data assimilation module that assimilates the sea surface temperature data from a medium-resolution imaging spectroradiometer and an advanced ultrahigh-resolution radiometer. EcoFish was then statistically validated, which resulted in high correlations for water temperature and salinity as well as low errors in comparison with in situ experimental data.


Introduction
The Baltic Sea is a shallow inland sea connected to the North Sea by narrow (4-28 km wide) straits. The topography of these straits, and in particular their low depth (around 5-6 m in most shallow points), impedes the free exchange of waters between the Baltic and the North Sea, causing significant water exchanges between these seas only during large infusions [1,2]. The Gulf of Gdańsk is the southern part of the Baltic Sea and is less affected by infusions, but it is exposed to influences from the land. The highly urbanized and industrialized coast has a huge impact on the environmental conditions of the Gulf in addition to the waters coming from the Vistula River [3,4]. An additional obstacle in the exchange of the Bay's waters with the open sea is the Hel Peninsula. It serves as a natural land barrier, marking the border between the Puck Bay and the Gdańsk Basin.
The morphometric conditions of the Gulf of Gdańsk favor the heterogeneity in salinity. Visible differences can be observed between the shallow coastal area hydrologically belonging to the surface layer of the Baltic Sea and the remaining deeper zone of the gulf [5,6]. However, smaller differences exist between the deep water region of the Gulf of Gdańsk and the open sea, where a typical for the Baltic Sea, layered water structure can be observed [7].
Comprehensive understanding and description of processes occurring in water is possible through the combined use of mathematical models, modern in situ research, and observational techniques in the form of satellite remote sensing [8][9][10][11].
In order to increase the intensity of knowledge transfer and the use of scientific potential by fishermen, and consequently contribute to the sustainable development of fisheries while increasing the protection of the Gulf's ecosystem, a three-dimensional, prognostic ecohydrodynamic model named EcoFish was built. The EcoFish model is configured for the Gulf of Gdańsk area and being developed under the project "FindFish Knowledge Transfer Platform-Numerical Forecasting System for the Marine Environment of the Gulf of Gdańsk for Fisheries" [12].The main goal of the project is to create a platform from which users (in particular fishermen and scientists) will be able to obtain knowledge and information about the physical and biological conditions of the Gulf of Gdańsk. The EcoFish model will be an essential tool in achieving this goal.
The most important modeled variables (especially water temperature, salinity, and dissolved oxygen concentration) that determine the state of the Gulf of Gdańsk ecosystem will also serve as input data for the Fish Module. Output information from this module, together with the data on fish preferences and expert knowledge, will allow the creation of maps of the most favorable environmental conditions for the occurrence of industrial pelagic fish in the Gulf of Gdańsk region, i.e., herring, sprat, and flounder.
The implementation of the satellite data assimilation module for sea surface temperature (SST), and also for chlorophyll a later in the biological part, in the EcoFish model will ensure that the obtained model results will be even more accurate (close to reality) and the model will correctly reflect the state of the Gulf of Gdańsk environment. In this study, the results from the EcoFish model were statistically validated.

Study Area
The domain of the EcoFish model includes the extended Gulf of Gdańsk (Figure 1), which is the southern part of the Gdańsk Deep area, located in the Gotland Basin. A straight line connecting Cape Rozewie with Cape Taran delimits the proper Gulf of Gdańsk. This line crosses the deepest parts of the Gdańsk Deep, with a maximum depth of 118 m. Along the coastal zone, there is a wide strip of shallows widening to the west of the mouth of the Vistula River. The slope of the bottom in the coastal zone is varied. The greatest decline occurs at the headland of the Hel Peninsula, where the bottom rapidly drops to a depth of 70 m [13]. The Bay of Puck is a shallow area of the Gulf of Gdańsk located in its western part [14,15]. Due to its geographical location and special hydrological conditions, the Bay  The EcoFish model consists of two active and two passive components. The active prognostic models are Parallel Ocean Program (POP) and Community Ice CodE (CICE). The passive components are responsible for providing atmospheric data fields and fresh water from catchment area. The main part of this system is the ocean model in which horizontal mixing is represented by biharmonic operator that is implemented by applying the Laplacian operator twice. It is performed in both viscosity and diffusion schemes; however, the mixing coefficients are different and equal to −3 × 10 14 and 0.3 × 10 14 , respectively (note: the sign was omitted). Vertical mixing used by authors in EcoFish system is called k-profile parametrization (KPP) [17]. Shear instability is parameterized in terms of the gradient Richardson number, while diffusivity (active tracers) and viscosity (momentum) are parameterized as diapycnal. Each active component has its own time step. The CICE time step is 10 min and it is equal to the coupling time step. Ocean model has typically two mods and time step is divided into two parts. The baroclinic part makes one step in 60 s. Because the model has linear free surface formulation, it is not needed to make substeps for the barotropic part. CESM is intended for use in global applications; thus, it was adapted for this purpose. Barotropic equation has been modified for possibility of assimilation of sea level at the boundaries [18]. Additionally, at the boundary area, salinity has been forced to have values from external model using restoring. The restoring time scale is not constant and depends on distance from the boundaries. Exactly at the boundaries, the restoring time scale is in the range of days and it increases to infinity toward the center of the domain at a distance of about 20 km (40 model cells). Similar models were already presented in [18,19].

Open Boundary
The results from the EcoFish model analyzed in this paper were spatially limited to the Gulf of Gdańsk area. However, the entire domain of the model is slightly larger and borders on the west and north with the open Baltic. Therefore, it is necessary to provide the model with boundary conditions. Apart from water temperature and salinity, it is necessary to prepare the data of sea surface height and barotropic components of sea currents. The data to the model boundary are provided by a 3D CEMBS [20,21] model with a horizontal resolution of 2 km. Since the exact determination of the data range on the border required trials with different settings, it was decided to prepare data on the entire FindFish domain, not only on its border. This allows the use of the same data regardless of the finally adopted settings. The fact that the results of the 3D CEMBS model are used as the source of boundary conditions means that 3D CEMBS's calculations must be performed prior to EcoFish.

Atmosphere Forcing
At the water-atmosphere border, the EcoFish model is fed with meteorological forcing. These data come from the UM (Unified Model) developed at the Interdisciplinary Modeling Center of the University of Warsaw (ICM UW). Parameters directly used as inputs are as follow: • 10 m wind speed, • 2 m air temperature, • specific humidity, • sea level pressure, • precipitation (convective and large-scale), • downward shortwave and longwave radiation Missing parameters are calculated by the atmospheric data module, which is an integral part of the EcoFish model. Air density and scattered and direct shortwave radiation in the visible and near-infrared range are determined that way.

River Discharge
Due to the fact that modeling of surface runoff requires the use of a hydrological model, The Soil & Water Assessment Tool (SWAT) software was used [22][23][24]. SWAT is a small watershed to river basin-scale model used to simulate the quality and quantity of surface and ground water and predict the environmental impact of land use, land management practices, and climate change. This adaptation of SWAT was developed under the WaterPUCK-Integrated Information and Prediction Service project [11]. Meteorological data (precipitation, wind, temperature, and atmospheric pressure) comprise a key element of any water balance model. The SWAT hydrological model is based on real-time observations (local meteorological station) as well as shorter weather forecasts (ICM UW website). The conversion of rainfall data into surface runoff is accomplished using the SCS (Soil Conservation Service) Curve Number procedure through the cumulative runoff volume and concentration time. SWAT covers sedimentation in surface and groundwater, and the transport model also includes snow cover. Additionally, the transport of nutrients and pesticides was taken into account for use at a later stage related to the launch of the biochemical part of the EcoFish model.
The SWAT model was created for six rivers (numbered 8-13) from the described domain (Table 2). For the remaining seven rivers (numbered 1-7), information on the runoff volume was taken from the Hydrological Predictions for the Environment (HYPE) model. It is a physics-based semidispersed catchment model that simulates the flow of water and substances as they travel from precipitation to discharge into the sea. We used historical time series from 1980-2010 for the geographic domain of Europe available as daily averages. Spatial resolution is determined by dividing the land area into catchments for which the HYPE data represent mean values at the estuary. The volumes for the years 2014 to 2020 have been calculated as the multiyear average over the available 30-year period.

Simulation Run
The EcoFish model was validated using the results of a seven-year simulation from 1 January 2014 to 31 December 2020, preceded by a two-year model spin-up. This run was calculated using the above-described model configuration with satellite data assimilation for surface temperature. Results were recorded four times a day as 6-h averages. The simulation results and model validation are presented in Section 4.
Additionally, in order to verify the correctness of the assimilation module itself, the same simulation was carried out without the assimilation of satellite data for the surface temperature. Both runs were then compared against satellite data.

Data Sets Used for Evaluation
Two in situ databases were used for the statistical analysis of the EcoFish model. Detailed descriptions are provided below.

ICES
The main in situ database that was used to validate the water temperature and salinity in the EcoFish model is the hydrochemical database provided online by the International Council for the Exploration of the Sea (ICES) via the https://ocean.ices.dk/HydChem/, accessed on 1 Febuary 2021 website. 17,902 measurements from 2014 to 2019 were used for comparison. Data for 2020 was not yet available in the database at the day of this validation.
Most of the data was from 2014 (6130 measurements) and 2015 (6848 measurements). From 2016 to 2019, 1055 to 1444 measurements per year were available. ICES measurements covered almost the entire domain ( Figure 2), except its eastern part, which covers areas along the Vistula Spit to the easterly coast of Gulf of Gdańsk. The region with the densest measurement grid was the strip along the 19 • E meridian.
ICES measurements were relatively uniformly distributed throughout the water column. In the surface layer from 0 to 5 m, 1082 measurements were available, in the layer from 5 to 30 m-5429, in the layer from 30 to 80 m-8987 and 2404 measurements for depths greater than 80 m.

Fishing Cruises
We had an additional database that was used to validate the model. It was a set of measurements taken under the FindFish project during fishing cruises. Valeport MIDAS CTD model 500 and GPS 73 Worldwide by GARMIN were used for this purpose. Recording of physicochemical data was carried out during fishing with towed and set gears. The following parameters were saved in the device memory: Data from the MIDAS CTD instrument were collected in the area north of the Vistula River mouth and in the strip in the open sea, parallel to the Hel Peninsula ( Figure 2). The measurements have a high temporal and spatial resolution. Therefore, they were averaged to match the model mesh. After this operation, we had 15,312 records. The highest density of measurements (8533 values) were from 30 to 60 m deep layer. It was closely related to the optimal fishing depth. The first data was collected on 22 May 2018 during pelagic fishing. To date, 422 hauls/releases have been made, of which 306 were harvested with towed gear and 116 with set gear. Most of the data was collected in 2020 (7293 measurements) and 2019 (5954 measurements). The remaining 2065 measurements were taken in 2018. The months with the highest amount of data collected were April (3150) and March (2395). Measurements were taken from the decks of 10 vessels, and the schedule of their cruises was determined on an ongoing basis in relation to the frequency of occurrence of fish species being the target of the catch. In addition, meteorological data were also recorded, including air temperature, wind strength and direction, cloud cover, and sea state, which were collected in the form of questionnaires completed by skippers.

Satellite Data Acquisition and Processing Module
The satellite sea surface temperature data (SST) used along with EcoFish model come from the SatBałtyk project database [8,9] and are based on measurements from a mediumresolution imaging spectroradiometer (MODIS AQUA) and an advanced ultrahigh-resolution radiometer (AVHRR). These data are calibrated to local conditions in the Baltic region and subjected to atmospheric correction and filtration. Raw maps are linked in space, geometrically corrected for changes in satellite position, and radiometrically corrected in case of numerical errors in data transmission. Data in the SatBałtyk system have a horizontal resolution of 1 km and cover the entire area of the Baltic Sea. Daily average values obtained from the combination of all satellite images available on a given day are used for assimilation. Assimilation takes place during the assimilation window with the peak of the window set to 12:00. This is taken into the account during the process of combining data using weighted averaging. The data acquisition management module automatically detects the presence of new data, downloads new files to a local storage, processes the data, and saves the result together with control files and metadata. Data processing consists of their interpolation onto the 575 m grid of EcoFish model and saving in the file format compliant with the model requirements. Thanks to the aforementioned control files, the module managing the entire system can monitor the status of satellite data on an ongoing basis and, if available, start assimilation.

Satellite Data Assimilation Module
The assimilation module is based on modified and extended CESM model [25,26] components that are an integral part of the EcoFish model. This allowed for smooth introduction of satellite information to computed data with each time step. In addition, it allowed the reuse of a number of settings and functionalities already available in the CESM, e.g., parameterization of the length of the assimilation window currently set to 24 h, the frequency of assimilation (each time step) and the modules for reading and processing data by the model. The method of assimilation selected in the EcoFish model takes as input the values of a given variable V mod obtained from the calculations of the model and satellite measurements V sat . Moreover, the method accepts a number of parameters that control its behavior and describe the data source. The most important of them are: • data_type-allows you to specify the frequency with which data for assimilation appears, e. With each calculation step, each assimilation module checks, on the basis of the data_type and data_inc parameters, whether new assimilation data should appear in a given step. If so, it loads the next file with assimilated data V sat . For the sake of example, let us assume that data appears every 24 h. Of course, these data should not be entered into the model at once, in one time step, as this would disrupt the continuity of equations and the equilibrium of the system. Hence, using the interp_freq parameter, one can choose that the data should be entered gradually, e.g., every 0.5 h. Using the restore_tau parameter, one can specify that the model should reach the assimilated values after a period of 12 h. Having these parameters, the module divides the current difference between V sat and V mod into the number of steps resulting from interp_freq that fall within the restore_tau period, i.e., in this case by 24.
Here, dV step is a partial increment introduced to the model in a given assimilation step (0.5 h in this example). The value of a given model variable depends on many factors, e.g., transport, radiation, biological processes, etc. Therefore, the difference dV calculated at the beginning between the model data and the measurement data must be constantly corrected to obtain the expected value at the end. Hence, it is updated at the frequency specified by the interp_inc parameter. The resulting value of the assimilated variable is calculated by adding the calculated increment to the model result: The figures below present the application of the satellite data assimilation module in the EcoFish model for two selected days from the simulation period. On the left side, one can see the satellite image for the surface temperature, in the middle, the model result with visible effect of assimilation, and on the right side, the results from the model version without the assimilation. The first scene ( Figure 3) was taken on 28 April 2019. On that day, due to the heavy cloud cover, the satellite image was severely restricted and provided information only about a small area within the domain ( Figure 3a). As a result of the assimilation module, the surface temperature in the vicinity of the Hel Peninsula and the eastern coast of the Gulf decreased, which can be seen in the Figure 3b. Figure 3c shows the surface temperature from the version of the model without the assimilation.

EcoFish Model Validation
The accuracy of the EcoFish model results, in the period from January 2014 to December 2020, for water temperature (Section 4.1.1) and salinity (Section 4.1.2) was checked by comparing them with available in situ observations from the ICES database and measurements made with the MIDAS CTD instrument during fishing cruises. This was conducted for the two versions of the EcoFish model-one with the assimilation module enabled (EcoFish+A) and one without assimilation (EcoFish−A). The most important statistical quantities were determined, such as Pearson's correlation coefficients (r), root mean square errors (RMSE), standard deviations (STD), and biases (in relation to means), which illustrate the tendency of the model to systematically overstate or underestimate the result.

Water Temperature
At first, we wanted to verify the impact of the assimilation of satellite SST module itself. Therefore, the surface temperature was validated for 2018 using two separate versions of the EcoFish model-one without (EcoFish−A) and one with satellite data assimilation (EcoFish+A). All other model parameters were identical. Results from both runs were then compared against assimilated satellite data. Table 3 provides a statistical summary of the temperature obtained from all three sources to give a better overview of its characteristics. Table 4 contains comparison of both model versions against satellite data used for assimilation. As one can observe, the differences between model and the satellite are larger in case of the model without assimilation (Table 4). Pearson's correlation coefficient (r) increased from 0.95 to 0.98. The root mean square error (RMSE) decreased by almost 1 • C from 2.31 • C to 1.45 • C, and the tendency of the EcoFish model to systematically underestimate the results by −1.12 • C for the nonassimilated version decreased to −0.56 • C for the assimilated version. This confirms that the assimilation mechanism itself was designed and implemented properly and it yields significant changes in obtained results. It is worth noticing that the assimilation does not result in 100% agreement between satellite and model data. This can be contributed to several factors. Most importantly, in order to save disk space, results saved by the model are averaged over the period of 6 h. Satellite data are also averaged over the assimilation period, but the assimilation is designed to have the peak alignment of the data precisely at 12:00. Apart from that, the model surface layer is much thicker than the surface layer measured by satellites. Because of that, it would be wrong to exactly replicate satellite results in the model. The assimilation module is parameterized to leave some level of freedom between model and satellite data. Last, one must take into the account that even though satellite data are assimilated with every time step, each water cell in the model is subjected to vertical and horizontal currents that dissipate the assimilated information. Given the above, obtained results are satisfactory and prove that the assimilation module works as expected.
In order to additionally emphasize the impact of assimilation on the improvement of the model results, a separate comparison of the temperature in the surface layer itself (from 0 to 5 m) with the data from the ICES database was made. For the surface, within the model domain, 1082 observations from the ICES database were available. The result of this analysis is presented in the Taylor diagram ( Figure 5) in the form of normalized statistical coefficients. Absolute statistical values are presented in tabular form (Table 5).  The ability of the EcoFish model (with an active SST satellite data assimilation module and without) to correctly project the real environment conditions was verified by comparing the model results for temperature with all available observations at all depths. The result of this validation is presented in the Taylor diagram [27] (Figure 6) in the form of normalized statistical coefficients. Absolute statistical values are presented in tabular form ( Table 6).  The Pearson correlation coefficient (r) calculated between the assimilated model (EcoFish+A) and the data from the ICES database was 0.94. When compared with data from fishing cruises using the MIDAS CTD instrument, it took the value of 0.87. The decline may be related to the varying density of in situ data from surface to bottom. ICES data was relatively homogeneously distributed in the water column, while the cruise data most often occurred at fishing depths, i.e., 30 to 60 m. Due to the use of assimilation, the model better reflects the temperature closer to the surface, as indicated by high correlations that slightly decrease with depth. Table 6 shows also that when comparing to ICES data, a nonassimilated version of the model (EcoFish−A) produced slightly better statistics than the assimilated version (EcoFish+A). This might be related to the database itself and should not be taken as a sign that assimilation worsened the results. It is opposite in the case of MIDAS CTD data in which statistics for EcoFish+A are better than for the version without assimilation. Summarizing, high correlations and low RMSEs were recieved from both versions of the model.
The mean square error (RMSE) for EcoFish+A is 1.33 • C when compared with ICES data and 1.83 • C when compared with cruise data. The model has a similar bias for both in situ databases of −0.36 • C and −0.34 • C, respectively. This means that the model tends to systematically underestimate the results slightly. The standard deviation for both databases used during the validation is similar and amounts to 3.66 • C for ICES data and 3.57 • C for data from the MIDAS CTD instrument (Table 6).
When analyzing the vertical profile (Figure 7), created while taking all observations from the ICES database into account, compared to the corresponding values from the EcoFish model (with and without assimilation), a high correlation can be seen, lasting from the surface up to about the 13th level of the model (to a depth of 80 m). Below the 13th level, the data are less correlated and the EcoFish model has a tendency to slightly underestimate the results. The reason for this is that the POP model has KPP implemented for surface layer only. Since it is designed for global issues, the vertical mixing scheme has been modified for better representation of surface layer in the Baltic Sea. Small changes of the interior shear mixing, suggested by Durski et al. [28], have been introduced. For better reproduction of the bottom layer, the dependence of the quadratic drag formula on thickness of the lowest model cells was applied (the logarithmic profile of roughness height). There is no turbulent closure though that will cover the bottom layer, which results in the underestimation of temperature ( Figure 7) and salinity (Figure 8) in the bottom layers.
The most important thing is that the model correctly reflects the temperature drop in the thermocline layer (on average from the 3rd to the 9th depth level). Not only is a high correlation visible, but close and overlapping ranges of the double standard deviations are also visible. Only below the 7th level (below the depth of 35 m) are the model data characterized by a smaller standard deviation than the observational data from the ICES database.
Therefore, it can be concluded that, despite the slight discrepancies between the model results and the observations (mostly at greater depths), the vertical mixing in the EcoFish model has been correctly implemented, and the model itself reflects real conditions well.

Salinity
Another physical variable derived from the EcoFish model that we validated is salinity. Salinity is a good parameter to check if the model is capable of handling water masses, as it does not undergo any transformations (gains and losses) in the marine environment. Salinity results from both versions (with and without assimilation) of the model were compared with available in situ observations from the ICES database using the same locations as for temperature ( Figure 2). The result of this comparison is presented in the Taylor diagram ( Figure 6) in the form of normalized statistical coefficients, as well as in the form of absolute values in table (Table 7). The Pearson correlation coefficient (r) calculated between the model data (EcoFish+A) and the observations from the ICES database took the value of 0.94, which is better than in the case of the nonassimilated version of the model. The root mean square error (RMSE) was 0.8 PSU; this is a satisfactory result, while having relatively large standard deviation of 1.27 PSU. The model has a low negative mean bias of −0.01 PSU, which could signify that it responds well to changes in salinity. However, this is largely due to the fact that at low and medium depths (from the surface to about 14th level), the model results for salinity are slightly higher than the ICES observational values, followed by a trend that reverses. At higher depths, the model begins to underestimate salinity from about 0.5 up to 2.0 PSU. This was demonstrated on a vertical salinity profile (Figure 8) created using all observations from the ICES database.
When analyzing the vertical profile, several characteristic zones in the water column can be seen. On the surface, we observe an increased standard deviation, both for model data and observations. It is the result of mixing fresh waters from river runoff with sea waters, causing increased salinity dynamics on the surface. Then, between 10 and 40 m deep (up to 60 m depending on the location), the isohaline layer stretches, which is noticeable both at a constant average salinity of about 7-8 PSU and the size of the standard deviation, which in this layer is the smallest throughout the entire water column. Below 35 m, the STD begins to rise gradually, reaching about 2 PSU by 70 m and remaining at this value to the bottom. Below the isohaline layer, there is a transition layer with a clearly delineated halocline, especially for the curve determined using ICES observations. Average salinity starts to increase systematically from about 55 m (11th level) and stabilizes only at a depth of 100 m (20th level). The curve determined using the model data also shows the presence of a halocline, but it is not so clearly marked. The model results do not have the same dynamics as in situ data, which is especially visible in a much smaller STD than in the case of observations. In the water layer from the 12th level to the bottom, the salinity obtained from the EcoFish model increases from about 9 PSU to 12 PSU, while the increase in salinity for the experimental data is greater and goes from 8 PSU to 13 PSU (Figure 8). The highest salinity values seen in the experimental data can be related to the infusion effect, while decreases can occur during periods of long-term stagnation.

EcoFish Model Simulation Results
In this section, we present the average monthly temperatures, salinity, sea surface height, and currents in the surface layer for EcoFish model with an active SST satellite data assimilation module. The averages reflect the simulation period from 2014 to 2020. The results for individual months have been transferred to the Appendix A section (Figures A1-A8 and Tables A1-A23). Additionally, figures with the differences between the monthly averages of parameters from the models with (EcoFish+A) and without (EcoFish−A) assimilation are presented.

Water Temperature
The model domain is characterized by a strong seasonal variability of the surface temperature ( Figure 9) for each year from the seven-year simulation (2014-2020). The greatest dynamics occurs in the southern part of the domain, which includes the southern part of the Gulf of Gdańsk and the Bay of Puck, which has the greatest number of low-depth coastal areas that react quickly to atmospheric forcing. The remaining regions are characterized by relatively lower seasonal variability of the surface temperature.
The mean water surface temperature for the entire analyzed model domain was 10.43 • C. The months with the lowest average temperature are generally February and March, and the coldest month in the simulation period was March 2018 with an average temperature of 2.05 • C. The warmest months are usually the summer ones, i.e., July and August with the record-breaking August 2018, in which the average water temperature in the surface layer reached 21.23 • C ( Figures A1 and A2 and Table A1).
Temperature extremes for single model cells were most common in shallow coastal regions. The most characteristic is the shallow inner Bay of Puck called Puck Lagoon (Figure 1), which is not only separated from the open sea by the Hel Peninsula, but in its eastern part, there is a characteristic bathymetric anomaly in the form of shallow water (marked with curved dashed line), which additionally limits the exchange of water masses with the outer Bay of Puck (and whole Gulf), influencing to local temperature extremes. The lowest temperatures were in January and February 2014, falling to −0.43 • C (Table A2), and the highest were in July 2014 and 2018, when the temperature exceeded 28 • C ( Table A3).
The lowest standard deviations ranging from 0.34 • C to 0.81 • C (0.53 • C on average) were obtained by the model for February and March. It is the period before spring when the surface layer is cooled down after winter, and both air temperature and sunlight do not yet reach high values to cause significant local changes. The largest deviations from the mean surface temperature appear in May and June (1.22 • C to 3.11 • C) with values exceeding 3 • C in May 2017 and 2018 (Table A4). The picture below ( Figure 10) presents the average monthly surface temperature for the period 2014 to 2020. It can be seen that there are four characteristic periods of temperature variation in the domain. The longest, five-month cold period, lasting from December to April, when the average surface temperature is low and relatively stable, ranging only from about 3 • C to 7 • C. Then, a four-month warm period lasting from June to September with average temperatures ranging from about 15 • C to 19 • C. There are also two transitional periods. The first transition (uptrend) period is May, when the temperature rises sharply from 5 • C in April to 15 • C recorded in June. The second transition period (downward) is October and November, when the water cools down quickly from an average of 17 • C in September to 6 • C in December (Table A1).
The reservoir that responds most quickly to external factors is the aforementioned Puck Lagoon. Thanks to its specific location and bottom topography, it has the greatest variability of the surface temperature and the minimum and maximum temperature achieved. In winter, local ice cover is observed, while in summer, it is exposed to toxic algae blooms stimulated by high temperatures and the deposition of inorganic phosphate from rivers. Looking at Figure 11, that shows the difference in the monthly average SST values between the model with (EcoFish+A) and without assimilation (EcoFish−A), it can be seen that the assimilation influences the increase of the temperature in the warm season (from May to August) and the decrease in the cold one (from October to January). It is related to the vertical resolution of the model, and more precisely to the thickness of the first (surface) layer, which is 5 m.
It is more than the actual surface layer for which the SST is measured from satellite instruments. A thicker layer responds slower to atmospheric forcing, as it has a higher heat capacity. Assimilation of the satellite data into the model adjusts the surface layer temperature to the atmospheric conditions faster; hence, it responds faster to the changes.

Salinity
During the seven-year-long simulation, the average monthly salinity in the surface layer was between 7.31 PSU and 7.76 PSU, which gives an average of 7.47 PSU for the entire time period (Figures 12 and 13). The annual course of salinity in the study area is usually established assuming lower values in the warm/summer months (minimum for April 2014) and higher values in the cold/winter season (maximum in February 2014) ( Figures A3 and A4 and Table A9).  Figure 14 shows that there are differences in surface salinity distribution when comparing models with and without SST assimilation. The model with assimilation (EcoFish+A) gives almost 1 PSU higher salinity in the south area of Gulf of Gdańsk from May to August. The reason for that is the increase in water circulation (Section 4.2.4 ) leading to more salty waters located in the deeper bottom layers being transported to the surface as a result of mass conservation law.  (Table A12), rarely exceeding 0.3 PSU. In the southern part of the domain, which includes the mouth of the Vistula River, salinity can fluctuate from about 2 PSU to even 8 PSU in the summer months ( Figure 13). The lowest salinity recorded in the model for the simulation period, amounting to 1.51 PSU, was recorded in March 2014 (Table A10). It was the result of the spring rise and runoff from the Vistula River. The highest of 8.70 PSU was obtained for January 2018 at open sea (Table A11). The dynamics of changes in salinity along the surface layer is influenced by a number of factors. Among others, it is the amount of river runoff, seasonal changes in the thermal structure of water or changing meteorological conditions. However, in the bottom layer, the salinity distribution over the year seems to be relatively uniform, strongly related to the bathymetry ( Figure 15). The highest salinity at the bottom (12.66 PSU on average) occurs in areas of great depth, in particular in the Gdańsk Deep. Possible fluctuations there, are no longer the result of cyclical processes with a seasonal frequency, but rather, of irregular events such as sea inflow. Only in 3 out of 84 considered months, the average monthly salinity at the bottom exceeded 13 PSU (Table A15)

Sea Surface Height
The lowest values of sea surface height determined in the model were −29.58 cm in November 2015 and −31.09 cm in October 2016 (Table A18). The highest, 60.19 cm and 77.21 cm, were obtained for March 2020 and January 2015, respectively (Table A19).
The average sea surface height was 1.50 cm ( Figure 16 and Table A17). Individual monthly averages are characterized by higher standard deviations (between 4.96 cm and 6.47 cm) in the months from October to February (Table A20). In the remaining months, the average standard deviations are usually smaller and reach values of about 3 cm (Figure 17).  The process of SST assimilation has minor but noticeable influence on sea level ( Figure 18). The increase of salinity in the coastal area ( Figure 14) results in higher density, which has direct impact on pressure distribution. Therefore, an increase in salinity causes a decrease in sea level, keeping the hydrostatic balance of the system.

Currents
The spatial distribution of sea currents inside the domain is much more characteristic and repeatable than in the case of sea surface height. The average current velocity in the surface layer for 2014-2020 was 6.73 cm·s −1 with an average standard deviation of 5.23 cm·s −1 (Figure 19 and Tables A21 and A23). The strongest currents were obtained for December 2016 and January 2015 and reached the speed of 104.45 cm·s −1 and 120.09 cm·s −1 (Table A22), respectively.  This is also an area of frequent coastal upwelling, which causes cold water masses to rise from the bottom to the surface. The process responsible for this phenomenon is the Ekman transport, associated with the persistence of the southeastern wind running along the Hel Peninsula. For example, the August 2015 current map ( Figure A8) shows a strong northwesterly current along the peninsula due to winds blowing in August, resulting in upwelling. The surface temperature map from 25 August 2015 ( Figure 21) shows a large horizontal gradient. Gradients observed in this region are often reaching up to 5 • C·km −1 [29]. Temperature has a direct effect on the density of water, and the density becomes lower as the temperature increases. As the result of such modification, the same wind stress gives stronger currents because of the second law of dynamics, which is presented as the difference between assimilated and nonassimilated sea currents ( Figure 22). Consequently, bigger sea currents increase water circulation in the coastal areas of Gulf of Gdańsk. Images of current roses in the surface layer for selected three characteristic regions within the domain (Figure 2) are presented below.
The VR region (Vistula River) covers the coastal and shallow, southern part of the Gulf of Gdańsk within the mouth of the Vistula River. Monthly averages of surface currents in this area rarely exceed 16 cm·s −1 (Figure 23). In this region, eastern currents constitute the dominant part of directions. In the months from December to April, they have 39% and a greater share of all directions. Such current direction causes that the water flowing out of the Vistula has a difficult outflow and is distributed to the east along the shore of the Gulf. The long-term presence of such currents limits the spread of river waters and reduces the zone of fresh and sea water mixing. In the GD region (Gdańsk Deep), which covers the deepwater area of the domain, located directly above the Gdańsk Deep, both the distribution of average velocities and directions is much more homogeneous than in the case of other regions ( Figure 24). This is due to the high variability of wind directions and velocities over this area. Moreover, due to the great depths occurring here, the bathymetry does not have such a strong influence on the distribution of currents as in the case of areas close to the coast. From November to February, statistically more often there are currents heading in the eastern, northeastern and southeastern directions (about 60% of cases), but in the remaining months, the situation is less diversified. For example, in the summer months (from June to September), the southeast, south and southwest are dominant directions. They account for 48.6% in June, 57% in July, 58.3% in August, and 53.6% in September. In this region, the average monthly current speeds are comparable to those obtained for the VR region. Currents with speeds in the range 4-16 cm·s −1 often constitute even 60-70% of all velocities here. Moreover, in each month, a small share of current velocities exceeding 16 cm·s −1 , and sometimes even 24 cm·s −1 , can be found. In the model results for HP (Hel Peninsula) region, we observe the highest average monthly current velocities in the entire domain ( Figure 25). Average values exceeding 24 cm·s −1 appear here every month and constitute from 3 to 22% of all speed ranges. The directions of surface currents are specific in this region. In each month, the dominant directions are northwest (prevailing in the warm months, from May to November) and southeast (in other months). Together, they constitute over 60% of all currents in each month. Such structure of surface currents enables rapid movement of water masses along the Hel Peninsula, mixing together the Gulf of Gdańsk and Baltic Proper waters. The surface current distribution in the vicinity of Hel Peninsula is mostly induced by dominating westerly winds [30], but also depends on the large-scale internal water cycle in the Baltic Sea [31].
The current velocities in the bottom layer along the VR and HP regions rarely exceed 4 cm·s −1 . Only in GD they have a higher share of over 10%, especially from October to February. The rose of currents for the Gdańsk Deep indicates the existence of a dominant northern bottom current ( Figure 26). This suggests that the waters at the bottom are most often pushed toward the Gotland Basin. On the monthly average maps ( Figure A9), from December to March, however, the dominance of the southern current can be observed. It moves the water masses in the shallow-water direction of the southern part of the Gulf of Gdańsk.

Discussion
This paper presents the hydrodynamic part of the three-dimensional numerical model EcoFish. The model has been used to simulate the hydrodynamics of the Gulf of Gdańsk. EcoFish has a satellite data assimilation module for SST, which uses data from the SatBałtyk project database, consisting of photos from a medium-resolution imaging spectroradiometer (MODIS AQUA) and an advanced ultrahigh-resolution radiometer (AVHRR). The task of this module is to assimilate the available surface temperature satellite data into the model domain in order to improve the simulation results, allowing for better determination of the dynamics of changes in physical parameters. The EcoFish model covers the South Baltic Sea and, more precisely, the entire Gulf of Gdańsk together with the Bay of Puck. The model domain is connected via an open boundary with the Baltic Sea from the west and the north.
The article presents statistical validation of the EcoFish model, which allowed to verify the correctness of the results obtained from it in terms of seasonal and spatial variability of the simulated water temperature and salinity. For this purpose, the available in situ observations from the ICES databases were used, along with the database created during fishing cruises carried out under the tasks of the FindFish Knowledge Transfer Platform. For the entire analyzed simulation period from January 2014 to December 2020, basic statistical values were determined, such as root mean square errors (RMSE), standard deviations (STD), and Pearson correlation coefficients (Section 4.1).
The validation showed that the EcoFish model results for water temperature were consistent with in situ observations. To confirm this, two experimental databases were used. Almost 18,000 measurements were available in the ICES database, distributed relatively evenly throughout the domain, except for the shallow region along the coast where no monitoring was carried out, or at least data from this region were not publicly available ( Figure 2). The correlation of the EcoFish model with these data (ICES) was 0.94 with the root mean square error (RMSE) of 1.33 • C. As a result of comparing the modeled temperature against the data from the database created during fishing cruises, a correlation coefficient of 0.87 were calculated. This is a satisfactory result, taking the strong concentration of cruise data in the belt from the mouth of the Vistula River in the northwest direction into account. Thus, the data come both from the area where there is mixing of river waters (from the Vistula River) with the waters of the Gulf and from the area where the strongest currents occur in the entire domain (the belt along the Hel Peninsula). It should also be noted that most of these data came from fishing depths, i.e., 30 to 60 m, where the impact of SST satellite data assimilation is no longer visible.
The correlation of the model results for salinity with the ICES data at the level of 0.94 and the low root mean square error of 0.8 PSU suggest that the model copes well with the transport of water masses. It also proves that the rivers in the model have been correctly implemented and that the outgoing freshwater is correctly mixed with the saltwater of the Gulf and distributed by currents in its area. Additionally, when analyzing the vertical profile (Figure 8), both the isohaline layer and the formation of a halocline at lower depth levels can be seen, which proves that the model correctly reflects the dynamics of salinity changes in the water column.
When analyzing the seven-year simulation period of the EcoFish model (from January 2014 to December 2020), it can be observed that the temperature of the Gulf of Gdańsk waters is subject to strong seasonal changes and depends mainly on changes in air temperature and solar radiation. They are also largely influenced by convection processes and wind-induced mixing. The changes in the water temperature of the Gulf also show the influence of the Vistula River, whose waters increase the temperature in the Gulf in spring and summer and lower it in autumn. The lowest values of surface water temperature occur in January and remain at the level of about 0.1 • C (Table A2). On the other hand, the lowest average values of surface water temperature in the entire domain occur in February (Table A1). During this month, the surface waters of the entire reservoir are characterized by a similar temperature, and the differences do not exceed 2.5 • C. In the following months, the temperature of surface waters increases (the fastest in the coastal zone). The highest spatial differentiation is observed in the model results for May and June (differences amounting to about 7 • C). The highest average surface water temperatures occur in August (Table A1).
The location of the Gulf of Gdańsk and its specific bottom topography favor the occurrence of salinity diversification. Significant differences in its distribution occur between the shallow coastal area and the deeper part of the Gulf, which resembles the waters with a layered structure typical of the Baltic Sea (with the presence of a halocline and a thermocline). The shallow-water coastal zone of the Gulf of Gdańsk is influenced by fresh waters entering it from rivers and other types of surface runoff. The Vistula River has the greatest impact on changes in salinity, with huge volumes of fresh water flowing out (average flow exceeding 1000 m 3 /s), causing local salinity drops below 7 PSU. Its influence is also noticeable in the surface layer of the deep-water part of the Gulf of Gdańsk, mainly in the spring season, when, due to currents, river waters mix with sea waters and are carried into the Gulf.
When analyzing the distribution of currents in the studied domain, a characteristic area can be distinguished, stretching along the Hel Peninsula. The strongest surface currents occur there, often exceeding 20 cm·s −1 . Two directions dominate there, depending on the season. Northwestern currents are mainly observed in the model results for the summer months, pushing the water from the Gdańsk Basin toward the open sea and are accompanied by the formation of coastal upwelling and downwelling. In the remaining months, southeastern currents predominate in this region, carrying waters toward the inner Gulf of Gdańsk. The distribution of surface currents around the mouth of the Vistula is also peculiar, where the most common is the eastern current, which distributes the waters flowing out of the Vistula along the shore of the Gulf. Its long-term presence limits the spread of the Vistula waters and reduces the zone of mixing fresh water with sea water. HP and VR regions are two specific regions along the Polish coastline where coastal up-and downwelling events occurs. This is induced by several factors, which the most important are the dominating westerly winds [30], bathymetry and vicinity of the coastal formation. In the HP region, in addition, the large-scale circulation of the Baltic Sea with the surface current pushing waters along the Hel Peninsula into the Gulf of Gdańsk [31] is responsible for the currents distribution.
The factor that has the greatest impact on changes in the sea surface height is wind. Certain areas can be distinguished with the greatest variation in the SSH. These are coastal areas, in particular around the Hel Peninsula, the southern coast stretching from the Bay of Puck, along the Vistula Spit, as well as the eastern shore of the Gulf of Gdańsk (Figures A5 and A6).

Conclusions
In this paper, we present a numerical model of the Gulf of Gdańsk EcoFish model with an active module of satellite data assimilation for surface temperature. This version of the EcoFish model is being developed and used within the framework of the project "FindFish Knowledge Transfer Platform-Numerical Forecasting System for the Marine Environment of the Gulf of Gdańsk for Fisheries". EcoFish is the basic element of the platform that provides fishermen and scientists with the current and forecast hydrodynamic, chemical, and biological conditions of the Gulf of Gdańsk. It also produces forecasts determining the most favorable environmental conditions for the occurrence of industrial pelagic fish in the South Baltic region. The aim of this research and development project is to increase the intensity of knowledge transfer and the use of scientific potential by fishermen, and consequently contribute to the sustainable development of sea fisheries and increase the protection of the Gulf of Gdańsk ecosystem.
To verify the correctness of the EcoFish model, a statistical analysis was carried out by comparing the model results with the in situ observations for the simulation period from January 2014 to December 2020. Satisfactory results were obtained and the compliance of the model results for water temperature and salinity with the available observations was confirmed. Correct mapping of the physical conditions inside the domain allows the model to be used for further simulations with an active part of the ecosystem. To do this, it is required to have a model that correctly simulates the physical conditions of mixing in a water body and heat exchange, controlling the heating and cooling of water masses. This is of great importance for the simulation of biochemical factors and the primary production process that will be conducted by the biochemical part of the EcoFish model.
The decision to validate the model for water temperature and salinity resulted from the fact that these two parameters serve as input data for the Fish Module, in which the Habitat Suitability Index maps are determined based on the environmental preferences of fish. The final product of the project will be sharing the FindFish platform as a website that will provide all the results and forecasts in operational mode. Acknowledgments: Calculations were carried out at the Academic Computer Centre in Gdańsk. We are grateful to the anonymous reviewers for valuable comments on earlier versions of the manuscript.

Conflicts of Interest:
The authors declare no conflict of interest.

Abbreviations
The following abbreviations are used in this manuscript:

ICES
International