Environmental Impact Assessment of Quarrying Plants: An IRF Kriging Solution to Noise and Airborne Dust Monitoring

Corresponding Author: Dario Lippiello, University of Roma Tre, Rome, Italy Email: dario.lippiello@uniroma3.it Abstract: A geostatistical approach is tested for the assessment of the environmental impact generated by industrial plants on the surrounding areas. To this end, a case study of a tuff quarry in the Italian region of Lazio is presented. The method utilised is Intrinsic Kriging with external drift which allows spatial maps of different hazardous agents like noise and airborne dust to be made through the integration of two different kinds of information. The initial data is taken from field measurements, which are carried out in accordance with the current regulations at various points in residential or built up areas near the site. Auxiliary data, integrated with the former as external drift, are obtained from specific propagation models respectively referred to sound propagation in outdoor environments by means of ISO 9613 norm and airborne dust propagation by means of ISC3 by US EPA approach. The use of a combination of these two information sources enabled both acoustic and airborne dust maps to be made which indicate not only the overall spatial trends of the acoustic and airborne climate of the area surrounding the site but also the local variations. Such results are generally difficult to achieve due to the limited data set and usually come at the price of increasing the number of measurements taken in the field. In order to quantify the contribution of introducing this integrated approach, scenarios are presented with and without external drift, using the same field samples as input in each. The results indicate a reduction in estimate variance for the cases analysed when using Intrinsic Kriging with external drift. Furthermore, it was observed from the analysis of cross validation and estimate variance that the use of external drift may also compensate for information loss, which may occur due to a reduction in the number of receptors in the monitor network for example, without affecting the quality of the final map. In the case study, on halving the number of measurements taken while introducing the suggested drift, the quality of the result remains similar.


