Defining a Standard Methodology to Obtain Optimum WRF Configuration for Operational Forecast : Application over the Port of Huelva ( Southern Spain )

In this contribution, we calibrate the meteorological model weather and research forecasting (WRF) for operational forecasting in the Port of Huelva managed by the Authority Port of Huelva. Meteorological forecasting will allow reducing the impact of the meteorological phenomena over weather sensitive activities in the region. Concretely, the meteorological modeling developed will be used to analyze meteorological hazard impacts and to improve the management of the local air quality. To achieve these goals, numerous sensitive analyses corresponding to different model options have been developed. These analyses consider different physical and dynamical options, the coupling of very high resolution physiographic database (topography and land uses), and data assimilation. Comparing experiments, results with observational measures provide us by the Spanish National Meteorology Agency (AEMET). During a representative period, the optimum WRF configuration for the region is obtained. Calibration has been focused on wind due to this is the main risk factor in the region. When the model is satisfactorily calibrated, WRF is evaluated using whole modeling years 2012 and 2013, working with very high horizontal resolution, up to 0.333 km of horizontal grid resolution. Results obtained from the evaluation indicate that the numerical weather prediction system developed has a confidence level of 70% for the temperature, 81% and 66% for the wind speed and wind direction respectively, and 90% for the relative humidity. Methodology designed defines the quality control assurance of high-accuracy forecasting services of Meteosim S.L.


Introduction
During the last years, the use of numerical weather prediction (NWP) models has increased significantly due to: higher accuracy of the models, easier access to community models, computational advances, use of atmospheric models by specialist from other disciplines as air quality or oceanographers, etc. [1].In the same way, the number of model users has grown exponentially in recent years, and more universities, research centers and companies use meteorological models to investigate and to offer applied services.And currently, there are a great amount of web pages providing meteorological forecasts based just on models.
Meteorological/NWP models have a wide range of options to set up: physical options, dynamical options, horizontal model resolution, number of vertical layers and density, domains architecture and nest down options, assimilation data, time-step, spin-up time, etc.It is a fundamental factor when configuring a model the selection of the parameterizations and options that are used [2] [3].And the best combination for one region is not necessarily applicable to another [4].Scientific community has drawn up guides and recommendations on the use of meteorological models [1] [5]- [7] that modelers may consider previously to simulate.In this sense, authors have defined and established an own standard methodology that could be applied in any region and allow to improve meteorological forecast for operational purposes.The methodology and the procedure define the quality assurance mechanism of services and products developed by Meteosim S.L. and it is based on experience acquired during the last 13 years and on the high expertise of the staff working in the company.
Furthermore, the reliability of the meteorological model decreases if not taken into account the particular characteristics of the area of its application.The limitation of the numerical models is mainly associated to the correct reproduction of the world that tries to simulate.It means that if the world that "feeds" the model is far from reality, the output of the model will also be.To solve this problem, meteorological models can be applied in a grid with very high resolution (horizontal grid spacing around 1 km), and physiographic high resolution data may be introduced in the model, allowing incorporating the physical characteristics of the region with greater accuracy and representativeness.In this sense, numerous studies have demonstrated how the incorporation of topography [8]- [12] and land uses [13]- [15] databases with higher resolution significantly improve the forecast performance.Topography and land cover influence on wind patterns, land surface atmosphere interactions, such as low level turbulence and transport, convection, precipitation, runoff, etc. [16]- [23].
In this paper, we focus our attention on improving meteorological operational forecasts of the Port of Huelva (Huelva), a region in Southern Spain.This work aims to investigate the optimum configuration of a meteorological model that allows to reduce the uncertainty and, therefore, to increase the confidence level of the forecasts.We have used WRF model (Section 2) to obtain the meteorological forecasts and we have defined a procedure to calibrate the model in a customized way for Huelva but it can be applied to any region.These forecasts will be used as an early warning system and will allow improving the risk management associated to daily activities, and particularly, to manage more efficiently the atmospheric pollution generated as consequence of aggregate handling and storage piles associated to the treatment of different solid materials, providing a better air quality for the region.
The study we present here is innovative for three reasons: 1) it defines a standard methodology which allows improving meteorological forecasts using numerical sensitive analysis and a pre-established procedure, and it could be reproduced in any region; 2) the study uses and configures a mesoscale meteorological model over the metropolitan area of Huelva, incorporating for first time in the region, new physiographic databases and very high resolution in a 0.333 km horizontal grid, which allow to consider the extreme complexity of the Port of Huelva managed by Authority Port of Huelva and its influence over meteorological patterns; and 3) finally, in comparison with experiences in other ports or maritime regions of around the world, the meteorological forecasting system developed increases the resolution and the accuracy [24]- [28].
A complete description of the methodology followed to conduct the optimum WRF configuration for operational forecast is presented in Section 2. A detailed analysis of the results obtained is presented in Section 3, and finally, conclusions are reported in Section 4.