Introduction
The environmental impact assessment is an important step towards verifying the environmental compatibility of industrial sites with sensitive surrounding areas both during the planning stages of the site and while carrying out the urban planning of the adjoining areas.
Noise and airborne dust emissions may be considered among the most critical hazards to be found in quarries due to the design and nature of this working environment. Not only does dust and noise impinges on the health and safety of workers but they also impair the quality of the air from an environmental perspective.
In this regard the EU Directive and the corresponding Italian legislation, require an initial forecast during the site's planning stages regarding both the acoustic and the airborne dust impact the site might have and regular verification that the conditions are being upheld throughout the working activity according to the nature of the plant. This forecast must demonstrate that the effects of the emissions from the planned industrial plant are compatible with the assigned land use for the area, respecting the limits stipulated by the relevant legislation and by the local authorities. More precisely the categories of land use and their relative limits, as for acoustic climate, are defined by these local authorities for the area under their control (see Table 1). Moreover, the Italian Decree n°152/2006 defines procedures to be adopted in order to control and reduce environmental impact due to airborne emissions from industrial activities. Among these, quarries represent one of the most pollutant as for airborne dust in general and PM 10 in particular.
The norm in force focuses its attention to those procedures to be implemented by the quarry management in order to acquire air quality permit. This procedure can be summarized in the phase of source identification and then their characterization by samples aimed to define the type of the airborne pollutant and its quantitative parameters such as emission rate and concentration.
The geographical boundaries of a given area destined for one particular land use are also indicated on maps drawn up for that purpose. With regard to extractive activities such as quarries, over and above evaluations and noise propagation forecasts, a further guarantee is required to confirm that standards regarding noise pollution and airborne PM 10 emissions are being respected. This entails carrying out a number of both acoustic measurements and PM 10 samples in locations stipulated by the local council and Regional Agency for the Protection of the Environment (ARPA), before and after the extraction activities have started. Moreover, prior to the completion of the site, the manager of the extractive activity must supply a plan of the provisions to be made (and times necessary for their realization), should the measurements indicate that the noise emissions or airborne concentration exceed the limits stipulated. In other terms, for those activities which require field measurements, the authorization for any subordinate work to begin must be granted by the time those measurements are carried out. In such cases, carrying out measurements becomes in effect a noise or air pollution test and complete authorization may only be granted once it has been verified that these measurements do in fact respect the stipulated limits.
The periodical nature of these tests is another characteristic which concerns the extractive sector. The activities are in constant evolution and sound or dust sources like haul trucks moving on unpaved roads may be mobile and thus generate both emissions and a climate which are variable in relation to the particular point which is being quarried at the time. Therefore the pollution tests must be carried out not only at the beginning of the activity, but will need to be constantly updated as the work advances and as the work face of the quarry changes.
Furthermore, while the work is underway it will be necessary to monitor the effects generated on certain receptors (generally those located in a built up area) where the effect of these agents is greater.
For the above reasons, in the extractive sector it may be necessary to include noise and air quality measurements from the field, which regularly check that emissions remain within the limits consented by the law, alongside forecasts based upon sound propagation or PM 10 dispersion models in the atmosphere.
Numerous studies have been carried out which examine the problem both of the acoustic (Neto et al., 2012) and airborne PM 10 impact (Bluvshtein et al., 2011) in relation to extractive activities. They cover the aspect of climate forecast or else that of the evaluation of the impact using specific receptors in the surrounding areas. Again in the same field one should consider studies on the manner of signal propagation from the characterization of specific quarry sources (Ozcan et al., 2012;Patak et al., 1999), as well as those about the quantification of attenuating factors due to effects of dispersion and atmospheric absorption (Cinar and Sensogut, 2009).
The research presented here shows, by means of a case study, that both the model based forecasting and the measurements recorded at the receptors may be integrated into a stochastic frame of reference through Kriging with external drift while gaining indisputable advantages. This technique (Matheron, 1973) is widely used in environmental sciences by means of many application studies (Baume et al., 2008;Van de Kassteele and Stein, 2006;Genikhovich et al., 2002), until more recent applications (Alfaro Degan et al., 2016). It allows the value of a regionalized variable to be estimated at a point in which it is unknown through the use of both its monitored data and data from an auxiliary variable, correlated to the primary and known throughout the dominion. Studies applied in various fields of environmental science, have demonstrated that the use of this auxiliary variable improves the quality of the estimate (Hudson and Wackernagel, 1994;Van de Kassteele et al., 2009).
The result in the case given is that of attaining maps which would be difficult to produce from either of the models on their own or else only from the measurements from the receptors, especially in the case of a restricted number of receptors. These maps show the global effects of sound and airborne propagation in the area surrounding the site and allow for an immediate comparison to be made with maps drawn up by the local authorities indicating the noise limits for each area. Consequently, any problematic areas become evident and points to be investigated are highlighted, such as a boundary between two areas with a different land use assigned to each. This is a considerable advantage with respect to maps obtained from forecast models, which must be compared with field samples, following which the points indicated by the authorities must then be checked with an official noise or airborne dust level test. On the other hand, the general spatial trend of the phenomenon would be obtained by using only measurements, which would be difficult to produce by estimation unless there were a considerable number of measurements recorded in the field. However, since the number of receptors, whose the measurements are required by law, may be reduced, it is particularly useful to adopt the method proposed here which enables the global spatial trend of the pollutant to be obtained without needing a particularly complex monitoring phase. Should periodical updates be required, which is typical in these cases, the advantage would be amplified.

Materials and Methods
The work is organised as follows. In the first phase field tests were carried out and then sampled data were post processed by means of especially written calculus codes. The final part was then oriented to data comparison and analysis.
The entire monitoring program was developed by means of the synergic contributions of two research groups: The Safety and Environment Laboratory and the Applied Physics Laboratory from The Department of Industrial and Mechanical Engineering of Roma Tre University.
In the first part acoustic samples were carried out by means of 30 receptors adjacent to the property in which the extraction plant is located (see Fig. 1). Each measurement was performed in order to guarantee that they were representative of typical conditions so that any events that may disturb the measurements might be determined (e.g. traffic or other manmade noise or dust emission from nearby sources). Subsequently, using the source model derived from ISO 9613 described in the following part, the machinery manufacturers' noise output specifications, as well as geomorphological information about the site, an assessment of the noise impact on the environment was carried out.
The same approach was applied with regard to airborne PM 10 concentration.
In the first stage, a sampling campaign, covering the same selected area, at the same 30 receptors, was carried out. Each gravimetric sample had a two hour duration in order to collect a significant dust mass on the filter, in relation to the selected flow rate. In particular the gravimetric sampler consisted in a sampling pump connected to a PM 10 terminal The flow rate was assessed at 1 m 3 /h. A 5 µm pore size PVC filter was equilibrated in a room for balance at 20°C ±1°C and 50%±5% RH and then weighed with analytical balance (Mod. Exacta series ABT 120-5 DM) with 0.01 mg sensitivity.
In the second stage, using the Gaussian model supplied by US EPA by means of ISC3 software, an evaluation of PM 10 concentration was defined. In this study ISC3 operates in short-term mode (ISCST3) so that the following meteorological data were recorded for each hour of simulation: flow vector (which is the opposite of wind direction), wind speed, temperature, Pasquill Gifford stability class, rural mixing height and urban mixing height. In order to guarantee maximum simplicity in the development of the models, the entire plant was modelled starting from the contribution of the sources present ( Fig. 2 and 3), whose combined effect at any given point of the domain is calculated as the sum of the sound energy or emission rate from the sources which are simultaneously active on site. When the sources were stationary they were allocated to their respective locations within the working area, whereas in the case of mobile sources, in accordance with the FAST methodology (Alfaro Degan et al., 2005), they were allocated to the geographic point which is the geometric median of the area in which the machinery operates.
The environmental impact assessment by both models was developed in the surroundings of the plant into to a square-shaped area, with sides measuring 310 m, in which the quarry is located at its southwestern part (see Fig. 4). This domain area was divided into a regular grid with 700 nodes spaced at intervals of 12 m.
Finally, by means of a specially written program, these data were utilised in the Intrinsic Kriging with external drift method in which the primary variable consists of the set of values of sound pressure and PM 10 concentration obtained through monitoring, whilst the secondary variable is again sound pressure or PM 10 concentration but calculated from the source models. On determining the acoustic or airborne dust climate, to facilitate interpretation of the data obtained, an iso level map was drawn up from the 700 nodes of the grid defined previously. To highlight the contribution made by the proposed estimation method, two different scenarios were implemented.
Firstly a comparison between homologous maps is examined; the map including external drift is compared with the map realized from the same set of measurements but without external drift. In both cases the estimation is carried out using the Intrinsic Kriging method.
Subsequently, in order to test the reliability of the method with a reduced number of receptors, estimates are carried out and corresponding maps are plotted, using a reduced set of receptors with half the number of measurements (i.e. 15 samples rather than 30). Again in this case both the estimates were implemented, both with and without external drift.
Each of the proposed scenarios is verified experimentally by means of validation or, in the case in which all the data are used in the estimate, by means of cross validation using the leave-one-out method. The calculation of the estimate variance and the other reference parameters completes the analysis and provides points for consideration in the discussion section.