Methodology
In the following sections, we show a more detailed description of the studied area features, the simulations domains periods analyzed as well as a more comprehensive explanation of the modeling approach and model evaluation.

Studied Area, Simulation Domains and Episodes Selected
The methodology defined to obtain the optimum WRF configuration has been applied over Huelva, in Southwestern Spain.Air quality levels achieved and the risk management in a complex harbor located very near of the population made remarkable the implementation of a very accuracy meteorological model in the zone.
Huelva had a population of 149,410 in 2010 and it is located along the Gulf of Cadiz coast in the mouth of the Odiel and Tinto Rivers.Huelva has a Subtropical-Mediterranean climate characterized by dry and hot summers and wet and mild winters.Temperatures higher than 40˚C are usually reproduced on summers.Highest wind velocity is reproduced during the sunset and episodes of severe gale (force 9 in the Beaufort scale) affects the region occasionally during the year.
In the city of Huelva and its metropolitan area coexist the core of the population of the province of Huelva, greenhouse zones, nature reserves (very near of Doñana Park) and one of the most important industrial poles in the South of Spain.The activity of the industrial sector is mainly characterized by the Port of Huelva, divided in two sectors: the inner port (near the city of Huelva) and the outer port, being this last the most important.Activity in the Port is associated with a high flow of loading and unloading operations and material handling.In this sense, the meteorology influence the atmospheric pollution generated by these processes and the activity, in itself, is conditioned by the meteorological conditions.
In • To validate the model we have considered the two full years, 2012 and 2013, and therefore, 730 WRF daily simulations have been conducted.• And to analyze the experiments of data assimilation we have considered an operational forecasting period, corresponding to the period compressed between 08/24/2015 and 10/24/2015.For data assimilation analysis a total of 300 WRF simulations have been conducted, corresponding to 60 days and 5 experiments.A total of 4150 WRF daily simulations have been realized.