Sampling Campaign
The area analyzed is located to the north of Rome, near the Municipality of Riano, which is an important industrial area in Lazio where Tuff is extracted.
An aerial view of the site is shown (Fig. 1).
The final product consists standard sized blocks for use in the building trade. The mineral deposit is divided into various areas and the rock is extracted by means of vertical cutting using the stone cutter (Fig. 2) followed by the block removing machine (Fig. 3) which cuts horizontally as well as crosswise, separating the blocks and allowing them to be removed from the deposit.
The monitoring network permitted the sampling of the Equivalent Continuous Sound Pressure Level (Leq) and PM 10 concentration across the 30 nodes of the grid shown in Fig. 4. Each location was observed for a period of thirty minutes (noise) and two hours (PM 10 ).  In each of the receptors indicated, the L eq is given by: Where: T o = Represents the observation time equal to 30 minutes P 0 = The sound pressure level reference value equal to 20 µPa p A (t) = The sound pressure level, weighted according to curve A, due to the active sources during the observation interval (Pa) At the same time, with regard to the gravimetric technique, each filter was weighted using an analytical balance with a precision of 10 µg, including field blanks. The weight of each sampled filter was recorded beside its corresponding tare weight. The concentration was finally assessed according to the NIOSH 0600 Standard as follows: Where: C = Airborne PM 10 concentration (µg/m 3 ) W 1 = The tare weight of filter before sampling (µg) W 2 = The post sampling weight of filter (µg) B 1 = The mean tare weight of blank filters (µg) B 2 = Mean post sampling weight of blank filters (µg) V = Volume as sampled at the nominal flow rate (m 3 )

Propagation Models: Noise and Airborne Dust
The secondary variable covering the full domain on a regular grid was identified in Sound Pressure and airborne concentration from deterministic propagation models.
As for noise propagation the model calculates the sound pressure level at the receptors, starting from the noise output data of the active sources, whilst considering the dispersion and reflection effects of the signal in the propagation phase. From standard ISO 9613-2 (1996) Attenuation of sound during propagation outdoors is: Where: L S,A,f = The predicted sound pressure level for the single frequency and at a certain distance from the sound source (dB(A)) D I = The directivity index of the source L W,Ai = The sound power level of the i-th source given by the machine manufactures (dB(A)) K 0 = The shape factor relative to the emission characteristics of the machine D s = The geometric divergence term defined as: where, d represents the distance from the source to the receptor (m).
The last term represents the set of factors that contribute to the absorption of sound energy: Where: D atm = The term related to absorption due to meteorological conditions (dB(A)) D g = The term related to absorption due to the ground (dB(A)) D b = That due to the presence of any barriers or obstacles (db(A)) In the case study presented here, the entire plant is modelled starting from the contribution of the four sources present whose combined effect at any given point of the domain is calculated as the sum of the sound energy from the sources that are simultaneously active on site. The sound level at a given point outside the plant is therefore calculated as follows: Where: L PA,E = The sound pressure level at a given point (dB(A)) N = Represents the number of active sources L Pai = Represents the sound pressure level at a given point generated by the i-th source obtained starting from Equation 3 (dB(A)) As for dust dispersion the secondary variable covering the full domain on a regular grid was identified in PM 10 concentration by means of a Gaussian model: Where: C = The hourly concentration at downwind distance (g/m 3 ) Q = Pollutant emission rate (g/s) K = A conversion factor w s = Wind speed (m/s) σ y = Standard deviation of lateral concentration distribution (m) σ z = Standard deviation of vertical concentration distribution (m) V = The Vertical term (include the images method) D = The Decay term (m) us = The mean wind speed at release high (m/s) y = The crosswind distance from source to receptor (m)