Modeling Approach
The regional and mesoscale meteorological model used for the study has been the Weather Research and Forecasting-Advanced Research (WRF-ARW) version 3.7 [29], developed by the National Center of Atmospheric Research (NCAR).It is a universally used community mesoscale model and a state-of-the-art atmospheric modeling system that is applicable for meteorological research, climate scenarios and numerical weather prediction.WRF is a fully compressible and non-hydrostatic model with terrain-following hydrostatic pressure coordinate.Different options that WRF offers can be combined in many different ways.WRF has different parameterizations for microphysics, radiation (long and short wave), cumulus, surface layer, planetary boundary layer (PBL) and land surface as physical options.To obtain WRF highest accuracy, it is essential to carry out a sensitive analysis of these different options by numerical experiments.In the same way, the definition of the simulation domains, spin up, vertical resolution or nesting architecture determine the accuracy, and therefore uncertainty, of WRF results [30] [31].
In this research, we have configured WRF model for operational forecast and we have considered a past period to calibrate and define optimum options and validate using a representative past period.The initial and boundary conditions for the operational configuration over domain d01 were supplied by the National Centers for Environmental Prediction (NCEP) Global Forecast System (GFS) with horizontal resolution 0.25˚ and updated every 6 hours.Otherwise, to calibrate and validate the model a new parent domain has been defined (called d00), covering all Iberian Peninsula, south of France and northwestern Africa.The parent domain d00 has been defined to adjust the coupling between global and mesoscale model.In this case, for the sensitive analysis initial and boundary conditions have been supplied by the NCEP/NCAR Climate Forecast System Reanalysis v2 (CFSv2) [31] with 0.5˚ of spatial resolution and 6h of temporal sampling.CFS is a better atmospheric representation than GFS because incorporates a higher amount of observations and measurements.In both cases, operational configuration and to calibrate and validate the model, two-way nesting was used for the external domains (d00, d01, d02 and d03) and one-way nesting for d04 innermost domain.
In d04 Large-Eddy-Simulation (LES) technique has been applied.This option replaces the traditional patterns of atmospheric boundary layer using the mathematical model of turbulence originally proposed by [32] and evolved into the latest versions of WRF.This technique allows to evaluate more effectively atmospheric turbulence and allows explicitly calculate all the physical and dynamical processes that characterize the microscale, to deal well cloud formation, turbulence, heat transfer, the exchange of heat or gases fluxes, etc.This technique is considered relevant when the horizontal resolution meteorological model works is below 500 m [7] [33].
Meteorological modeling system developed works operationally in a computing cluster owned by Meteosim S.L. with 27 nodes and more than 260 cores.The operational system developed is updated every 6 hours and offers forecasting for the next 96 hours.

Sensitive Analysis and Calibration
As we have commented previously, the main goal of this paper has been to define a standard methodology to achieve the optimum WRF configuration in every region of the world.To do this, we define experiments using different option combinations.But there are multiple and numerous combinations of WRF options, being not feasible to analyze all of them.For this reason, we have considered next sensitive analysis: The first four-kind of experiments have been applied over a past period of analysis, whilst the fifth kind of experiment has been evaluated in an operational period.
In all WRF sensitivity experiments, we have used a ceteris paribus experimental approach [34].This approach is based into modify only one configuration option and holding all else constant.This approach has always been followed except when exist some WRF restrictions.

Physical Options
First of all, we have realized experiments modifying physical options.A total of 18 experiments have been evaluated progressively, as Table 1 shows.Two of them by varying microphysics schemes, four of them by varying radiation scheme options, three by varying cumulus scheme and eight of them by varying PBL and surface layer schemes (at the same time due to model restrictions).The first numerical experiment corresponds to the default WRF options (defined as INI experiment).Secondly, microphysics experiments (MPH) are analyzed.Microphysic scheme allow to predict water phase transitions in the atmosphere and to consider snow and hail.Thirdly, longwave (LWR) and shortwave (SWR) experiments are analyzed.These schemes define radiation parameters depending on cloud cover, location, gases and aerosols in the atmosphere, time of the year, etc. Fourthly, cumulus experiments (CUM) are analyzed.Cumulus parameterization is used to predict the collective effects of convective clouds at smaller scales as a function of larger-scale processes and conditions.And finally, experiments of PBL and surface layer have been carried out (PBL experiments).PBL and surface layer schemes define boundary layer fluxes (heat, moisture, momentum) and the vertical diffusion process.

Dynamical Options
Usually, sensitivity analysis does not include a specific treatment of dynamical options and are focused on physical options, parameterizations and horizontal resolution.In this research, we have wanted to investigate some dynamical options that they are recommended for very high resolution following recommendations and best practices in WRF [7].We have focused on damping and diffusion options.These options could improve model-top reflection of mountain waves, remove poorly resolved structures and reduce noise at model scales similar to grid-spacing.In Table 2, we show the different experiments analyzed.INI experiment corresponds to default WRF dynamical configuration.Experiment defined as DIN1 incorporates horizontal diffusion 6 th order for all levels in the model.This diffusion can help to reduce noise at model scale, in particular relevant when weak wind occurs.A modification of the horizontal diffusion 6 th order factor defines the second dynamical experiment called DIN2, from 0.12 to 0.36, adapting [54] recommendations.Experiment called DIN3 incorporates Rayleigh option (in comparison to INI experiment) to damping gravity waves produced by topography in upper levels, by adding an implicit gravity wave damping layer near the model top.And finally, DIN4 experiment considers Rayleigh option and the modification of the horizontal diffusion 6 th order factor.

Number of Vertical Levels
There are numerous papers which demonstrate that increasing the number of vertical levels is related to an improvement of the accuracy of the forecasts [56]- [59].For this reason, we have started with the WRF default configuration and we have increased this number to define experiments.Vertical structure of WRF default configuration (version 7) includes 30 vertical layers covering the whole troposphere and resolution decreasing slowly with height in order to allow low-level flow details to be captured.The first 9 levels are inside atmospheric boundary layer (below 1500 m), with the first level at approximately 57 meters, and the domain top is about 50 hPa.
Our experiments are focused to increase the number of layers inside the atmospheric boundary layer.And we have defined two experiments based on WRF working with different verticals levels:  In the three cases, the outer vertical layers of the atmospheric boundary layers are coincident among them, and the difference is the density in the atmospheric boundary layer.In Figure 2, we show a comparison between levels for every experiment and their distribution in sigma coordinates.

Physiographic Model Database
The terrestrial data sets for the WRF model are built using NCEP geographical data.These consist in global data sets for soil categories, land-use, terrain height, annual mean deep soil temperature, monthly vegetation fraction, monthly albedo, maximum snow albedo and slopes.The highest horizontal resolution available in WRF model is 30 arc-seconds, approximately 1 km at the equator, and it is called "GTOPO30".WRF uses land-use categories from US Geological Survey (USGS) 24-category data with a highest resolution also of 1 km.
When very high resolution is required or the region is very complex from a topographical point of view, the accuracy of these databases is not enough, and it is necessary to couple higher resolution databases.In this research, we have evaluated different topographical and land use database information.
In the case of the topography, we have adapted terrain data information belonging to the National Aeronautics and Space Administration (NASA) Earth data (http://gdex.cr.usgs.gov/gdex).The NASA SRTM (Shuttle Radar and Topography Mission) 3'' (approximately 90 m at the equator) [60] is one of the topography database used (hereinafter referred to as SRTM).And the other topographical database considered is The ASTER (Advanced Spaceborne Thermal Emission and Reflection Radiometer) global DEMv2 [61].This information is available at 1'' (approximately 30 m at the equator).
Regarding land uses, two databases have been analyzed: Corine Land Class 2006 (CLC2006) and Climate Change Initiative Land Cover (CCI-LC).The first one has a 100 m resolution and it has been developed for the European Environment Agency (EEA) covering most countries of Europe and Turkey.Whilst, CCI-LC has 300 m as resolution and it has been prepared for the European Spatial Agency covering the entire world.Higher resolution physiographic database has been coupled in WRF following [62] [63] approaches.These databases have been adapted to d03 and d04 domains.
Two experiments has been defined combining topographical and land use databases.• One experiment, defined as HRP1, considering ASTER and CLC2006 information; • And other experiment, defined as HRP2, considering SRTM and CCI-LC databases.
In Figure 3, we show the comparison between topographic and land use information coupled in WRF using different databases over d04 domain.
Using ASTER topographical database we have applied methods described in [15] [64] to evaluate the topographical complexity of the region, obtaining that the mean elevation is 66 m with a standard deviation of 2 m, being therefore a region with a low complexity from a topographical point of view.Regarding land use, domain d03 is 40% covered by water, 7% coniferous forest, 7% non-irrigated arable land, 6% agricultural areas and only a 1% urban areas.

Nudging Options and Data Assimilation
Last experiments realized have been focused on the use of three-dimensional variational data assimilation (3DVAR).The Gridpoint Statistical Interpolation (GSI) system developed by the Development Testbed Center (DTC) has been applied.This system is based on Spectral Statistical Interpolation (SSI) [65] and it is prepared to be coupled into WRF.Our methodology considers nudging model towards observations (defined as observational or station nudging) and analysis (defined as analysis or grid nudging) that they are used for dynamical initialization.Grid nudging is recommended for coarse resolution, whilst observation nudging is recommended for fine scale [66], and both of them can be combined.In this contribution, we have tested both nudging options applied over different domains and considering meteorological observations from: metars, radiosoundings and monitoring stations (data defined as global data); irradiance information from EUMETSAT satellites (data defined as satellite data); and from local meteorological stations managed by AEMET (Spanish National Meteorology Agency).This information has been coupled and combined into WRF defining different testing experiments as we show in Table 3.

Data Sets and Model Evaluation
As we have commented to calibrate and validate WRF model, we have run the model during different periods and we have compared modeled and observed values.Modeled values and observed values from local meteorological stations are compared using different statistics.Values from the local meteorological stations have been considered as observed values.Considered stations are managed by the Spanish National Meteorological Agency (AEMET) and correspond to the stations of Huelva Ronda Este (37.28˚N, 6.91˚W, 19 m a.s.l),CartayaPemares (37.22˚N, 7.08˚W, 15 m a.s.l) and El Arenosillo (37.10˚N, 6.74˚W, 41m a.s.l).
To evaluate the model performance, four statistics have been selected among the large amount of methodologies that can be applied [21] [67]: the Mean Bias (MB), the Mean Absolute Gross Error (MAGE), the Root-Mean-Square Error (RMSE) and the Index of Agreement (IOA).These statistics provide information on how  uncertain a model is, regarding to the observations [5] and according to them a benchmark is given following [68] [69] suggestions.
MB is an evaluation of the data tendency; positive values mean that the simulated values are overestimating the observed ones.Similarly, negative values indicated an underestimation of the simulated values over the real ones.MAGE is used to measure the closeness of the modeled and observed values.IOA provides a measure of the match between the departure of each prediction from the observed mean and the departure of each observation from the observed mean.This is a statistic to evaluate with a single value the goodness of fit of a modeling system with respect to the observations [70].IOA has theoretical range of 0 to 1, with a value of 1 suggesting perfect agreement.RMSE is calculated as the square root of the mean squared difference in modeled and observed values.It is commonly used as a measure of the overall model performance.
Regarding wind direction, statistics have to be considered very carefully due to circular nature of wind direction.For this reason in the case of the wind direction a modification of the traditional formula of MB and MAGE has been applied [3] [71] [72].
The evaluation performed is focused on the inner domains, d03 and d04, since the final aim of this study is to find the best model setup for high resolution domains.Statistical evaluation of the meteorological data is achieved by comparing the modeled parameters to the meteorological station observations of temperature at 2 m, wind speed at 10 m, wind direction at 10 m and relative humidity at 2 m.The statistics used for each meteorological parameter and its benchmarks are shown in Table 4. Wind speed and wind direction are calculated considering calms below 1 ms −1 , as wind direction is not reliable for lower speeds [30] [73].The statistics have been calculated from hourly data of the model and observations.In this paper, the attention is focused on wind performance due to its influence over air quality management and risk management in the zone.The benchmark value of 30˚ is a valid reference value for meteorologically simple areas (locations with low topography complexity and/or land use variance and which meteorology depends basically on the synoptic scale).Typically for the Iberian Peninsula and meteorologically complex areas this reference value is significantly higher [30] [71] [74]- [76].For this reason, statistical showed in Table 4 are complemented with the Directional Accuracy (DACC) statistic [15] [77].This statistic is a parameter that quantifies the percentage of occasions that the atmospheric model uncertainty for this variable is less than 30˚ (Equation (1)) being D i the difference between modeled and observed values.
DACC is applied into the sensitive analysis and model evaluation.For the model, complete evaluation also has been defined parallel statistics for the rest of variables: Temperature Accuracy (TACC, Equation ( 2)), Wind Speed Accuracy (WACC, Equation ( 3)) and Relative Humidity Accuracy (RHACC, Equation ( 4)).Values obtained for these statistics have been considered by the authors as the numerical representation of the accuracy, confidence level or reliability of the meteorological modeling system.

Results and Discussion
In the next sections, we present the results of the sensitive analysis and a complete two years evaluation using the proposed statistics to compare modeled and observed values.

Sensitive Analysis and Calibration
Sensitive analysis results for physical and dynamical options and physiographic model databases are showed in Table 5. Modeled and observed values have been compared for every meteorological parameter: temperature, wind speed, wind direction and relative humidity, and the one that showed best results for the maximum meteorological parameters has been selected as "best case".When results are very similar among them or not conclusive differences are obtained, wind speed and wind direction has been the key variable to decide.Results obtained by the experiments of number of vertical levels (VER), dynamical options (DIN) and physiographical database (HRP) should be compared individually with the INI experiment (based on default WRF configuration).In the case of physical options: microphysics (MPH) experiments should be compared with INI; long wave radiation (LWR) with MPH best case; short wave radiation (SWR) with LWR best case; Cumulus (CUM) with SWR best case; and planetary boundary layer (PBL) should be compared with CUM best case.
Results in Table 5 showed that INI experiment provides high accuracy, fulfilling most scientific recommendations in relation to the model uncertainty for different atmospheric variables.The meteorological variables analyzed have their greatest sensitivity settings vertical levels and in the scheme of planetary boundary layer as we discussed below: • Experiments with different vertical levels show a slight improvement in temperature, relative humidity, speed and wind direction in comparison with INI experiment.Both experiments, VER1 and VER2, significantly increase the computational time (around 30%), assessing the increase of computational cost versus the improvement in the quality of the predictions, VER1 scheme is selected as the optimal for carrying operational weather forecasting.• Experiments with different microphysical schemes show a slight improvement in the wind, while the other variables remain practically constant, and for this reason the SBU-Lin scheme is chosen as optimal for the performance of operational weather forecasting.• Experiments with different long wave radiation schemes show slight improvement in estimating wind velocity using the scheme RRTMG replacing RRTM, while variables temperature, wind direction and relative humidity remain practically constant.In turn, it has been observed that the FLG scheme is not entirely stable, providing computing errors in the 3% of cases, which rules out their use as operational prediction scheme.• Experiments with different shortwave radiation schemes do not show significant improvements in the prediction of weather variables, looking ever so slightly enhanced temperature.In addition, the experiment using the FLG scheme very significantly worse predictions, probably because such scheme long wave radiation is selected RRTMG rather FLG.Dudhia scheme has been selected as optimal for the performance of operational weather forecasting.• Experiments with different cumulus parameterization show improvement in the wind speed, at the same time worsening in estimating the wind direction, while the variables temperature and relative humidity are not significantly modified.Kain-Fritsch scheme is chosen as optimal for the performance of operational weather forecasting.• Experiments with different planetary boundary layer schemes show improvements in temperature using the QNSE scheme and wind speed using ACM2 or MYNN3.While in the case of the wind direction is MYJ and YSU which provides better results.The relative humidity does not seem very sensible to PBL schemes.Schemes as ACM2 or QNSE work best in particular situations (convection conditions or stability, respectively) but we focus our attention on minimize the uncertainty in any meteorological situation.Therefore, Yonsei University scheme is chosen as optimal for the performance of operational weather • Experiments with different dynamic configuration show a slight improvement in predictions of temperature, wind speed and direction, while relative humidity does not show significant differences.In view of the results, the options used in the experiment DIN4 for conducting operational weather forecasting are chosen.• Experiments incorporating very high resolution physiographical databases improve significantly wind speed and relative humidity, while the temperature and wind direction remain practically constant.This result is entirely consistent with expectations: firstly because the topography is not complex in the area, was not expected a significant improvement in wind direction; and secondly, because incorporating CLC2006 land uses provides a higher roughness length, which implies a lower wind speed, associated with the fact that the model in its default configuration tends to overestimate the wind, means that the application of CLC2006 tends to reduce this overestimation.In view of the results, we consider the use of ASTER and CLC2006 physiographic databases for operational weather forecasting.

Data Assimilation
In Table 6, statistical evaluation of the different data assimilation numerical experiments is showed.In this case, as we comment previously, the analysis period corresponds to an operational period of two months.Moreover, initial and boundary conditions are provided by GFS, whilst for the sensitive analysis CFSv2 is used.We can observe that the experiment INI using WRF default configuration offers high accuracy, similar to that obtained from calibration experiments.The nudging experiment that optimizes the accuracy of weather forecasts is NUD5, which corresponds to make grid nudging over d01 (9 km) using satellite data, and observational nudging over d02 (3 km) and d03 (1 km) using local measurement stations.

Optimum WRF Configuration
Once sensitive analysis has been concluded we can summarize the WRF configuration that minimizes the uncertainty and optimize the confidence level of forecasts for operational purposes in the region of Huelva.In Table 7, we show the physical and dynamical options selected, best physiographic databases incorporated, number of vertical levels configured and nudging applied.Options that we do not evaluate but we incorporate in the optimum configuration are: initialization (GFS 0.25˚) and planetary boundary layer scheme (LES) for d04 domain.

Statistical Model Evaluation
In Table 8, the comparison between modeled and observed values from WRF simulations and meteorological local stations is showed.As we comment before, we have considered the period compressed between 01/01/ 2012 and 12/31/2013.In this evaluation, we also incorporate a local meteorological station inside the Port of Huelva (37.20˚N, 6.93˚W) allowing to compare the performance obtained using high resolution (1 km) and very high resolution 0.333 km) and different PBL parameterization, YSU for d03 and LES for d04.
Results show a high accuracy in the prediction of temperature, with MAGE lower than 2˚C.Temperature is underestimated in all stations, being in Huelva Ronda Este where highest underestimation is reproduced.Regarding the comparison between d03 and d04 results over Port of Huelva, we can observe a significant reduction of the accuracy, from 79% (d03) to 57% (d04).In the case of wind speed we obtain different results, when d04 outputs is applied the accuracy increases from 61% to 84%.In the other stations, wind speed accuracy is enough high, with values between 77% and 85%, and RMSE is lower than 2 ms −1 .Wind speed is overestimated in Huelva Ronda Este, Cartaya Pemares and El Arenosillo.Wind direction is also slightly overestimated and the values obtained are very closer to the recommendation value of 30˚.Directional accuracy obtained is compressed between 64% and 67%, and a significant reduction of the accuracy is observed when d04 is applied, from 64% to 47%.Finally, relative humidity prediction shows a very high accuracy, with values between 86% and 92%.