Estimation Method: Intrinsic Kriging with External Drift
According to a probabilistic approach to modelling the physical mechanism of sound propagation and PM 10 , dispersion, the n observations z(x) are interpreted as a realization of a random function Z(x) at location x. In this framework Intrinsic Kriging with external drift is the method used to estimate the Sound Pressure Level (SPL) and Concentration (C) in the selected area. This approach allows the estimate of a random primary variable Z(x) known only at a small set of points in the domain while taking into account not only the n available data of the primary variable itself, but also considering another variable s(x) exhaustively known in the same domain. As variable s(x) should be strongly correlated with the primary variable Z(x), the former may be used to provide additional information regarding the spatial structure of the latter.
Considering a linear estimator we may obtain the value of the variable as a linear combination of N data values: in which i λ coefficients are to be determined, together with µ i Lagrange multipliers, as results of the following Kriging system (Equation 9, 10 and 11): Equation 9 and 10 are respectively imposed to ensure an unbiased estimation for primary variable and secondary one, while Equation 8 is the optimization condition that implies the minimization of estimation variance.
In particular, dealing with a non stationary phenomenon the conditional expectation of the primary variable may be expressed as a deterministic trend function: It varies smoothly at local scale and it is expressed as: where, the f l functions are all the monomials of degree less than k (degree of the drift) and the a l 's are unknown coefficients. The total number of monomials, denoted L + 1, is equal to (K+1)(K+2)/2 in R 2 .
In many cases the spatial variable exhibits a trend that is not always simple or easy to interpret and, in such cases, the theory of intrinsic random functions of order k (IRF-k) has proven to be an efficient methodology (Kitanidis, 1985;Chiles and Delfiner, 1999).
The basic idea is that, instead of working with the experimental data, one works with combinations of the experimental data that filter out any drift in the same way that finite differences filter a trend in time series analysis. Such linear combinations are called Allowable Linear Combinations (ALC) or generalized increments of order k when they filter out a drift up to order k.
From an applied point of view the most interesting aspect of IRF-k theory is that the generalized covariance of order k characterizes an IRF-k in the same sense that a covariance function characterizes a stationary random function: thus dealing with non stationary cases, means the K(h) function may not be defined for a given h distance by taking pairs of sampled points.
Under these conditions, the generalized covariance is determined as the result of an automatic fitting process in which the Generalized Covariance is expressed as a given a priori model whose coefficients are chosen according to the Minimum Norm Quadratic Unbiased Estimation (Rao Radhakrishna, 1971) then developed for big data set (Stein, 1986).
In particular the Generalized Covariance is expressed as: in which the three parameters C 0 , b 2m+1 , b 2m are not negative and assessed using the minimum norm quadratic unbiased estimation method. For the purpose, some calculus computer code is developed and run.
In practice once a K order is assigned, the program allows one to define the set of coefficients minimizing the norm of the matrix system, thus providing an unbiased estimate.
The only notable specification refers to the SPL variable itself: since the estimation method is based on a linear interpolation of values (Equation 8), the whole procedure is developed considering the regionalized variable as sound pressure (p a ) and only once the estimation is carried out, the conversion from sound pressure to its corresponding level is carried out. As for PM 10 concentration the conversion is not required.