Conclusions
A standard methodology to select the optimum meteorological configuration in any region has been defined.This procedure defines the quality control assurance of Meteosim S.L. to develop high-accuracy forecasting services.We can summarize the methodology in ten steps as presented below.
Quality Control Decalogue of Meteosim S.L. for meteorological forecasting services: 1) To analyze previously experiences about atmospheric modeling in the region of interest or near regions.
2) To analyze previously experiences about atmospheric modeling for similar purposes and sectoral activities in any region in the world (for example experiences in ports or air quality purposes).
3) To select an appropriate meteorological model to adapt to the region of interest, it should be considered the options included in the model that it could be tested and evaluated.
4) To define modeling architecture of domains, horizontal grid resolution and vertical levels depending on synoptic, mesoscale and local meteorological patterns in the region.Also, it is important to consider the aim of the forecast and the final utility (as air quality or weather risk management).8) To incorporate to the operational forecasting system 3DVAR or 4DVAR data assimilation combining grid and observational nudging.9) To validate the weather forecasting system developed during a minimum period of 1 year comparing modeled and observed values using referred statistical parameters.And to compare the metrics obtained with benchmark values.10) To obtain the confidence level or reliability value for every meteorological variable of interest.
In this case, we have applied this methodology in Huelva, a coastal region in Southern Spain.Results obtained after the application of the methodology have defined the operational system to predict meteorological phenomena for the next days in the Port of Huelva.Modeling results will be used as early warning system and it will be a tool that will allow reducing the effects of meteorological hazards and associated risks for the air quality in and out of Port of Huelva.A total of 4150 WRF daily simulations have been realized distributed among numerical experiments of sensitive analysis and model evaluation.Sensitivity analysis developed using different set-up WRF options has concluded that those experiments with different vertical levels are the most influential in improving weather forecasts, followed by the physical configuration (mainly the atmospheric boundary layer), and ultimately the dynamic configuration.Compared to the results of physical and dynamical experiments, the prediction of weather variables is more sensitive to the addition of high resolution physiographical databases.Uncertainty modeling after that sensitive analysis is reduced into a 19% in the case of temperature, a 34% for wind speed, a 2% for wind direction, and a 9% for relative humidity.Statistical evaluation realized using a comprehensive representative period (two full years, 2012 and 2013) has showed that the operational meteorological system designed and developed offers a confidence level or reliability compressed between a 66% and a 79% for temperature, 61% -84% for wind speed, 64% -66% for wind direction and 86% -92% for relative humidity.These values have been obtained using all official measurement stations in the region of interest.
Results obtained from highest grid resolution model output (d04) show that wind speed is improved significantly regarding the same results obtained using d03 output, increasing the accuracy from 61% to 84%, and slightly in the case of the relative humidity.In contrast, accuracy in temperature and wind direction is significantly reduced, from 79% to 57% and from 64% to 47% respectively.We need to evaluate other model options over d04 and use more meteorological measurement information to obtain a better performance evaluation for this very high resolution.