Results
The acoustic and airborne climate were monitored in the area shown in Fig. 4 in which red symbols represent sampled points.  The acoustic climate was modelled by means of the sources shown in Table 2 which indicates the machinery considered in the case study with the relative noise emission specifications expressed in dB(A). These sound power levels along with the data relative to the attenuation of the signal during its propagation constitute the input data of the model.
The forecast is shown in Fig. 5 by means of an SPL map showing the isobels. This is drawn up from the data obtained from a regular grid with nodes spaced 12 m apart in a square-shaped domain with sides of about 300 m long. The quarry plant is located in bottom part of the map and in the southwestern corner of the area.
For quick reference, in Fig. 6, the resulting data from the monitoring campaign are compared with the values obtained for the same points by means of implementing the described model. Figure 6 a scatterplot is given between the equivalent sound pressure level found at the points monitored (abscissa) and the corresponding values obtained from the model (ordinate). Both series are expressed in dB(A).
In order to employ the estimate using IRF Kriging, an ad hoc algorithm was utilised which could calculate the estimate from a grid of defined points, starting from the available dataset. A further option of external drift allowed auxiliary information to be taken into consideration.    To implement the estimate and to the compare the case, two scenarios were considered: the first with all the 30 measurements (Fig. 7) and the second which utilized 15 survey measurements chosen at random from among the 30 measurements available. Each input dataset was then considered both with and without external drift.
Each result is indicated with its relative SPL map showing the isobels and the quality of the result obtained tested by means of validation (in the cases where further measurements were available), or otherwise using the Leave One Out method cross validation.
For each scenario, Tables 3 and 4 show the estimate parameters, i.e., the degree of K of the intrinsic function and the coefficients which characterise the model of generalised covariance, K(h) defined in equation 14. More specifically, Table 3 indicates the parameters relative to the estimate carried out using 30 measurements, whilst Table 4 shows the estimate when only 15 measurements are utilised.
Lastly, with reference to the populations of the four scenarios, Table 5 indicates the linear correlation R and the mean relative error (given in dB). The bottom row shows the same parameters by the ISO model.
PM 10 concentration was assessed by means of the ISC3 model by American Environmental Protection Agency (EPA).
In particular, the emission factors were defined in accordance with data supplied by means of AP 42 document (USA-EPA, 2004;2006) and then corrected to the local conditions. Each source was characterized with a specific emission rate whose total contribution determined the average hourly concentration all over the selected domain. The total airborne climate was, in fact, assessed as sum of the contributions generated by single sources.
In Table 6 these data are summarized in relation to machine and corresponding model.    Figure 15 the PM 10 isobel map is drawn up. The selected domain is the same square shaped area of 310 meter side monitored in the previous case. The airborne concentration is expressed in µg/m 3 .
Using the same methodological approach shown for acoustic impact, some estimates were carried out using different sets of sampled data. In particular two cases were highlighted using the whole set (30 samples) and 15 samples chosen at random among the former.
After the output concentration map of ISC3 model, a simple estimation by means of Kriging was carried out and the corresponding map was drawn up. Finally the code was applied by considering the data from sampling campaign as the primary variable and the data from the model as an auxiliary one. The output is shown in Fig. 16.
For a further detailed analysis data referred to the application of the IRF theory are summarized in Table 7 whereas in Table 8, the linear correlation R and the mean error (given in µg/m 3 ) are summarized for each scenario.

Discussion
The comparison of the acoustic and airborne maps allows a number of considerations to be made. Firstly, it should be noticed that in the case study of Intrinsic Kriging without drift (Fig. 7), and especially for a subset with a limited number of measurements (Fig. 11), the application returns an estimate with low variability on a small scale in which the generalised covariance constitutes the only nugget effect. This, however, prevents local effects from being perceived. The introduction of drift, on the other hand, represents an improvement in this regard. It is possible to appreciate the single contributions of the two components which are combined to produce to the final estimate; on the one hand the data derived from the sample and the auxiliary information from the model on the other hand. The way in which the external drift contributes to the estimate may be seen in Fig. 9, 13 and 16. In any cases in fact, a trend is seen which recalls the variable structure of the phenomenon and is the result to the integration of the auxiliary information provided by the estimate model with the sample data. For a direct comparison one should note the map in Fig. 5 and 15 which show respectively the acoustic and airborne climate obtained by means of the use of the propagation models. Moreover, the acoustic model utilised in the case study, as may be observed in Fig. 6, is characterised by a generalised underestimation of the sound pressure values. On the other hand the PM 10 dispersion model over estimates sampled values by a 15-20% factor. This aspect stems from the decision to consider the planning stage of the site when implementing the model. Thus no measurement is recorded in order to characterise the sound source (which are modelled according to the noise emission data supplied by the machinery manufacturer) nor the airborne dust ones. This factor has contributed to a certain divergence between the forecast data and the measured data. But this condition has not affected the applicability of the method, however and has allowed its performance to be examined under unfavourable conditions constituting a good stress test. The calibration of the study variable (i.e. the sound pressure level) is guaranteed by the information supplied by the sampling campaign. Therefore the systematic underestimation introduced by the acoustic model or the over estimation introduced by the airborne propagation one, is corrected and is no longer observed in the estimates, especially when the number of measurements are increased (see Fig. 10 or Table 8).
By and large, it is possible to affirm that, all inputs being equal, the introduction of external drift contributes to an improvement in the quality of the estimate and may be observed both from the synthetic indicators in Table 5 and 8.
As evidence of this, it may be appreciated in Table 5 and 8 that the correlation index increases with the reduction of mean error. In this regard, the greater density of the points along the diagonal in Fig. 10 with respect to the corresponding Fig. 8, rather than the trend in Fig. 14 with respect to the trend in Fig. 12, goes to confirm this tendency. A further point for consideration is a comparison of the data relative to the estimate using 15 measurements and drift, with respect to those using all 30 measurements without the auxiliary help of the model (Table 5 rows 2 and 3; Table 8 rows 2 and 3). The data show that a drastic reduction in the number of sampling measurements in the input (in this case equal to 50% of the total) is not translated into a drastic reduction in the quality of the estimate, which on the contrary, according to the correlation factor and the general information given in Table 5 and 8, remains more or less unchanged.
Lastly, on observing the case using all the available measurements and with drift, an improvement in the estimate may be noted with respect to that obtained by means of the model. In this regard one need only compare the data relative to the correlation coefficient which increases from 0.73 to 0.77 and from 0.62 to 0.70.