•
Physical options: 18 experiments modifying physical options in relation to WRF default configuration have been realized.Two experiments to analyze optimum microphysics scheme; four experiments to analyze longwave and shortwave radiation scheme; three experiments to analyze cumulus; and eight experiments to analyze optimum planetary boundary layer scheme.• Dynamical options: 4 experiments modifying horizontal diffusion 6 th order option and factor, and damping option have been developed.• Number of vertical levels: 2 experiments modifying the number of vertical levels and their density in the lower layers have been realized.• Physiographical database: 2 experiments modifying land use and topography databases with different resolution have been realized.• Nudging options: and 5 experiments applying grid and observational nudging have been applied.

Figure 3 .
Figure 3.Comparison between different topography and land use databases coupled in WRF.(a) (c) (e) Corresponds to topography database: GTOPO, ASTER and SRTM respectively; (b) (d) (f) Corresponds to land use database: USGS, CLC2006 and CCI-LC respectively.

5 )
If high resolution is required, it should obtain local physiographic databases and couple them into the meteorological model selected.6) To define numerical experiments to analyze the sensitivity to the model to changes in parameterizations and schemes (PBL, cumulus, microphysics, radiation, etc.), vertical levels and physiographical databases coupled.7) To analyze the numerical experiments results comparing modeled values and observed values provided by local official meteorological stations during a representative period.

•
One experiment considering 36 vertical levels, with 15 vertical layers below 1500 m and with the first level at approximately 16 m; • And other experiment considering 42 vertical levels, with 21 vertical layers below 1500 m and with the first level at approximately 8 m.

Table 1 .
Numerical experiments developed corresponding to physics options.

Table 2 .
Numerical experiments developed corresponding to dynamic options.
Figure 2. Sigma level distribution for every experiment defined.

Table 3 .
Numerical experiments developed corresponding to nudging options.

Table 4 .
Statistics used for model evaluation and benchmark values.

Table 5 .
Statistical evaluation of the different physical-dynamical-physiographic numerical experiments carried out for the temperature, wind speed, wind direction and relative humidity using d03 model output.In bold the optimum value for each statistical.Corresponds to the reference benchmark value used.Experiment abbreviations are defined in Table1, Table2 and Table 3. *

Table 6 .
Statistical evaluation of the different data assimilation numerical experiments carried out for the temperature, wind speed, wind direction and relative humidity using d03 model output.In bold the optimum value for each statistical.* Corresponds to the reference benchmark value used.Experiment abbreviations are defined in

Table 7 .
Configuration options selected as optimum for meteorological forecast over the coastal region of Huelva.

Table 8 .
Comparison between modelled and observed values from WRF simulations and meteorological local stations.