Conclusion
On the basis of the above it may be concluded that the Intrinsic Kriging with external drift methodology of estimation constitutes a feasible instrument to evaluate the acoustic and airborne dust impact generated by quarry related activities has on the surrounding environment.
More specifically, this methodology allows knowledge of the physical phenomenon of sound propagation or dust dispersion to be integrated with data derived from monitoring, thus integrating a stochastic approach with a deterministic model.
Taking into account results regarding these two hazardous agents, it may be observed that the adoption of the external drift contributes to improve the quality of the estimation in both cases determining a growth of the correlation factor R as shown in Table 5 and 8. Moreover in the case of airborne dust the external drift also contributes determining the reduction of the K-order in the IRF method (Table 7) while the same could not be observed for noise (Table 3 and 4). This outcome seems to suggest a more significant contribution supplied in the case of airborne dust nevertheless a more consistent data set would be required to confirm this hypothesis.
The proposed method is directly applicable in many cases since both sources of information are often available. On integrating this information, even in the presence of few receptors positioned in the area surrounding the extractive site (i.e., with a restricted set of measurements), the results show that it is possible to produce a better map with respect to maps carried out without external drift, i.e. either those made using only information from measurements or those made only from the model. Furthermore, if the number of measurements is the same, the proposed method can plot maps which not only show the general trend, but which also better indicate the local variations.
Research is currently underway to examine how this model responds in other case studies in order to test the way it behaves as the correlation varies between the physical drift of the model and the field measurements.
In conclusion, a possible application is the comparison with the limits set by the authorities for various zones of different kinds of land use and appears particularly useful at points located along the boundaries of these zones where other measurements may otherwise be necessary. Lastly, as sound or airborne dust sources and machinery may be mobile and generate a variable climate, the evaluation of the environmental impact should be updated frequently. If the new configuration should have a smaller number of measurements than before due to fewer sensitive receptors, then it has been demonstrated that, also in this case, the use of Intrinsic Kriging with external drift not only allows a map to be realized but that the quality of the map is comparable to that of previous maps utilizing more receptors but without external drift.

Funding Information
This research was developed by the Scientific Research Fund of Roma Tre University, Engineering Dept.

Author's Contributions
Dario Lippiello: Designed the research plan, participated in all experiments, coordinated the data analysis and contributed to the writing of the paper.
Guido Alfaro Degan: Organized the study, participated in all experimets and contributed to the writing of the paper.
Mario Pinzari: Designed the research plan with Lippiello and critically reviewed the manuscript.

Ethics
All of the other authors have read and approved the manuscript and there are no ethical issues involved.