WO2020025984A1 - Method of use of a lidar device and operatively associated lidar data processing unit for providing real-time monitoring of meteorological parameters - Google Patents

Abstract

Method performed by a lidar device and operatively associated lidar data processing unit, the lidar device being configured to scan aerosol layers of the atmosphere by emitting pulses at a wavelength (λ) and receiving return signals, which are processed by the data processing unit that executes a series of algorithms and uses innovative techniques to provide the Technical Effect of real-time monitoring of meteorological parameters and generate a quantitative/qualitative report of atmospheric layering, detect the Planetary Boundary Layer Height, measure wind speed and derive human visibility. The lidar can operate in a vertical, slant or horizontal direction, in a 2- dimensional or 3 -dimensional operating mode; When measuring in a slant direction, it is continuously re-calibrated during sequential movement at predetermined slant range steps through setting the extinction coefficient αaer(λ,R) and backscattering coefficient βaer(λ,R) to new calibrating values that retain reference of corresponding values of the immediately precedent step.

Description

Method of use of a lidar device and operatively associated Iidar data processing unit for providing real-time monitoring of meteorological parameters

THE FIELD OF THE ART

The invention is directed to a method being performed by a lidar device and an operatively associated lidar data processing unit with the lidar device being configured to scan the aerosol layers of the atmosphere by emitting a plurality of pulses at a wavelength (l) and receiving a plurality of lidar return signals, the latter being processed by the operatively associated lidar data processing unit with a scope of providing real-time monitoring of meteorological parameters and generating a quantitative and a qualitative report of atmospheric layering, detection of the Planetary Boundary Layer Height (PBLH), wind measurement and derivation of human visibility.

THE BACKGROUND OF THE INVENTION

It is known to use a lidar device to measure transmission, extinction end back- scattering in the atmosphere, in order to detect gases or particles therein and determine their concentration and distance from the site of the lidar. By way of e xample V.E. Derr, "Estimation of the extinction coefficient of clouds from multiwavelength lidar backscatter measurements", Appl. Opt., Vol. 19, No. 14, refers to investigation of clouds, The aforementioned publication describes a lidar which transmits a laser beam via a transmission lens sys tern. The transmitted laser beam is absorbed and scattered in the atmosphere. The scattered light is measured by a receiving lens system which has an aperture angle approximately equal to that of the transmitted beam. The extinction coefficient of the atmosphere and consequently the density of admixed gases and particles can be determined from the intensity of the scattered light, the known scattering cross-sections of various gases and particles and the spectral adsorption, since the laser wavelength is known.

Monitoring of weather parameters and acquiring reliable weather forecast becomes increasingly important in the present era of greenhouse effect and climatic change. Weather forecasting is the application of science and technology to predict the conditions of the atmosphere by collecting quantitative data about its current state at a given place and time and using meteorology to project how the atmosphere will change. Weather forecasting now relies on computer-based models that take many atmospheric factors into account; yet, a high inaccuracy is associated with weather forecasting that is, amongst other things, due to the localized nature of measurements and assumptions made in relation to extending their validity over a broader region and disregarding the chaotic nature of the atmosphere. Furthermore, inaccuracy is associated with the error involved in measuring the initial conditions and to an incomplete understanding of atmospheric processes.

Aviation safety, i.e. flight and ground safety at airports is an example of a technical field wherein monitoring of meteorological conditions and deriving visibility is an application of critical importance. As reported by the Commission for Instruments and Methods of Observation (CIMO-WMO-No. 807, 2008), estimation of the atmospheric visibility is a complicated issue that is because of the psychological and physical condition of the observer and has to do mainly with the atmospheric extinction coefficient compared with molecular and aerosol particles in solid or liquid form in the atmosphere. The extinction phenomenon comes mainly from scattering and less from absorption of light. The visibility estimation varies and is based upon individual perception, the light source characteristics and the transmission factor. A major flight safety issue is therefore to obtain visibility measurements, especially from the pilot’s point of view and to be able to have accurate measurements of some meteorological parameters, not only above the airport, but also in Slant Range (SR) directions to provide visibility information to the pilot.

During the last years many operative techniques have been developed to estimate the atmospheric visibility and provide real time meteorological conditions at airports. For instance, the forward scatterometers use the in situ forward atmospheric scattering technique with light beams transmitting and receiving units to estimate visibility ranges and the transmissometers, which measure the extinction coefficient through forward scattering but the transmitter and receiver are spaced apart, usually at distances 10-75 m, the scintillometers, which measure the sensible heat flux, the ceilometers, which measure the aerosol backscatter coefficient, the distrometers which measure precipitation, the telephoto- meters/cameras, which measure the daylight and object atmospheric contrast visibility lidars, which measure the laser beam extinction and/or backscattering, the cellular networks which can detect atmospheric conditions by cellular power transmittance attenuation and the nephelometers which measure the aerosol scattering over a wide angle 0-180°. Radars are able to provide weather conditions parameters like rain, giant aerosols and cloud properties. Finally, visibility may be estimated by the human eye based on previously known distances along predefined locations around airports.

These above mentioned techniques applicable in the measurement of visibility and analogous techniques applicable in measurements of other meteorological parameters, have many disadvantages, as they are based mainly on in situ data, presuming homogeneous atmospheres, whilst, when slant range measurements are considered in“real” atmospheres, these techniques become less operational, due to large approximations considered, especially, under unstable atmospheric conditions.

Additionally, automated algorithms for atmospheric layering detection and noise subtraction have been proposed for several years now. The detection was made typically using the slope method which also is disadvantageous because of the presumption of a homogeneous atmosphere. Up to date, there has been disclosed no means for comparing aerosol backscatter and extinction coefficients derived from vertical and slant range (multi-angle) measurements in non-homogeneous atmospheres, except using synthetic lidar data in homogeneous atmospheres.

Further, whilst prior art may provide at increased costs, less accurate measurements of weather parameters, already manifested, such as the values of visibility or wind speed, at a particular location and a particular time, prior art is deficient in that it fails to provide continuous real-time accurate data and in that it cannot provide prediction of changes of these values in the course of near future. Accordingly the algorithms used in processing such lidar data of the prior art are not capable of providing such perspectives.

The Planetary Boundary Layer and the Height thereof (PBLH) stands for one of the most important parameters of meteorology and atmospheric phenomena appearance and can be used to provide reliable data for weather forecasting. According to the geographical and geodesic location, the sun’s radiance and the season, a variability of weather phenomena and growth of PBLH may occur. Many techniques have been developed in order to predict and detect this height safely, with the most accurate one to be the expensive, radiosonde method, with in situ measurements, that naturally cannot be endlessly repeated to provide a near real time estimation of the PBLH. In a specific location PBLH cannot also be identified continuously by satellite means that may only estimate the PBLH at a specific location periodically as their predetermined path is arranged to pass above this specific location. No integrated system of providing actual measurements of weather parameters as manifested that is not based on estimates and does not necessitate faulty assumptions can be found in the prior art. Accordingly no capacity is available in the prior art for continuous real time forecasting of weather that could be based on accurate continuous real time measurements of atmospheric layering and PBLH. Furthermore, there has been no set of algorithms proposed for processing measurements and providing reliable results on current values of meteorological parameters and of the trends thereof in the course of time.

Having this in mind, the object of the invention is to provide the technical effect of accurate real time measurement of meteorological parameters with a lidar device that affords to provide this effect at low cost in association with a series of novel algorithms that are being processed and of concurrently applied novel techniques by the lidar data processing unit operatively associated with the lidar device, wherein the lidar device may advantageously provide the above data also by measurements effected with a single beam and at a single wavelength.

An object of the invention is to provide a continuous, real time monitoring of an all-inclusive spectrum of meteorological parameters that can be processed with the proposed series of novel algorithms that are being processed and of concurrently applied novel techniques by the lidar data processing unit to provide through use of varying combinations of the applied algorithms future forecasting of meteorological parameters that could prove decisive in affording best practices in handling eminent severely adverse weather changes or handling anthropogenic hazards such as a fire or the dispersion of a toxic layer emanating from an industrial accident.

An object of the invention also is to provide the lidar device with an array of additional sensors that may further maximize the continuous, real time monitoring of an all-inclusive spectrum of meteorological parameters, such additional sensors including temperature and humidity sensors that may provide data for the determination of dew point that is an important feature to be recorded in a variety of applications.

Accordingly an object of the invention is to provide the real time all-inclusive spectrum of meteorological parameters data and means of processing these data obtained by the invention to a plurality of end users and a plurality of related authorities being involved, including applied Meteorology for current weather and forecasting, civil protection for pollution monitoring and early warning of ash or other harmful substances and aerosols approaching a populated area or for the movement of aerial masses and the detection and monitoring of swarms of mosquitos and birds transitioning from an area to another, thereby addressing civil health issues, for biomass burning thereby providing early alarm for fire handling and securing forest, animal and human life, for Aviation Safety through providing Visibility and full report of aviation sensitive meteorological parameters to tower controllers and to pilots of aircrafts, as well as to scientific bodies, such as universities and research institutions for the advancement of related scientific research.

It also is an object of the invention to provide the real time all-inclusive spectrum of meteorological parameters data and means of processing these data for Space Applications particularly related to Observation of the Earth or other planets, through Satellite systems with easily applicable space settings, even without the need of Doppler devices, but merely with the single wavelength lidar device of the invention, thereby advantageously reducing the overall load of satellites dedicated in performing such monitoring.

It is another object of the invention to provide measurements of the aforementioned all- inclusive spectrum of meteorological parameters with the lidar device pointing at the vertical, horizontal or slant direction and maximize the accuracy of data obtained through recalibrating the lidar device in each incremental step of slant measurements, such recalibrating and associated processing algorithms and techniques of the invention thereby enhancing the accuracy of data obtained relating to detection of atmospheric layering in any given direction. A further object of the present invention is to provide suitable algorithms and techniques, in order to provide detection of the PBLH, of an extensively wider area and not merely of a specific location, with a 3 dimensional lidar in vertical or slant pointing.

It is also an object of the invention to provide with these algorithms the technical effect of secure and precise measurements irrespectively of weather conditions, these algorithms being adjustable by the user, in order to work automatically at any location and at any time of the year with the lidar device pointing in either vertical or near vertical or slant or horizontal direction.

It is a further object of the invention to provide accurate real time wind speed measurements through tracing of a real time layering and distribution of atmospheric layers to provide results that are of critical importance to airport tower controllers, atmospheric scientists and remote sensing community, who may find these algorithms extremely helpful at their area of expertise.

It is an object of the invention to provide the proposed set of algorithms and techniques for estimation of the atmospheric layering, retrieval of the PBLH and provision of data for the abovementioned all-inclusive spectrum of meteorological parameters through the actual structure of aerosol backscattering coefficient p_aer(k,R) and/or aerosol extinction coefficient ot_aer( R), and/or Range squared Corrected lidar Signal (RCS) in 3 dimensions and through creation of different kinds of“mathematical tools” based on the backscattering lidar signal, the concordance of which is based upon the approximation of PBLH retrieval, according to the Bayes theorem.

SUMMARY OF THE INVENTION

The invention proposes a method being performed by a lidar device located on ground, maritime or space environment and a lidar data processing unit operatively associated with said lidar device to provide real-time monitoring of meteorological parameters through detection of atmospheric layering including detection of the Planetary Boundary Layer Height (PBLH), said lidar device being configured to scan the aerosol layers of the atmosphere by emitting a plurality of pulses at a wavelength (l) and receiving a plurality of lidar return signals, each of said signals providing return signal parameters (p) comprising a Range Squared Corrected lidar Signal (RCS) that is the received power P’(k,R) after atmospheric and electronic noise background (BG) correction, a range (R) dependent variable extinction coefficient a_aer(k,R) and a range (R) dependent variable backscattering coefficient a_er(k,R), characterized in that:

said lidar device being selectively operational in a vertical, slant or horizontal direction, in a two-dimensional or 3 -dimensional operating mode;

said lidar device being calibrated through setting said extinction coefficient a_aer^,R) and said backscattering coefficient p_aer(k,R) to a zero value for non-aerosol presence and detection of molecular layer whilst emitting in a vertical direction and is continuously re-calibrated during sequential movement at predetermined slant range steps through setting said extinction coefficient a_aer(k,R) and said backscattering coefficient p_aer(k,R) to new calibrating values that retain as reference the values of a_aer(k,R) and/or or p_aer(/.,R) of the immediately precedent · step; said lidar device operating with a lidar ratio C(k,R)= a_aer(k,R)/p_aer(k,R), said lidar ratio (C) being set at a predetermined value that corresponds to the type of atmospheric layer being detected, said predetermined value being deducted through the following steps:

obtaining average values of a_aer(k,R) and corresponding range limits from an established diagram of a_aCr()sR) versus visibility depicting measurements employing a wavelength of 550 nm;

calculating a ratio for each pair of adjacent atmospheric conditions, including a ratio for sky crystal clear / sky clear, a ratio of sky clear / light haze, a ratio of light haze / haze, a ratio for haze / thin fog, a ratio for thin fog / light fog, a ratio of light fog / moderate fog and cloud;

applying said ratios for pairs of adjacent atmospheric conditions, said average values of a_acr(k,R) and corresponding range limits in an established diagram of a_aer(k,R) and a_er(k,R) versus wavelength to obtain values of a_aci(k,R) and p_aer(k,R) at the wavelength (l) of operation of said lidar device;

said lidar data processing unit operatively associated with said lidar device being configured to process said lidar return signals to generate a quantitative and a qualitative report of atmospheric layering, detection of the (PBLH), wind measurement and estimation of visibility.

Preferred embodiments of the invention present solutions through use of the abovementioned method for providing accurate real time monitoring of meteorological parameters including atmospheric layering and distribution, provision of the PBLH, wind speed and visibility measurement.

BRIEF DESCRIPTION OF THE DRAWINGS

The invention will be best understood by those skilled in the art by reference to the accompanying drawings in which:

Fig. 1 shows an illustrative diagram of the lidar device of the invention being installed at ground level and adapted to emit signals and receive signal responses in a vertical, horizontal or slant direction with a typical runway of an airport being depicted adjacently to the lidar device.

Figs. 2a-2c show well established and widely acknowledged diagrams as published in R. M. Measures,“Laser remote sensing. Fundamentals and Applications,” Krieger, Sys No 9247, MEA 621.3678 (1992), and in particular:

Fig. 2a is a plot of Aerosol volume backscattering coefficient (in rrf'sr ¹) as a function of wavelength for different types of clouds and haze.

Fig. 2b is a plot of Aerosol extinction coefficient (in m ¹) as a function of wavelength l (pm), and

Fig. 2c shows a diagram of the Variation of atmospheric extinction coefficient in km ¹ with visibility range Rv, at the wavelength of 550 nm and e=0.02.

Fig. 3 shows results of the application of“DENOISING 1/2” algorithms of the invention on values of a_aer (k,R) retrieved at 355 nm.

Fig. 4a shows a diagram of altitude versus time depicting an actual RCS signal acquired with the lidar device of the invention.

Fig. 4b shows a diagram of the type of atmospheric layer versus range based on the signal depicted in Fig. 4a, wherein“VISIBILITY” and“WEATHER PHENOMENA” algorithms of the invention have been applied with the scaling from darker to lighter color indicating variation of the type of atmospheric layers being scanned.

Fig. 4c shows another diagram of altitude versus time depicting an actual RCS signal acquired with the lidar device of the invention.

Fig 4d shows a diagram of the type of atmospheric layer versus range based on the signal depicted in Fig. 4c, wherein“VISIBILITY” and“WEATHER PHENOMENA” algorithms of the invention have been applied with the scaling from darker to lighter color indicating variation of the type of atmospheric layers being scanned.

Fig. 5a is a diagram of altitude (height) versus time depicting an actual RCS signal acquired with the lidar device of the invention at l064nm from the LRSU NTUA (Laser Remote Sensing Unit of the National Technical University of Athens) data at 12:29:40 UTC of 01/02/2016.

Fig. 5b is a diagram presenting p_aer(k,R) values versus range (R) at the time of receipt of data from the RCS signal of Fig. 5a.

Fig 5c is a diagram presenting beginnings and endings of atmospheric layers versus range (R) at the time of receipt of data from the RCS signal of Fig. 5a as provided by MUSTI-L/D algorithm of the invention.

Fig 5d shows a black and white visualization of the type of atmospheric layer versus range (R) at the time of receipt of data from the RCS signal of Fig. 5a.

Fig 6a shows an RCS signal acquired at l064nm from LRSU NTUA data (1/6/2014 at 06:00 UTC).

Fig 6b is a zoomed view of Fig. 6a for the time period of 07:44:30 till 08:07:50 UTC where MUSTI-PBLH gave PBLH estimation.

Fig. 6c is an indicative record of p_aer(k,R) at that period of time (07:46:10 UTC).

Fig 7a shows an RCS signal acquired at l064nm from LRSU NTUA data (31/10/2011).

Fig 7b is a zoomed view of Fig. 7a for the time period of 14:33:10 till 14:56:30 UTC where MUSTI-PBLH gave PBLH estimation.

Fig 7c is an indicative record of p_aer(k,R) (14:43:10 UTC).

Fig 8a is a comparison of the temporal evolution of the PBLH as derived using the extended Kalman filter (EKF) technique (left-hand side images) and “MUSTI”, “DECIS” and “WEATHER PHENOMENA-PBLH” algorithms (right-hand side images), for different meteorological conditions: clear sky case, cloud case, dust, Etesian flow and sea breeze case. The arrows denote the recording timeframe, with the curved dashed line (left side) denotes the PBLH produced by EKF and the white dashed lines (right side) denote the averaged PBLH as produced by “MUSTI”, “DECIS” and “WEATHER PHENOMENA-PBLH” algorithms.

Fig 8b is a comparison of the the temporal evolution of the PBLH as derived using the extended Kalman filter (EKF) technique (left-hand side images) and“MUSTI”,“DECIS” and “WEATHER PHENOMENA-PBLH” algorithms (right-hand side images), for different meteorological conditions: Sea breeze case and cloud case. The arrows denote the recording timeframe, with curved dashed line (left side) showing the PBLH produced by EKF and the white dashed lines (right side) denoting the averaged PBLH as derived by “MUSTI”, “DECIS” and“WEATHER PHENOMENA-PBLH” algorithms. Fig. 9a shows a 3D scanning lidar performing the scans for a so called FASTPLAN technique (FAST way for production of PLANetary boundary layer height), by scanning uniformly, at steps of 15° or 7.5° or even 5°, in order to acquire a scan of a wide area of the sky for PBLH retrieval.

Fig. 9b presents the FASTPLAN technique in slant measurements and vertical columns presentation of p_aCr(k,R) (and/or a_aer(k,R), RCS).

Fig. 9c is a“WEATHER PHENOMENA-PBLH” presentation in 3D, using FASTPLAN technique. PBLH presented in white line.

Fig 10 is an RCS signal acquired at l064nm from LRSU NTUA data (31/10/2011) with a zoomed file for the time 14:34:50 UTC where “WEATHER PHENOMENA-PBLH” produced a PBLH estimation of 1605m appearing with the white line.

Fig. 11 is a flowchart of “DECIS”,“MUSTI” and“WEATHER PHENOMENA-PBLH” algorithms.

Fig. l2a shows lidar signals which are being depicted as follows: the older one with a darker line and the newer one with the lighter line in a series of recorded data files - signals.

Fig. l2b shows examples of similarity‘frames’ captured where a speed measurement is being carried out.

Fig. 12c is a presentation of slant lidar measurements where (Up) the darker the color stands for the upwards wind and lighter, downwards wind direction and (Down) darker the color for the east wind directions and lighter for the west wind direction. The rest of the air mass missing between continuously measurements applies for 3rd dimension wind movement, by using high repetition frequency of beam in pulsed signals. The speed of the 3rd dimension movement of the air mass could be calculated through the rate of air mass loss in that (+ or -) direction.

Fig. l3a is a presentation of the“WEATHER PHENOMENA-WIND” algorithm of four sequential lidar data recording files every lOs.

Fig. 13b shows the tracking of a certain atmospheric layer beginning at about l500m from the equipment and travelling until l750m in 30s, thereby presenting a layer tracking (wind) speed of 8.3m/s.

Fig. 14 is the flowchart of the “3D STEPPING TECHNIQUE”, “DENOISING 1/2”,

“WEATHER PHENOMENA” and“VISIBILITY” algorithms. The l— . jl sign denotes the repetition of the above steps from vertical to the horizontal pointing of the lidar device, for the RVR (Runway Visual Range) measurement.

DESCRIPTION OF PREFERRED EMBODIMENTS OF THE INVENTION

The invention will be hereinafter described by reference to the illustrative embodiments presented in the accompanying drawings.

The invention relates to a method being performed by a lidar device and an operatively associated lidar data processing unit, wherein the lidar device and operatively associated lidar data processing unit is adapted to provide real-time monitoring of meteorological parameters through detection of atmospheric layering including detection of the Planetary Boundary Layer Height (PBLH). The lidar device is configured to scan the aerosol layers of the atmosphere by emitting a plurality of pulses at a wavelength (l) and receiving a plurality of lidar return signals, each of these signals providing return signal parameters (p) comprising a Range Squared Corrected lidar Signal (RCS) that is the received power P'(k,R) after atmospheric and electronic noise background (BG) correction, a range (R) dependent variable extinction coefficient a_aer(k,R) and a range (R) dependent variable backscattering coefficient p_aer(k,R). A series of algorithms are proposed that use the data obtained by the lidar device to provide processing of the lidar return signals to generate a quantitative and a qualitative report of atmospheric layering, detection of the (PBLH), wind measurement and estimation of visibility.

The performance of this lidar system and algorithms has been extensively tested and was found to provide advantageously precise real-time results that provide solutions in problems being encountered in this technical field by systems of the prior art.

Backscattering lidar data from a multiwavelength Elastic-Raman lidar of LRSU (Laser Remote Sensing Unit) of the NTUA (National Technical University of Athens), as well as a 3 dimensional lidar have been used in the aforementioned testing procedure. The advanced stationary lO-wavelength elastic-Raman-DIAL lidar system was located in Athens (37.96°N, 23.78°E, 220 m) at the premises of NTUA. This system is configured with a capacity to perform, simultaneous independent measurements of RCS, a_aer(k,R) and p_aerU,R) (at 355, 532 and 1064 nm). The system was also capable of measuring the water vapor in the troposphere with H₂0 Raman channel at 407nm. The LRSU Raman lidar is part of the ACTRIS Network and performs regular validations of the CALIOP lidar on board the CALIPSO space-borne platform.

The algorithms were run at 355 nm and the results obtained were compared to those obtained simultaneously at 1064 nm, under varying weather conditions. The algorithms of the invention are proposed to be incorporated to the data processing codes of a 3D scanning lidar system. This system preferably is a multiwavelength with 355nm co-polar, 355nm cross- polar and 387nm nitrogen Raman channels.

The tests have been carried out with the lidar emitting, mainly, at 355nm because of the eye safety conditions, provided on that region of the spectrum. Particularly, eye safety reaches up to 1 km from the lidar system according to the EU standard on laser safety EN 60825-1 :2007. The data obtained were transposed to the visual spectrum (at 550 nm), following the World Meteorological Organization (WMO) and the International Civil Aviation Organization (ICAO) rules of daytime visibility.

THE LIDAR EQUATION

A lidar device is known to operate in accordance with the following equation:

where, P(k,R) is the received power, l is the laser wavelength, R is distance, P_0L is the power of the transmitted laser pulse beam, RF is a reference altitude for which a molecular atmosphere is assumed, A₀ is the diameter of the receiver’s telescope, AR is the spatial resolution of the lidar, x^) is the geometrical form factor, q(R) is the lidar opto-electronic efficiency, while a_aer^,R) and a_Ray ) R) are the extinction coefficients for Mie and Rayleigh scattering, accordingly.

According to Klett [J. Klett,“Stable analytical inversion for processing lidar returns,” Appl. Opt. 20, 211-220 (1981) and“Lidar inversion with variable backscatter extinction ratios,” Appl. Opt. 24, 1638-1643 (1985)], the lidar equation can be solved to provide the aerosol extinction coefficient (a_aerU,R)) and aerosol backscattering coefficient (b_3bG(l^)), assuming a relationship between these two coefficients, the so-called lidar ratio C(>.,R)= a_aer^,R)/ p_aer(^R) :

S(R) º Inti* (A, R) R² ] = \n[RCS] (6)

Again, RCS is the Range Squared Corrected lidar Signal, and P’(7,R) is the received power after atmospheric and electronic noise background (BG) correction:

F (A, R) = R(l, R)-BG (7)

In case of vertical or near vertical measurements, a reference height is being set to designate the height at which the atmosphere is purely molecular. In this case, it is possible to retrieve the values of p_{aer( A},R) and/or a_aer(7,R) from that height to the ground. These measurements are, thus, able to provide a clear view of the aerosol load located along the line of sight (LOS) of the lidar beam.

The Aerosol Optical Depth (AOD) marked as (t) is given by Equation (8) below:

^,o,a, (^ ^R) = ^TC_t ^{+ T} _C, + - + ^T _C, (⁹)

where, R' is the distance at which a_aer^,R) has been measured with a range resolution of 1 bin, that was equivalent to 7.5 m in the Lidar device being used in the measurements carried out to illustrate operational performance of the present invention. The total AOD (xtotai) is the sum of xcn, which is the individual optical depth (n=l,2,3,...) of the (n) different atmospheric layered areas of measuring.

The lidar device of the invention is selectively operational in a vertical, slant or horizontal direction, and can be located in the ground, maritime or space environment, in a two- dimensional or 3 -dimensional operating mode and in cooperation with the operatively associated lidar data processing unit and the algorithms processed thereby provide a substantially improved accuracy in monitoring meteorological parameters through detection of atmospheric layering including detection of the Planetary Boundary Layer Height (PBLH).

An important characteristic of the invention is the employment of a lidar ratio C(L,R)= a_aer(7,R)/p_aer(7,R) that is being set at a predetermined value that corresponds to the type of atmospheric layer being detected. A classification of these variable atmospheric layers and variable values of a, b and C is presented in Table 1 herein below, this Table being constructed through use and elaboration of the diagrams presented in Figs. 2a-2c that have been derived from the well established and widely acknowledged publication in R. M. Measures,“Laser remote sensing. Fundamentals and Applications,” Krieger, Sys No 9247, MEA 621.3678 (1992). The elaboration of the diagrams is conducted as follows: The average values of a_aer(k,R) and corresponding range limits for each atmospheric condition depicted in the diagram, namely for sky crystal clear, sky clear, light haze, haze, thin fog, light fog, moderate fog/cloud, are being received from Fig. 2c that presents variation of a^(k,R) versus visibility depicting measurements employing a wavelength of 550 nm. It is noted that the wavelength of 550 nm is a median value of the wavelength range (400-700 nm) applicable to human vision. It is well known that the human visual spectral range spans from 400 nm to 700 nm, peaking at 555 nm for daytime and 507 nm for night conditions, according to CIMO (Commission for Instruments and Methods of Observation) for relative luminous efficiency of monochromatic radiation.

Next a ratio is being calculated for each pair of adjacent atmospheric conditions, including a ratio for sky crystal clear/sky clear, a ratio of sky clear/light haze, a ratio of light haze/haze, a ratio for haze/thin fog, a ratio for thin fog/light fog and a ratio of light fog/moderate fog or cloud, such ratios indicating the difference in strength of adjacent atmospheric conditions. It is herein noted that the values of these diagrams of Fig. 2a for b_3a-(l,II) and of Fig. 2b for a_aer(k,R) are calculations of a computer program called MAPP (Modular Atmospheric Propagation Program) developed by SRI (Stanford Research Institute), for maritime measurements, because of the global nature of its application.

Then Figs. 2a, 2b are used, wherein“high altitude haze” is considered equivalent to the condition designated as “Light Haze” in Fig. 2c and “cumulous cloud” is considered equivalent to the condition designated as“Moderate Fog” in Fig. 2c. Further care is taken to convert the logarithmic base values and dimensions of km ¹ of Fig. 2c to the dimensions of m ¹ as a_aer^,R) is shown in Fig. 2a.

Having made the aforementioned analogies, the hereinabove strength ratios for pairs of adjacent atmospheric conditions and average values of a_aer(k,R) and corresponding range limits are being applied in Fig. 2a showing values of a_aer(k,R) versus wavelength and the lines/curves of Fig. 2a to go from the wavelength of 550 nm to a desired wavelength value, that is the value of the wavelength (l) of operation of the lidar device available for conducting measurements. In the case of the invention, measurements were carried out at the wavelength of 355 nm and therefore Table 1 herein below presents approximate ranges of values of a_aer(^R) and p_aCr(?-,R) at the wavelength l=355 nm.

TABLE 1

THE WEATHER PHENOMENA (W.PHEN) ALGORITHM

The above Table 1 is hereinafter being used by the algorithm designated“WEATHER

PHENOMENA” to provide a qualitative profile of layering of the atmosphere.

The“Weather Phenomena” algorithm is adapted to run with lidar wavelength set at 355 nm for generating a qualitative report of atmospheric layering, being adapted to check the retrieved values, through the respective C(k,R), of p_aer(L,R) and a_aer(k,R) providing a qualitative report of atmospheric layers as follows:

i) For every 1 to 2 space bins the algorithm checks the retrieved values of p_aer(k,R) and through the respective C(k,R), the a_aer(k,R),

ii) If p_aer(k,R) > 2*1 O ³ m ¹ sr ¹, then the outcome is Cl =“Cloud, no visibility”, with C(k,R) = 10 sr and a_aer(k,R) > 2*l0 ² m ¹,

iii) If 2*1 O ⁴ < p_aer(k,R) < 2*l0 ³ m ¹ sr ¹, then the outcome is MF=“Moderate Fog, max visibility 200 m”, with C(k,R) = 10 sr and 2*l0 ³ < a_aer(k,R)< 2*l0 ² m ¹,

iv) If 8*10 ⁵ < p_aer(k,R) < 2* 10^-4 m ¹ sr ¹, then the outcome is LF=“Light Fog, visibility 200 m to 1 km”, with mean C(k,R) = 8.125 sr (6.25 < C^,R) < 10 sr) and 5*1 O ⁴ < a_aer(k,R) < 2*l0 ³ m ¹, v) If 4*10 ⁵ < _aer(^R) < 8*10 ⁵ m ¹ sr ¹, then the outcome is TF=“Thin Fog, visibility 1 to 2 km”, with mean C(k,R) = 6.875 sr (6.25 < C(k,R) < 7.5 sr)

1

vi) If 2*10 ⁵ < p_aer(k,R) < 4* 10 ⁵ m ¹ sr ¹, then the outcome is HZ=“Haze, visibility 2 to 4 km”, with mean C(k,R) = 6.25 sr (5 < C(k,R) < 7.5 sr)

vii) If 8*1 O ⁷ < _aer(k,R) < 2* 1 O ⁵ m ¹ sr ¹, then the outcome is LH=“Light Haze, visibility 4 to 6 km”, with mean C(k,R) = 27.5 sr (5 < C(k,R)< 50 sr)

viii) If 2*1 O ⁷ < P_aer(k,R) < 8* 1 O ⁷ m ¹ sr ¹, then the outcome is SC-“Sky Clear, visibility 6 to

ix) Next the“Weather Phenomena” (W.PHEN) algorithm calculates the AOD according to Equations (8) and (9) from the lidar position towards its pointing direction.

x) If P_aerCW < 2* 10 ⁷ m ¹ sr ¹ (d_aer(k,R) < 1 * 10 ⁵ m ¹), no color discrimination is being made and the visibility goes higher than 10 km.

xi) A check of the AOD total value as derived from equation (9) in comparison with a_aer(^R) is calculated and visibility versus R is being constructed.

In accordance to the above data, visibility can be calculated and a categorization of atmospheric conditions-layering is being made at slant and/or vertical pointing lidar direction. So, with a 3D scanning lidar system, a profile of fully atmospheric layering at different angles, azimuth and vertical, can be provided in near-real time. Figs. 4a-4d illustrate examples of use of the“weather phenomena” and“visibility” algorithms of the invention. Whilst the abovementioned limits of values are related to measurements being made with the lidar device operating at 355 nm, an analogous classification may be elaborated for varying wavelengths of operation. It is a characteristic of this W.PHEN algorithm that it can provide a qualitative identification of prevailing atmospheric conditions after accurate measurements and not just estimations and further providing a 3D profile and classification of atmospheric conditions.

THE 3D STEPPING TECHNIQUE

There has been a lot of questioning on whether aerosol backscatter p_aer(k,R) and extinction a_aer(\R) coefficients can be derived from operational lidar slant range measurements without any large error estimations, resulting from assumptions like atmospheric homogeneity conditions presumed, in a large scale (e.g. typically over a few hundred meters) and on automatic retrieval of these values. Such questioning is by way of example provided by W. Vieeze, J. Oblanas, R.T.H. Collis,“Slant range visibility measurement for aircraft landing operations” (SRI, February 1972) or by R. H. Kohl,“Discussion of the interpretation problem encountered in single wavelength lidar transmissometers”, J. Appl. Meteor. 17, 1034 (1978).

Scanning the atmospheric volume by a 3-dimensional (3D) scanning lidar as proposed in the present invention, a_aer(k,R) and/or p_aer(/.,R) can be retrieved by making measurements at sequential incremental steps of the order of 1° or less, moving downwardly from the vertical towards the horizontal direction.

In this way, RF at Equations (2) and (3) denotes the distance at which we abstract 1-2 slant - height bins (1 bin is equal to the spatial resolution AR of the lidar device) for every new measurement and we retain the last value of a_aer(k,R) and/or or p_aer(7,R), at which RF was previously taken, as the new calibrating values for a_aer(k,R) and/or p_aCr(k,R). So, if RF, a_aer(7,R), Paer(^ R) are the last known values and RF-_{I W}, aaer(l,Rnew), b_3bί(l^_he\n) are the new ones, then we have:

R_F-„_{e, < Rr} (10)

(where N_bin=l , 2, 3, etc.) (11)

In this way, from each sequential lidar measurement the corresponding a_aer(7,R) and/or P_aer(^R) values are retrieved, taking into account the previous (last-known) values of these coefficients and by reducing RF at the same time, there is a continuous“calibration” of these values, from the previous to the next“lower” position of Ri (i=l, 2, 3, etc.). This stepping procedure is illustrated in Fig. 1.

In this way, it is proposed by the invention, for every height“drop” towards the horizontal position, to derive the values of a_aer^,R) and/or b_3eG(l^) from the vertical position to the final horizontal one, thus avoiding the ambiguities related to assumptions, like presuming homogeneous atmospheres, which produces large errors, especially over long distances.

However, it is herein acknowledged that the presently proposed technique of continuous recalibration also assumes a certain atmospheric“homogeneity” based on the a_aer(/-,R) or b_3eG(l^) values, between successive equidistant atmospheric layers (by moving at the incremental steps of the order of 1°, each successive atmospheric layer having a depth of 100 m over a horizontal distance of about 5.9 km. The limit of 1° is user adjustable and more precise measurements may be received if a lower angular value is used, especially depending on the instability of the weather conditions. It is herein noted that the more unstable or “heavy” the weather conditions are, the lower should be the limit of incremental steps of measurements.

Methods like the slope technique or the ratio method, as proposed by V. A. Kovalev, W. E. Eichinger,“Elastic Lidar,” Wiley Interscience (2004) and by R. H. Kohl,“Discussion of the interpretation problem encountered in single wavelength lidar transmissometers”, J. Appl. Meteor. 17, 1034 (1978), respectively, mainly use the assumption of homogeneous atmospheres over longer distances, to provide slant range atmospheric parameters. However in this way, the error on the derived parameters increases due to the instability of the prevailing atmospheric conditions, as a function of distance R. Other methods like the multiangle one, assume that a horizontal constant value of b_3et(l) is present at any distinct height, which is not the case for“real” atmospheres.

Moreover, the Optical Depth Solution according to V.A. Kovalev assumes that the aerosol backscatter to extinction coefficient is constant and the optical depth must be estimated by other independent measurements, in the vertical direction, like from a solar radiometer, or by a Raman lidar as discussed by A. Ansmann, M. Riebesell, U. Wandinger, C. Weitkamp, E. Voss, W. Lahmann, W. Michaelis,“Combined Raman elastic-backscatter lidar for vertical profiling of moisture, aerosol extinction, backscatter, and lidar ratio,” Appl. Phys. B55, 18 (1992). Although this method seems to work well under different atmospheric conditions, the problem becomes noticeable when trying to make slant-horizontal range measurements, where this value cannot be retained as constant and cannot be found, especially in horizontal measurements and in longer ranges, where the lidar signal may become too noisy.

Furthermore, the aforementioned technique of slant measurements of the invention maybe compared to the Boundary Point Solution (BPS) as described by G. Pappalardo et al. A. [“Aerosol lidar intercomparison in the framework of the EARLINET project. 3. Raman lidar algorithm for aerosol extinction, backscatter, and lidar ratio”, Appl. Opt. 43, 5370 (2004)]. The BPS solution assumes that the aerosol backscatter to extinction coefficient is constant and range-independent and sets the extinction coefficient as a known value at a specific range (boundary conditions). The same technique was used by Klett referred to hereinabove and this setting of initial conditions is the only similarity with the technique proposed in the present invention wherein, after the first vertical lidar measurements and the use of DENOISING algorithms, the reference (maximum) height is calculated, such as the values of a_aer(/_^,R) and/or b _aer^,R) bcCOfflG Z6G0. Then, it IS Gcisicr to retrieve t_aer(^,R) c Tld/qG aer (W for lower heights from the lidar signal.

So, with the technique of slant measurements of the invention and using the above mentioned equations 2 and 3, we expand this way of thinking, for slant and horizontal range measurements with different lidar ratio values, each time, taken from the “Weather Phenomena” algorithm for each measurement of interest, dependent on the type of the desired detection (e.g. thin fog, moderate fog, haze, etc.).

Regarding the visibility slant range measurements, according to WMO and ICAO rules for daytime visibility and by using the empirical and well qualified, Koschmieder law (at 550 nm), the visibility Vis (in km) versus the atmospheric extinction a (in km ¹) is given by:

Vis= 3/a_aer (X,R) (for 8=0.05 at 550 nm) (12)

Vis = 3.9\2]a_aer{X,R) (for 8=0.02 at 550 nm) (13)

where, e is a pure number, showing a contrast threshold, as a difference of the self-luminance of any object and the general luminance of the area viewed from a standing position. As reported by W. Vieeze, J. Oblanas, R.T.H. Collis,“Slant range visibility measurement for aircraft landing operations” (SRI, February 1972), the contrast threshold of 0.02 gives superior results to that of the value of 0.055.

Therefore, the estimated error of the retrieved Vis values directly depends on the accuracy of the retrieved values of a_aer(Z R) and/or p_aer^,R). When the classic Klett’ s inversion method is implemented, where a constant C(>_-5R) value is used as input (constrained to the mean aerosol optical depth (AOD) obtained from a nearby sun photometer, the a_aer^,R) values are retrieved with an uncertainty of the order of 20-30%. The Raman technique (during nighttime) retrieved the a_aer(7,R) vertical profiles with uncertainties of ~5-l5%.

The error estimation for the technique of slant measurements of the invention at the initial (vertical) lidar pointing measurement, concerning the retrieval of a_aer^,R) and/or p_aer(/-,R) depends only on the error estimation of the method used (e.g. Klett or Raman). In the case of a non-homogeneous atmosphere as assumed in the present invention and using Klett’ s method, we assume that we can obtain the values of C( ,R) and the ranges of a_aCr( ,R) for different types of aerosols (thin fog, haze, light haze, etc.) from Table 1 hereinabove. In this case we have overcome the problem of the C( ,R) uncertainty for these types of atmospheric conditions and we are therefore capable to obtain the profile of a_aer(Z-,R) and/or p_aei ,R) without any error related to the assumed value of C (Z.,R).

Then, if the calculation of visibility is of interest, we can derive its uncertainty (AVis) in the value thereof by taking the derivative of the visibility as derived from Equations 12 or 13: AVis = Vis [lsa_aer{X, R)/ cc_aer (l, /?)] (14)

where Aa_aer denotes the estimated error value of extinction coefficient. TABLE 2

The above Table presents various methods and techniques used to derive/estimate the atmospheric extinction and/or visibility and their associated uncertainty, compared to the 3D Stepping Technique of the invention that is applicable to either homogeneous or non- homogeneous atmosphere. This advantage and the fact that the error on the a_aer( R) and/or P_aer( R) estimation (as well as on the Vis) of the technique of the invention is approximately zero makes this technique more attractive in operational-commercial lidars and generally, in atmospheric visibility devices. It is an advantage of the 3D Stepping Technique that it may use 3D Lidar measurements in operational environment and perform slant and horizontal measurements without any theoretical assumptions with perfect accuracy for atmospheric parameters measurements.

DENOISING-l/2 ALGORITHMS

At this point, once we have calculated the RCS from equation (6), a denoising filter is applied to effectively reduce the noise of the lidar return signal and to automatically estimate the reference calibration height. Although various denoising techniques may be used to this end, two filters of the invention are used and proposed herein below.

Using the first filter, DENOISE- 1 and for distances typically less than R=5 km (R is a user adjustable parameter according to atmospheric conditions and power transmittance), a median spatial filter is been applied by means of ±2 range bins, as a first step. Then, at the second step, the values (at bin number n, counting from the beginning of the lidar signal) of a_aer (Z.,R_n) and/or p_aer(/.,R_n) are being checked, if they are equal or greater than a threshold (k= 1.1- 1.15 for _aer (^R) values and k=l.15-1.2 for a_aer (^R) values) from bin to bin (Eq. 13) holds also for the a_aer (l, R_n) values).

UW// ^A₊₁) ³ or /UA,/U// W ³* (¹⁵)

If one of the conditions in Equation 15 is valid, then, a_aer^,R_n) and/or p_aer( ,R_n) retain their initial values; the corresponding values at the next bin n+l (e.g., a_aer (Z.,R _n+i) and/or a_er( ,R_n+i)) are been tested (step three), if Equation 15 is valid considering k=l,5-l,55 for p_aer(/.,R) values and k=l .55-1.6 for a_aer(/-,R) values, from bin to bin. If the latter is the case then the values of a_aer(^R_n+i) and/or p_aCr(/.,R_n+i) are been set to zero.

Otherwise, they take the value of: U^ / W+/U^A ,)]/2 (16)

In any other case, there is no value change for a_aer (L,R_n) and/or p_aer(k,R_n). For instance for greater distances (i.e. R>5 km) the same procedure is being followed, but with ±3 bins. Afterwards (step four), a subtraction is performed of the negative values of a_acr(k,R) and/or aer(k,R) through the whole signal, because these values are clearly artificial and as that they are no longer needed. If:

< W = 0 and/or U« = 0 (17)

At step five, we set a threshold of a_aer (k,R_n)>60 km ¹ and/or b_3a (k,R_n)>6 sr 'km ¹ or _aer ( ,R_n)>0.6 sr 'km ¹ (typically, the a_aer(^R) values are one to two orders of magnitude higher than the values of _aer(k,R_n)) (V. A. Kovalev, W. E. Eichinger, “Elastic Lidar,” Wiley Interscience (2004)). In this case the signal is typically considered too noisy, and those values are set equal to zero.

Afterwards, at step six, the DENOISE-2 algorithm can be applied. We first divide the measurement range R in equal“distance packages”, each one ranging from 300-500 m up to 1000 m. Then, at step seven, the length of the N-th signal“string” (LSTN) (starting to count preferably from the max range (R_max), chosen by the user, towards the lidar position) within a set of“distance package”, is calculated from the following Equation.

where, å is the sum, n represents the corresponding range bin (from range b_min to b_max), LST_n is the smallest“string” length of the signal of its bin, d, is the smallest difference and A_n represents the corresponding P'(k,R_n) or a_aer(k,R_n) and/or p_aer(k,R_n) values. Finally,

is the full length of the signal, as shown in Fig. 3, throughout the whole lidar signal (typically up to R_max=7.5-l0 km from the lidar position). At step eight, a ratio calculation is made according to the following Equation:

(LST_N)/(LST_N < 1 + KG" (20)

where, sequential LSTs are calculated and compared to a constant value of ~l (e.g., 1±10 ^h).

If Equation (20) is valid then (step nine) the values of P'(k,R), a_acr^,R) and/or p_aCr(/sR) corresponding to LSTN+I, are no longer useful and are therefore set equal to zero. The idea behind this assumption is that, when large distances (e.g. > 8-10 km) are considered (i.e. the signal to noise ratio-SNR becomes smaller than 1-2 the received lidar signal becomes extremely noisy as acknowledged by Heese, B., Flentje, H., Althausen, D., Ansmann, A., and Frey, S., in“Ceilometer lidar comparison: backscatter coefficient retrieval and signal-to-noise ratio determination”, Atmos. Meas. Tech. 3, 1763 (2010).. Therefore, the LST values between adjacent equal LSTs become very close to each other and then noise“governs” the lidar signal, so the retrieval of a_aer(k,R) and/or p_aer(k,R) becomes useless and the associated values of a_aer(/.,R) and/or p_aer(/ ,R) to the LSTN+I are now of no interest.

If Equation (20) is not valid anymore, then the procedure stops, and the respective values of An and of LSTN+I are retained. The same procedure shown at steps seven, eight and nine can be re-applied by setting the start of the signal at the lidar position, and N is, then, counted from that point towards R_max. However, it is herein proposed to start N from Rmax towards Rmin (the location of the lidar position), because in this manner more reliable results are being retrieved, especially under noisy conditions.

The application of the herein proposed DENOISING algorithms shows that in most real- operational conditions cases (more than 50 cases have been tested), these algorithms manages to retain (in more than 85% of the studied cases) the useful (a_aer^,R) and/or b aer (k,R)) information values“hided” inside the noisy signal. These algorithms are found to work better in noisy signals (lower SNR values) and can be used in combination or separately. An example of application of both algorithms is presented in Fig. 3 where it is noticeable that useful _(Xacr(^R) information is been kept for further processing.

The algorithms can be applied, already, from the first vertical measurements, to the P'( ,R) values, using the BPS method (near-end or far-end) and automatically select R_ref, following the above mentioned procedure, where a_aer(L,R) and/or b_3a-(l^)) are close or equal to zero as observed from the longest distance of the lidar position. In this way an automatic selection of R_max is achieved for each lidar signal, which is recalculated each time according to the 3D stepping technique of the invention, by subtracting 1-2 range bins for every next measurement.

The value of P '^,R) remaining at that range bin (R_n= Rref preferably equal to zero) is then being chosen to calibrate the next slant measurement, with the difference of using P'(k,R) values instead of a_aer (k,R) and/or b_3a- (k,R). In this case, k, at step two, takes the value of 1.3 and at step three takes the value of 1.5. Then the algorithms proposed can be rerun using these P^'(7,R) values with similar results to those produced by using the a_aer (7,R) and/or b_3eG (k,R). In any case, the recommended values of k in this paper are a good approximation for “real” atmospheric conditions at mid-latitude regions and at the wavelength of 355 nm.

The diagram illustrated in Fig. 3 shows results of the application of the above DENOISING algorithms of the invention on the values of a_aer (^R) retrieved at 355 nm by the LRSU- NTUA (03-10-2016 at 08:02:10 UTC). The lighter grey line denote the retrieved a_aer( ,R) and the darker line is the denoised value of a_aer (k,R) after the application of the“DENOISING 1/2” algorithms.“Distance packages” are chosen, typically equal to 997.5 m for each LST, for a range resolution of 1 bin (7.5 m). It is a benefit of the DENOISING 1/2 algorithms that they afford identifying and subtracting the noisy part of a signal whilst keeping the real existing measured value at the same time, without“hearting” the valuable data contained therein.

THE VISIBILITY (VIS) ALGORITHM

According to ICAO and WMO rules, the visibility is defined as the greatest distance at which a) a black object of suitable dimensions, situated near ground, can be seen and recognized when observed against a bright background, and b) the greatest distance at which lights in the vicinity of 1000 candelas can be seen and identified against an unlit background (Aerodrome Meteorological Observation and Forecast Study Group (AMOFSG), AMOFSG/10-SN No. 11 (2013)). However, these two distances have different values in air of a given extinction coefficient, and the latter abovementioned distance (b) varies with the background illumination. Moreover, the former distance (a) is known as Meteorological Optical Range (MOR).

Nowadays, visibility measurements are provided to the airport authorities through different measurement devices and techniques. More specifically:

1) Forward scatterometers : they are forward scattering devices that measure the atmospheric scattering coefficient at the same location using light beams transmitting and receiving units; they are, usually, deployed along the runways and are able to provide visibility in a distance greater than 10 km of the device’s location. They also provide Runway Visual Range (RVR). Their disadvantages are that they assume homogenous atmospheric conditions (so they need constant calibrations) and measure scattering coefficient at ground level and not in slant ranges (SR);

2) Transmissiometers: they are also forward scattering devices that provide measurements of the atmospheric extinction coefficient (a). Their transmitter and receiver units are spaced at a certain distance and can provide visibility and RVR data at distances greater than 10 km of the device’s location. Their disadvantages are that they also assume homogenous atmospheric conditions and are also located at ground level, thus not being able to provide visibility data along the slant direction along which an aircraft approaches for landing.

3) Scintillometers these devices measure the sensible heat flux over long distances and so, in some way they can measure the atmospheric visibility; again, they are limited in ground measurements and cannot provide slant range or even vertical visibility measurements.

4) Ceilometers : these devices, employing low power laser beams, work on the same principle as the lidars and can use 3-dimensional scanning techniques. They measure the atmospheric backscatter coefficient (b) and can provide cloud base and ceiling, as well as vertical visibility. They provide the vertical aerosol layering and a value equivalent to the visibility; only, recently, the 3D scanning ones are able to provide atmospheric layering and visibility in a slant range for the approaching aircrafts.

5) Distrometers [3]: These devices use a laser beam, in situ, to measure the precipitation and also the presence of haze, fog etc. Through these measurements and based on the WMO and ICAO rules, they can provide“locally” visibility data.

6) Telephotometers / Cameras: These devices use the contrast of daylight gathered from the sky in comparison with that of the runway. They work only in foggy days and locally at the

Airport sites and specifically along the runways.

7) Nephelometers: These devices use a light source (LED or laser) to measure, in situ, the aerosol backscatter coefficient, and then, the atmospheric visibility can be retrieved using certain assumptions. Again, this technique provides measurements only in situ, whilst slant range visibility measurements are not available.

8) Cellular networks : They can be used by means of atmospheric condition detection by attenuation in a line of transmittance between transmitter and receiver at the range of 20 GHz to 38 GHz (until the range of 120 GHz).

9) Weather Radars : These devices can provide the atmospheric visibility, under cloudy/foggy conditions, based on the backscatter coefficient from hydrometeors, through radar reflectivity data.

10) Last, but not least, the old traditional way of measuring visibility, visually, by a trained aircraft controller or an experienced weather observer-meteorologist. This old fashioned way is indeed a safe way of revealing the actual visibility at an airport; however, this method is not automated, can be done only by trained personnel and can lead to large uncertainties.

In accordance with the present invention, following use of the abovementioned denoising filters and the choice of the slant range measurement technique, the VISIBILITY Algorithm is used to provide visibility measurement both to the airport tower controller in the direction of an aircraft approaching for landing and to the pilot of the approaching aircraft in the direction of the line of approach of the aircraft to the runway, according to WMO and ICAO rules of daytime visibility, at the airports. Using the values of a_aer(U,R) and p_aer(L,R) as provided in Table 1, and equations (8), (9) the AOD (Aerosol Optical Depth) and the total AOD can be calculated both from the lidar-tower point of view towards the maximum range of the lidar signal and from the maximum range of lidar’s signal, towards the lidar device near the runway (pilot’s point of view).

In addition, the daytime visibility (Vis), according to WMO and ICAO rules, can be estimated from the empirical and well qualified, Koschmieder law (at 550 nm), which connects the visibility (in km) with the atmospheric extinction coefficient a_aeilA,R) (in km ¹), as provided in Equations. (12) and (13).

By using the visibility versus the atmospheric extinction coefficient at 550 nm (cf. Fig. 2c), we are able to have a good approximation of the above required values at daytime, at 355 nm as described hereinabove.

In particular, after the physical calibration of the lidar device, as well as the application of the 3D stepping technique, the denoising of the lidar signals and the application of the Weather Phenomena Algorithm and its C_n^,R) values, the VISIBILITY algorithm is being applied. In order to acquire the visibility values (Equations 10 and 11) the Visibility algorithm calculates the value of r_totai using Equations (8) and (9) (favorably using the trapezoid method) from the lidar’s position to its pointing direction (typically 15° from horizontal). It checks at what range (Ri) of the lidar’s position the Xtotai(0,Ri) > 2*l0 ², which according to Fig. 2b is the value of the extinction coefficient needed for the detection of cumulus cloud, with no visibility through it and corresponds to fiacr (¾,R) > 2*10 ³ nT¹ sr ¹, with C(X,R) = 10 sr and oi_aeifA R) > 2 * 10 ² m ¹.

Then, the Visibility algorithm stops working at that range (R_t) and retains that range as the total visibility range at that time, according to WMO and ICAO rules for daytime visibility. As for the Visibility algorithm application, the AOD starts to be measured from the maximum range (Rmax) or from the range at which the lidar signal has useful information containing detection of the aerosols (R₂), towards the lidar’s position.

From that automatic calculated or manually selected range (Rmax or R₂ - user defined) and when Ttotai(Rmax, R3) > 2*1 O ² or r_totai(R2, R3) > 2*l0 ² (for the reason mentioned before), the visibility is equal to that range (R₄) where:

(21)

or R₄ = R₂ - A (user defined) (22)

counting from the range of Rmax or R₂ (user defined). In both cases, if r_totai(0,Rmax) £ 2* 10 ² through the whole distance (Rmax), then the visibility is equal to Rmax. The outcome for Visibility accordingly, is that:“The visibility of tower controller is Ri (of Rmax)” and“The visibility of the pilot from R_max (or R₂), is R₄ (or Rmax)”.

Another way to accomplish visibility measurements is by using the below mentioned Equation (23). Thus, taking available lidar measurements that have been made at certain locations at specific wavelengths, categorized by season and atmospheric composition by using the Absorption Angstrom Exponent (AAE) (which determines the wavelength dependence of absorption of aerosols) and the Extinction Angstrom Exponent (EAE) values (which is an indicator of the particle size of the atmosphere, together with significant meteorological conditions, like the relative humidity (RH), temperature, rain, etc., and the geographical location, one is able to acquire visibility measurements. The idea is that by performing lidar measurements at different wavelengths over different locations, one would be able to acquire the visibility and meteorological conditions with good approximation, at a single wavelength (i.e. at 355 nm) (cf. Equations 23, 24, 25), deployable at any airport (location):

a(2,)/«(A) = /?(W(A) = (23)

in order to be able to acquire:

/?(4) = */?(A₂) (24)

a(l_i ) = 1 * a(l₂) (25)

where, li, l₂ are different wavelengths, for example li=355 nm, where measurements have been taken and l₂=550 nm, at which the visibility and categorization of the atmospheric conditions have been described by Measures; L is a constant that is defined by different values, according to the calculations made at different seasons, chemical compositions of the atmosphere, at any geographical location using Equations (23), (24) and (25). Ambiguities about aerosol composition or mixtures that result at some AAE values are able to be reduced by clustering, as proposed by Russel et al. [“Absorption Angstrom exponent in AERONET and related data as an indicator of aerosol composition,” Atmos. Chem. Phys. 10, 1155 (2010)].

Figs. 4a-4d illustrate examples of use of the“weather phenomena” and“visibility” algorithms of the invention.

In Figs. 4a-4b a case study of the application of the“weather phenomena” and“visibility” algorithms is shown, where the airport tower visibility from the lidar’ s vertical pointing direction is found to be 523 m and the pilot’s visibility from a 3000 m height (distance) is of 1928 m. Fig. 4a in particular shows a diagram of altitude versus time depicting an actual RCS signal acquired with the lidar device of the invention at 1064 nm by the LRSU-NTUA (16- 05-2011). Further, Fig. 4b shows a diagram of the type of atmospheric layer versus range, wherein“VISIBILITY” and“WEATHER PHENOMENA” algorithms of the invention have been applied with the scaling from darker to lighter color indicating variation of the type of atmospheric layers being scanned with the 3D lidar pointing to the vertical (355 nm) at 11:30:10 UTC, starting at 417.5 m above ground level up to 3000 m. The Black and white visualization in this diagram (a visualization in color is produced in practice) provides visualization of meteorological conditions with the darker color indicating“strong” blurriness of the atmosphere starting from Cumulus Cloud, Moderate Fog, light Fog, Thin Fog, Haze, whereas lighter or no color (white) indicates a“Very Clear” or“Exceptionally Clear Sky”. Figs. 4c-4d present, with diagrams analogous to those shown in Figs. 4a-4b, a further case study of the application of the“weather phenomena” and“visibility” algorithms of the invention, wherein an actual RCS signal is acquired with the lidar device of the invention at 1064 nm by the LRSU-NTUA (31-10-2011 at 15:29:50 UTC ). In this case, the pilot’s visibility from the vertical height of 3000 m is 1370 m and tower’s visibility is 1512 m because of the strong backscattering (probably cloud) at this height with no visibility through it.

It is a benefit of the Visibility algorithm of the invention that it provides visibility measured in 3D through accurate measurements made by the 3D Stepping Technique, both for the Tower controller and for the Pilot point of view, including Slant Range Visibility (SVR) and Horizontal Visibility (HOV) with accurate Lidar measurements. ATMOSPHERIC LAYERING AND DISTRIBUTION

THE MUSTI-L/D ALGORITHM

The algorithm Multiple Space-Time filters (MUSTI) is able to provide aerosol Layering (L) and their spatial Distribution (D) based on the retrieved b_3a(l,K) coefficient (and/or a_aer^,R), RCS) using a 3D scanning lidar. When L and D is retrieved, one is able to measure the wind velocity in 3D by continuously measuring the p_aei ,R) coefficient (and/or aaer(L.R), RCS) profiles using a technique similar to the one proposed by Tomas and Rocadenbosch [“Wind retrieval from multiangle backscatter lidar profiles through anisotropic aerosol structures”, J. Geophys. Res. Atmos., 120, 7758-7776 (2015)]. Alternatively the“SIBESMEA” algorithm proposed herein below can be used. The algorithm MUSTI L/D is adapted to run for generating a quantitative report of atmospheric layering. This algorithm is also necessary in the processing of“MUSTI-PBLH” algorithm that is presented herein below.

Basically,“MUSTI-L/D” Algorithm uses a combination of the backscattering signal and its coefficients retrieved (p_aer^,R) and/or a_aer(/,R), RCS) values spectrum, certain space arrangement settings (‘space filters’) and first and second derivatives of p_aer(^R) and/or a_aCr(L,R), RCS and changes thereof through time. So, with the use of“space filters” being proposed and the effect of changes in the course of time, in the backscattering signal and its coefficients p_aer(L,R) and/or a_aer(L,R), and/or RCS, thereupon, the invention proposes means of capturing varying atmospheric layers.

In“MUSTI-L/D” algorithm, p_aer(L,R) (and/or a_aer( ,R), RCS) in the vertical or slant range (R) after range calibration of our lidar system, calculation of RCS (cf. Eq. 6) and denoising of the lidar signals using the DENOISING algorithms described herein above. The whole process can be made upon the whole received lidar signals (Eqs. 1, 2, 3, 4, 5). In this work, we chose to use settings appropriate for b_3eG(l,B) (Eq. 3). However, the way of thinking of the algorithms remains the same, with different setting values, for different kinds of measurement types, as analyzed below.

In accordance to the invention, the calculated values of b_3e,(l, ) (Eq. 3) are being“cut” into smaller“space filters” in order to create distinct conditions for processing. By way of example, 12 distinct space filters (G1 to G269) from 7.5m (1 bin) to 2017.5m (269 bins) are being used.

G3 = 3*bin=3* 7.5( ) (26)

where, G3, represents the“space filter” of a package of 3 bins.

In general, a plurality (n) of space filters Gl, G2, G3, ....Gn, is being created, each space filter representing a predetermined integer multiple of bins, wherein each bin is identical to the resolution of the lidar device being used.

Afterwards, we check if the calculated b^l,B) in any of these filters is greater than 8* 10 ⁷nT ¹ sr ¹ and C^,R) = 50sr (a_aer( ,R)>4* 10 ⁵m '), because a value greater than this, means that atmospheric phenomena like haze, fog, etc., are present, at the wavelength of 355nm, according to Table 1 herein above. Then, we check if the values of b_3eG(l,B) in any of these space filters from a package of bins to the next package of equal number of bins, are greater than 10% (percentage produced statistically in a great number of measurements), in order to find out if the difference between the b_3et(l,B) values, is within normal oscillation of the coefficient of the same layer or has to do with some important change that denotes different layering:

/U^ / ^W, > 0.1 (27) or

/U^ ₃₊₁ //U^,¾₃ > 0.1 (28)

where, p_aer( ,,R)_G3 is the value of p_aer(k,R) of the space filter of the package of 3 bins and p_aer(k,R)_G3_+i is the next value of p_aer(k,R) of the next‘space filter’ of the package of 3 bins also.

At the next step, the algorithm calculates the first order derivatives of all space filters, at any distance, (dp_aer(k,R)G3/dR) and thereafter calculation is made in order to produce the second order derivatives thereof (d²p_aer( ,R)G3/dR²).

M3DR = {d( _ae ,R)₀₃)/dR} (29)

M3DDR = {d²(R_aer( ,R)_G3)/dR²} (30)

M3DRP=M3DR> 0.212 (31)

M3DRN=M3DR< -0.212 (32)

M3DDRJ^M3DDR> 0 (33)

M3DDRN=M3DDR<0 (34)

where, M3DR is a matrix with elements, the values of the first order derivative of the package of 3 bins, at range R, from matrices and M3DDR is the matrix with elements, the values of the second order derivative, of the same package of bins, at R, from matrices, M3DRP is the matrix of the positive values of M3DR, M3DRN is the matrix of the negative values of M3DR, M3DDRP is the matrix of the positive values of M3DDR, M3DDRN is the matrix of the negative values of M3DDR.

Herein the term“agreement” is designated to declare existing values at the elements of the corresponding matrices, as a function of range, with values other than zero. If so, the new element of the new productive matrix, takes the value of 1 , when terms of“true” conditions are met, or, 0, when“false” conditions are encountered. All these matrices are being used in all kinds of space packages. By way of example, M169DR is the matrix with elements of the first order derivatives of the package of 169 bins, at range (R) and M55DDR is the second order derivative of the package of 55 bins at range R, etc.

Then, we check if these first order derivatives have the same sign (positive or negative) for a number of continuing space filters like G2, G3 and G5. So, if the I^st order derivatives with the same sign correspond to values greater than an angle change of 12° (statistically produced for estimating serious change of layering, i.e. a condition that cannot be associated to a layering normal oscillation-user dependable), we store these values at distance (R_x). A “beginning” of a layer is produced by positive values of first order derivatives (>12° angle change of values of p_aer^,R) - recommended) and positive values of 2^nd order derivatives, plus 10% change of averaged p_aer(k,R) values, from the next package of equal number of bins. Accordingly, an“ending” of a layer is produced by negative values of first order derivatives (>12° angle change of values of p_aer(k,R) - recommended) and negative values of second order derivatives, plus 10% change of averaged p_aer^,R) values from the next package of equal number of bins.

Fig. 5a presents a diagram of altitude (height) versus time depicting an actual RCS signal acquired with the lidar device of the invention at 1064 nm from LRSU NTUA data at 12:29:40 UTC on 01/02/2016. Accordingly, Fig. 5b presents p_aer(k,R) values versus range (R) at the time of receipt of data from the RCS signal of Fig. 5a. Then, Fig 5c presents beginnings and endings and Fig. 5d types of atmospheric layers versus range (R) at the time of receipt of data from the RCS signal of Fig. 5a as provided by MUSTI-L/D algorithm of the invention in use at a vertically oriented 3D lidar, operating at 355nm, starting at 417.5m above ground and terminating at more than 3000 m above ground level. Fig. 5c indicates‘beginnings’ of an atmospheric layer with a darker color and indicates endings thereof with a lighter color, in agreement with Fig. 5b.

Fig 5d shows a black and white visualization of the type of atmospheric layer versus range (R) at the time of receipt of data from the RCS signal of Fig. 5a. Use of the“Weather Phenomena” algorithm is made to characterize the type of layering and therefore a darker color indicates“strong” blurriness of the atmosphere starting from cumulus cloud, Moderate Fog, light Fog, Thin Fog, Haze. No color (white) indicates“Very Clear” or“Exceptionally Clear Sky”.

It is noted that two or more continuing identically shaded (darker or lighter) lines, without any different color line in between, indicates that p_aeit7i,R) (and/or a_aer^,R), RCS) continues to go higher (darker line) or lower (lighter line) in quantity (value), probably in a different kind of layer. This visualization provides a convenient understanding of atmospheric layering and distribution (L and D).

Common layering altitudes - heights are being tracked among equal and different continuing space filters and an average is being produced for these altitudes, in terms of 10 to 15 range - height bins (75 to 107.5m spatial tolerance). In this way we are able to produce, in total, atmospheric layering (MUSTI-L) and we are able to have a clear understanding of the distribution (MUSTI-D) of p_aer^,R) (and/or a_aer(7,R), RCS) in a 3 dimensional way surrounding a 3D lidar system or an area, for example of an airport.

MUSTI-L/D IN PROCESS

The algorithms were developed and a compound of matrices was created and used for the outcomes to be produced. All the matrices below present the corresponding valuable information with range (R) - height, accordingly for slant or vertical measurements, in particular:

(i) At first we take a recording file of backscattering (p_aer^,R) (and/or a_aer(/ ,R), RCS - user depending, but the value settings here responds to p_aer(^R)) with R, (±) 30° from vertical pointing, for grounded lidar system and vertical, or near vertical, measurements, accordingly (Eqs. 3 and/or 4), in order for the lidar system to reach molecular heights backscattering signal at the receiver. We start these measurements from a typical height of 2l0m above our lidar system, due to overlap issues.

(ii) Then we create two different matrices, B and C, where the elements of B are separated in number of total bins and has the values of p_aer^,R) as elements, that correspond to each bin (lbin=7.5m spatial resolution) and the elements of C has the actual range (R) in number of bins (or calculated in m):

B = {fi_aer{ ,R)} (35)

C= {£} (36)

(iii) Then, we create space filters by smoothing every 2, 3, 5, 9, 15, 27, 41, ... 269 bins the values of the elements of matrix B and we put these calculations into elements of new created matrices, like M2, M3, M5... M269. In this work we created 12 space filters (user dependable) and it is noted that, space filters has to do with the creation of ‘values of measure’ - conditions of measurement, for our processing algorithm: M3 - {G3,R} (37)

(iv) Moreover, we check, for every element of these matrices (M2, M3...) if any change of P_aer^,R) values from element to element of each m matrix, is greater than 10% and store these values at the respected R_x. It is noted that this step has to do with the monitoring of the values of p_aer( R) which denotes if the change is according to normal oscillation of the values of p_aei(^.,R) at a specific layer or has to do with the incident of a meteorological condition change, like the‘beginning’ or the‘ending’ of a cloud, haze, type of fog or other layer. We keep track of these values by means of R_x, for farther processing, to matrices H2, H3, H5 ... H269, respectively:

HS = {MS,[fi_aeJ ,R)_Gi / b_ae l,K)_aM > 0.1]} (38)

or

H3 = {M3,[p_aer{ ,R)_CM lfi_aer( ,R)_Gl > 0.\]} (39)

(v) We check which values of b_3eG(l^) of the matrices H3, H5 ... H269 are greater than

> 4*l0 ⁵ m ¹), where, according to Table 1 hereinabove, is the beginning of the meteorological conditions at 355 nm, upon which this algorithm is based. If these values of the elements of the matrices, are less than 8*l0 ⁷ m ¹ sr

', they are excluded from matrices H3, H5, ... H269. At this step (v) of the algorithm we can have the same results by setting values different than the abovementioned values of (P_aer^,R)>8*l0 ⁷ m ¹ sr ¹, C^,R) = 50 sr and a_aer(9sR) > 4*l0 ⁵ m ¹ at different wavelengths of our lidar system making use of the correlation process using Fig. 2a described hereinabove.

The values of p_aer^,R), less than 8* l0 ⁷ m ¹ sr ¹ and C(XR) = 50 sr (a_aer^,R) > 4* l0 ⁵ m ¹), seem to denote‘Sky Clear’ with a visibility of 6-10 km;‘Very Clear’ or‘Exceptionally Clear Sky’ with visibility greater than 10 km (at 355 nm). Otherwise, if a_aer(^R) is used, Eq. 2 must be taken into account and proceed to next steps of the algorithm according to correlated data obtained from Fig. 2a and keeping the corresponding value of C( R).

(vi) At this step, we create matrices of first order derivatives for every space filter from (M) matrices like M2DR, M3DR, M5DR ... M269DR (Eqs. 29, 31, 32) which values correspond to differences in angle change greater than 12° (statistically produced value avoiding normal oscillations of a specific layer and according to specific space filter, the difference can be varied in a way where the longer the space filter, the less the required change angle (> 10°) - user adjustable).

(vii) Afterwards, we separate these results of the last matrices to positive and negative first order derivatives. In this process, we put the positive first order derivatives at M2DRP, M3DRP, M5DRP ... M269DRP matrices and the negative ones at M2DRN, M3DRN, M5DRN ... M269DRN (Eqs. 31, 32).

(viii) We create different kind of matrices where positive and negative first order derivatives are divided, at DP2, DP3 ... DP269 for positive first order derivatives which are in agreement with H matrices, by means of values other than zero, at the respected space package of bins and DN2, DN3 ... DN269, which stand for negative first order derivatives (Eqs. 33, 34) with the same way of thinking:

DPS = {MSDRP,HS} (40)

DNS = {MSDRNHS} (41) where, DP3 (used commonly for“beginnings” of a layer) is the matrix that includes M3DRP and M3DDRP and DN3 (used commonly for“endings” of a layer) is the matrix that includes M3DRN and M3DDRN matrices which are in agreement with H matrices, by means of values other than zero at the respected space package of bins.

(ix) We create new matrices of second order derivatives for every space filter at the whole range (R_x) of (M) matrices named M2DDR, M3DDR ... M269DDR (Eq. 30).

(x) We separate the created M2DDR, M3DDR ... M269DDR at those with positive values, named M2DDRP, M3DDRP ... M269DDRP and at those with negative values, named M2DDRN, M3DDRN ... M269DDRN (Eqs. 33, 34), which are in agreement with H matrices, according to the same way of thinking like in the previous step.

(xi) We create the idea of“beginnings” and“endings” in matrices by creating M2S, M3S ... M269S matrices for ‘beginnings’ and include, DP3 ... DP269 as well as M2DDRP, M3DDRP ... M269DDRP matrices, in agreement, and‘endings’ who stands for M2F, M3F ... M269F and includes, DN2, DN3 ... DN269 as well as M2DDRN, M3DDRN ... M269DDRN matrices, in agreement:

MS = {DPi i^'DDRP_j (42)

MFi = {DM,M3DDRI (43)

(xii) At this point, all possible‘beginnings’ (M2S, M3S ... M269S matrices) and‘endings’ (M2F, M3F ... M269F) have been created. In the process of better clarification of a ‘beginning’ or an‘ending’, similar - close space filters products, in terms of range or height, must be in agreement. For example, in order for M3S to exist, M2S must exist at that range area of 3 bins and most probable M5S to, at the same range area, otherwise an oscillation of the layer probably exists, without denoting a different type of layer. The choice of the exact range (R_x) at which the algorithm denotes a‘beginning’ or an‘ending’ of a layer, is being made by the products of the space filter, of the same area of interest (area of a few meters), which corresponds to the greatest range (R_x)-distance from the location of the lidar system, for an‘ending’ (recommended) and the one with the smaller range (R_x) value-distance from the lidar system, for a‘beginning’ (recommended) is chosen. To better arrange these settings, the user must be familiar with the location of measurements, wherein for cold and clear layered atmosphere the products of the space filters should agree in terms of close range vicinity of nearby space filters as described, but for more disturbed weather, like the one existing in hotter climates and larger atmospheric layers, a wider range of space filters must agree in terms of existence. The above can also be decided automatically with the assistance of the DECIS (DECISion) algorithm that is presented and analyzed herein below.

(xiii) At this step we create new matrices named LAY matrices (LAY2, LAY3, LAY5 ... LAY269) and NLAY matrices (NLAY2, NLAY3, NLAY5 ... NLAY269) which stand for layering with LAY matrices denoting a layer with a“convex” configuration and NLAY denotes space between layers with a“concave” configuration. So, the LAY2, LAY3, LAY5 ... LAY269 matrices include the agreement, in terms of range, of DP2, DP3, ... DP269 and M2DDRN, M3DDRN ... M269DDRN matrices and accordingly, NLAY3, NLAY5 ... NLAY269 matrices, include the agreement of DN2, DN3,... DN269 and M2DDRP, M3DDRP ... M269DDRP matrices. It is typical for MS and MF elements of the matrices to surround the new produced elements, in terms of range, of LAY and NLAY matrices. In this way an isolation of existing layers is being made:

LA B = {D ,M2DDRN_S (44) NLA B = [DA3, M^,DDRP_f (45)

(xiv) In case of continuing MS matrix elements without any MF matrix element being encountered in between, in terms of range, it is highly recommended to use the MS element, that is most proximate to the location of the lidar system (at least in the vicinity of 10 to 15 range bins - user dependable), as the‘beginning’ of the present layer. Accordingly for the existing elements of the MF matrix without any MS matrix element being encountered in between, the MF matrix element with the larger distance (again in the vicinity of 10 to 15 range bins - user dependable) is recommended to be used as the‘ending’ of the present layer.

(xv) At this step, all the‘beginnings’ and‘endings’ that have been left, are been tested, if more than one range points (R_x) of a‘beginning’ or an‘ending’ accordingly, existing in close vicinity and with a spatial tolerance of 10 to 15 bins (depending on the accuracy the user needs). For values of LAY and NLAY matrix elements, different than zero, an average is being calculated in order to produce a clearer‘beginning’ or‘ending’ for the user and keeping track of this final layering at a final new matrix, that is called LAYER:

LA YER {LAY 2, LA 73,· -,LA 7269,

NLA Y 2, NLA 73,· · ·, NLA 7269}

(xvi) Finally, the results of this LAYER matrix are presented in different colors for ‘beginnings’ and‘endings’ in order to have a good visualization of the layering (L) of the atmosphere, as well as of the Distribution (D) of _aer(k,R), to lidar’s pointing direction (Fig. 5c).

(xvii) Layers, can also be colored as well, through the production of the optical depth (x) (Eq. 8) of each layer and a categorization can be made, according to Table 1 at any wavelength that the lidar system is transmitting (Figs. 2a, 2b, 4b). In this way, a layering can be made also, with the cooperation of the“Weather Phenomena” algorithm as illustrated in Figs. 5d. It is also practical for depolarization capability lidars to color the different kinds of aerosol at any distance and a profile of the aerosol layering depth of each kind, can be made, for reasons of being handy to use and visualize.

Obvious benefit of the MUSTI-L/D Algorithm is its Ability to perform quantitative Layering and Distribution of aerial masses also through a profile accurately and also in 3D.

THE PBLH

According to Stull ([An Introduction to Boundary Layer Meteorology”, Kluwer Academic Publishers, Dordrecht, (2013), p. 2) the Planetary Boundary Layer (PBL) is the lower part of the troposphere that is most influenced by the presence of the Earth and responds to surface forcing within a time scale of 1 hour or less. The PBL Height (PBLH) is defined as the highest part of the PBL, just below the so-called entrainment zone.

The usual systems used to retrieve the PBLH are radar backscattering [J. Klett,“Stable analytical inversion for processing lidar returns,” Appl. Opt. 20, 211-220 (1981)], lidar backscattering [S. H. Melfi, J. D. Spinhime, S. H. Chou, et al.,“Lidar observations of vertically organized convection in the planetary boundary layer over the ocean”, J. Climate Appl. Meteor., 24, 806-821 (1985)], ceilometers, High Resolution Doppler Lidars (HDRL), sodars, radio acoustic sounding systems (RASS) [International measurement confederation, Proceedings of the 2^nd Symposium on“Metrological assurance for environmental control”, Nova science publishers (1988)], radar wind profilers and radiosondes.

According to [R. B. Stull, “An Introduction to Boundary Layer Meteorology”, Kluwer Academic Publishers (1988) p. 21], atmospheric aerosols, which are mostly present at the lower troposphere are the best tracers of atmospheric motion, and thus, can be used to study the PBL structure. To trace atmospheric aerosols the best wavelength region lies between 1- 3 pm. According to the same paper the use of a multi-wavelength lidar would provide better results regarding aerosol tracing.

THE MUSTI-PBLH ALGORITHM

Taking into account the results of MUSTI-L/D algorithm and its LAYER matrix, usage is thereafter made of the so called‘time filters’ of the invention. In this case the goal is to find the PBLH. Basically, we keep track of the‘ending’ heights of layers through time, because PBLH is presented like an‘ending’ height at the backscattering signal and during a period of time of 15 minutes or more, but less than an hour in general, we are trying to locate the most common‘ending’ height and the outcome of“MUSTI-PBLH” algorithm is this height.

To begin with, we already have our lidar oriented in, favorably, vertical or near vertical position (it is recommended in less than 30° from vertical), for purposes of capability of the lidar system to reach molecular heights at the receiver. We continue measuring for at least 15 to 20’ (user dependable), in this vertical or near vertical orientation and acquire at least 15 or more lidar data files. It is time for the creation of the‘time filters’ of the MUSTI algorithms.

These filters have to do with the 15 or more data files that we are keeping track of. In other words, we keep track of the layering and distribution at the LAYER matrix, but this time, we check only the produced‘endings’. In this algorithm we treat PBLH as a layer‘ending’ with no‘beginning’ in every specific small time frame of 15 or more minutes, for vertical pointing measurements. The ideal would be half an hour of recording time of p_aer^,R), but in any case less than an hour, since it is within this period of time that PBLH is expected to have a small or medium change, especially in medium weather conditions. The idea is to take advantage of that unique characteristic of PBLH that looks like a common‘ending’ in the time scaling “MUSTI-L/D” algorithm.

Then, we check the heights of those‘endings’ and if they exceed a common PLBH area for that time of the year and for the location of measurement (colder or hotter climates with lower or higher commonly PBLHs, accordingly - user provided). We add (and subtract accordingly), some height from that area of commonly PBLHs (user adjustable) in order to produce a wider area of height search and in this way, we exclude larger height‘endings’ in order to exclude the possibility of capturing higher cloud‘endings’ rather than PBLH and we also exclude lower layer‘endings’, mistaken as PBLH (user dependable). It is noted that the common PBLH provided by the user is not mandatory for the rest of MUSTI-PBLH algorithm to produce its result concerning PBLH.

Afterwards, we check which is the most common‘ending’ layer in all times series of 15 to 20’ (maximum an hour) of the 15 or more data files. This is done because it is expected that the present PBLH will not be extremely changed (more than 300 to 500m) in such a small time period. Even in a case of disturbed weather conditions, where a violent change occurs, it will have a scaling that will be easily observed and captured by the tracked records of measurement, which have a time scale of a maximum of around one minute. Then, an average height of the measurements of most common‘endings’ in terms of 10 to 15 height bins difference (user adjustable) will be given as PBLH, even if the latter is dynamically and violently changing. If more than one common ‘ending’ is found with the same higher frequency of occurrence, several more files must be added to that search and more than 15 data files might be needed to add to the processing in order to have a common“ending” with at least double occurrences in comparison with the second most frequently occurring common‘ending’ (recommended). Most of the layers could be clouds and a cloud physically ends, or dramatically change their ‘endings’ in this short time period or even better, when the recording time is close to 30 or even 45’ (user adjustable). So, in this way, the clouds should not produce common‘ending’ heights in that search process and, thus, they are excluded as elements, from the PPBLH (Possible PBLH) matrix.

Of course, it is possible for any layer or cloud, in a dynamic or smooth atmosphere, to have an‘ending’, in some or all time frames that we are searching. It is extremely difficult though to have it steady, at the same range, when tracking space filters begins with a spatial resolution of 1 bin (=7.5m). In any way and in order to avoid that small change by tracking a common layer or cloud‘ending’, as PBLH, except from the adjustment of the recording time, we keep track of the LAYER matrix and the corresponding‘beginnings’ to the remaining ‘endings’, if these exist, by taking into account LAY and NLAY matrices of each recording file.

This is done, in order to decide whether to subtract these‘endings’ from a final PPBL matrix of Possible PBL heights. After every‘beginning’, an‘ending’ is normally expected, even if that denotes the desirable PBLH. To accomplish that automatically, we can also keep track of a number of packages that already produced PBLH before and their‘endings’, 2, 3, 4, ... (user adjustable), in order to extend the observation time and correct faulty PBLH retrievals, keeping in mind that a PBLH could have a small change in height but usually not more than 300 to 500 m in an hour, according to Stull ([An Introduction to Boundary Layer

Meteorology”, Kluwer Academic Publishers, Dordrecht, (2013)]).

So, for a disturbed weather, 2 or maximum 3 previous PBLH retrievals might be of use (recommended), but under milder weather conditions, 3 or even 4 PBLH retrievals and more, could be of valuable use. The reason is that in disturbed weather the PBLH might change dramatically each time and in easier layered weather conditions, more recording time might be needed to spot clearly, the exact PBLH. The maximum recording time that we could take into account is recommended to be an hour or less. The DECIS algorithm that will be presented herein below also takes part in that automatic decision, by answering the questioned‘present weather conditions’, i.e. the factor that is mainly responsible for the PBLH growth through time.

The algorithm can be restricted by a minimum and maximum altitude of commonly PBLH area, of 1000 to even 5000 m, depending on the geographical position (latitude) of the measurement and the most common PBLH values through previous years, at the corresponding time the measurements were taking place. The exact PBLH value will be delivered, after an average of most probable common ‘endings’, with slightly different heights, that have been left and have been chosen for the PPBLH matrix (Fig. 11).

MUSTI-PBLH IN PROCESS

At this point, provided that a clear vision of the layers of MUSTI-L/D algorithm is supplied, further processing of the lidar signals proceeds as follows:

(xviii) After performance of steps (i) to (xvii) of MUSTI L/D algorithm for more than fifteen (or more user adjustable) data files of p_aCr(L,R), which correspond to 15’ or more of recording time, the data obtained are put into an MLAYER matrix, which is a combination of the LAYER matrices, of all recorded data files. At this step we make use of the‘time filters’, for the first time in the MUSTI algorithms. We keep track of the MLAYER matrix which includes LAYER1 to LAYER15 or more (user adjustable) matrices, separated in 15 or more columns by 15 or more (accordingly) time data files: MLAYElh {LA YE It LA YER.,· -,LA YER5) (47)

(xix) In collaboration with the DECIS algorithm described herein below, another type of usage of the‘space - time filters’ of the“MUSTI-PBLH” algorithm is made. This step is used in different ways by means of different selections of ‘space filters’ and their product matrices agreement, according to the blurriness of the atmosphere (DECIS decision) and then, new ALAYER1, ALAYER2 ... ALAYER15 matrices (for 15 data files) are being produced accordingly (cf. steps from (i) to (xix)). For example, in a‘Clear Atmosphere’ (produced by DECIS or by WEATHER PHENOMENA algorithms), different‘space filters’ and their products (spatial derivatives, etc.) have to be in agreement for a‘beginning’ or an ‘ending’ to be calculated, like the‘space filters’ M2, M3, M5 and their products. For a ‘Blurry Atmosphere’ this agreement will not be enough and, thus, M7, M9, Ml 5 and M27 ‘space filters’ and their products, depend on the “blurriness” given from DECIS (or WEATHER PHENOMENA), have to be in agreement for a‘beginning’ or an‘ending’ to be calculated this time. The“blurriness”, given by DECIS, has to do with the wideness and/or compactness, in terms of height, of the existing atmospheric layers and the appropriate agreement of the corresponding‘space filters’ and their products for a‘beginning’ or an ‘ending’ to be calculated.

The bigger and widener the layers are, the longer‘space filters’ need to be taken into account. The thinner the layers are, the lower the total blurriness becomes and the smaller‘space filters’ are needed to be taken into account this time. The“blurriness” of the atmosphere can be countable according to DECIS algorithm (and/or WEATHER PHENOMENA algorithm), but in any case, DECIS cannot decide whether there is a layer, a‘beginning’ or an‘ending’ of a layer and the exact R that they occur. This is done by MUSTI-L/D algorithms with the assistance of DECIS.

(xx) A new matrix called NLAYER is created, which is a matrix where track is kept only of the‘endings’ with range.

(xxi) At this point, the user inserts commonly occurring PBLH area for the site under study, according to time of the day and season. A wider space in terms of height is thereafter created (user adjustable). Any records in the NLAYER matrix outside of these height limits can be excluded or being used as secondary, possible PBLH, if the latter is not provided by the model (user dependable).

(xxii) Optionally, especially for disturbed atmospheres - DECIS algorithm decision is provided. At this step, we check the remaining‘endings’ of the ML A YER matrix and we also check the ALAYER matrices to find their corresponding ‘beginnings’, which can be collected in another BLAYER matrix. A good tactic is always to extract‘endings’, when their corresponding‘beginnings’ are close to an area of 2 to 3 range bin heights (3 * 7.5m = 22.5m). An‘ending’ has a corresponding‘beginning’ at the lower height. We suppose here that there are no clouds or other layers thinner than 20m, across a wide area of the sky. In any case, if an atmospheric layer thinner than 20m exists, then this could probably be a thin cloud or a thin air mass travelling and could be subtracted. A new matrix called XLAYER is being created with the remaining calculated‘endings’.

BLAYER= {ALAYER,· ·, ALA YER 5, MLA YEI§ (48)

(xxiii) Optionally, in case of stable weather conditions prevailing with strong and stable layering in the course of time, we can make use of the previously produced PBLHs (2, 3 or even 4) in order to estimate the‘endings’. A wider XLAYER matrix is being produced, as it is being reinforced with previously produced PBLHs. (xxiv) At this step, we check the XLAYER matrix for common‘endings’ in terms of 10 to 15 range bins height (user adjustable) in all columns of our recording time. Common‘endings’ are being denoted at this step and an average is being calculated to those commonly remaining ‘endings’ in terms of 10 to 15 range bins height. These averaged common ‘endings’ to that entire recording time are kept at a new PPBLH matrix, as possible PBLHs.

(xxv) At the PPBLH matrix we check for the most common‘ending’ to be chosen as a final possible PBLH value. If more than one PBLH value is being produced with the same frequency of occurrence and their difference exceeds 300 m, then more data recording files are needed until the most common one appears at least with double frequency of occurrences than the next most common calculated height ‘ending’ (user adjustable). This type of processing occurs, if the recording time is less than 30 or 40 min. If the recording time already exceeds these time frames, then, less recording time (e.g. 10 min) is suggested to be selected so that a PBLH can be calculated. This is done, because at an overcast weather for example, cloud‘endings’ may appear at the top of the PBL and the PBLH might be changing with time, according to the prevailing weather conditions. Therefore, over the limit of 30 to 40min, if more than one PBLH values are produced with the same frequency, it is suggested to limit this recording time down to 15-20’ (user dependable) as a first step, but always less than 40 min, in order for recorded data files from‘clearer skies’ over the lidar system, to produce the true PBLH value with a greater possibility of appearance. In case of a PBLH difference less than 300m (normal oscillation of PBLH of a nearly disturbed atmosphere) and recording times less than 30min, an averaged PBLH value may be taken (user dependable).

(xxvi) If the most common element of PPBLH appears less times than half the number of the recording files, then similar recordings to our PPBLH matrix (within 10 to 15 bins) are taken into account and an average value of these recorded‘endings’ are limited to very few or even to one most common PPBLH element. The algorithm could also return to step (xxiii) with the results of PPBLH matrix and continue its processing with these new remaining data. This latter procedure can be carried out 2 or 3 more times (user dependable), if no PBLH is produced, according to conditions of the most frequent PBLH (step (xxv)). Otherwise, the algorithm proceeds to the next step.

(xxvii) The final, most common PPBLH recording is being given as PBLH. If the height of other PPBLH recordings presents the same possibility at that height, then the algorithm returns to step (xxvi) and is repeated till the last step, taking into account more or less data files, according to step (xxvi). This procedure could be used for 2 to 3 times or until a PBLH value is produced (user dependable), giving a‘flag’ of that condition to the user.

(xxviii) If at this step, still no PBLH can be calculated, then the closest PPBLH to the most commonly PBLH area, as denoted at the beginning by the user, in a range between 100 to 300m (user adjustable), could be delivered as PBLH value, giving, again a‘flag’ of that condition to the user.

(xxix) If none of the above conditions are met (steps (xxiv) to (xxviii)), then, no PBLH can be calculated and a‘No PBLH found’ alert will be given to the user. This can happen, for example, when disturbed weather conditions exist at various heights and, thus, no PBLH can be calculated within the range of the lidar system.

The way that the ‘products’ are been calculated by MUSTI and the differences of the ‘measurement values’ of the space filters, seems to limit down any fatal miscalculations of these algorithms. We found a successful PBLH calculation rate in more than 85% of the cases we studied within a range height of 100 to l50m (taking into consideration the variability of the PBLHs during the recording time, especially for disturbed atmospheres) (Fig. 6). The evaluation method of the PBLH calculations is based on the relevance between the PBLH‘derived’ by the RCS color plot (cf. Fig. 6) and the one produced by the“MUSTI- PBLH” algorithm, for same time period. We note here that the backscattering at 1064 nm (cf. Fig. 7a) seems to be one of the best ways of aerosol tracing and one of the most successful ways to detect and visualize aerosol presence and, thus, the PBLH, as a function of time. It is clear that:

b_Ioia1{l,K) = b„_o1{l,K) ₊b_aecl,K) (49)

and at 355nm, P_totai^₃₅₅,R) is different than _totai^io64,R) at 1064 nm and according to Mie scattering, aerosols are easier to detect at 1064 than at 355 nm, using Eq. 1. However, the way MUSTI works with derivatives and percentage of change of p_aer^,R), seems to do the exact same job at any wavelength including that of 355nm.

In any way, MUSTI and its input parameters (i.e. a 10% change in the values of p_aer(L,R)) are user adjustable according to the wavelength used and the behavior of p_aer(L.R) and accordingly of a_aCiiT,R) or RCS (Eqs. 2, 6), based on Fig. 1.

Additionally, if the lidar system has an aerosol depolarization capability, then MUSTI can be enhanced with an additional module, providing the ability of characterizing the aerosols by their shape, at any height or distance at the received channels. As MUSTI and DECIS are able to be adjusted to work at any wavelength.

Benefits of the MUSTI-PBLH Algorithm lie in the ability thereof to acquire PBL Height with a simple low cost Lidar device, in the small time that needs to be spent in processing (a few seconds only), the need of minimum time recording, the PBLH retrieval in near real time, the elimination of demand of former PBLHs, contrary to the EKF method wherein such demand exists. Further MUSTI-PBLH Algorithm does not necessitate or depend on the user’s ability to estimate the PBLH through elaboration of data deriving from an imaging of the backscattering signal or to set initial conditions of processing as required by the EKF method, wherein such user’s estimates, especially in combination with a likely change of initial conditions in the course of time (e.g. appearance of clouds) would lead to loosing track of the PBLH and providing faulty results.

THE DECIS Algorithm

As MUSTI may find difficulties, sometimes, to provide the true PBLH between the first few hundred meters to ~l200m above ground, due to the intense backscatter signals from low altitudes and the high decrease in the values of p_aer^,R) (and/or a_aer(^R), RCS) with range, the DECIS algorithm is being introduced to provide assistance and limit the search area as well as the possible faults of PBLH retrieval.

DECIS looks at the quality of the atmosphere through the laser beam of the LIDAR, by dividing the signal with space and looks on the intensity measurement of the backscattering signal, trying to understand where the atmosphere is close to or clearly molecular and give possible PBLHs as outcomes. Those outcomes can always be inputted as most probable PBLHs to“MUSTI-PBLH” to help with“MUSTI-PBLH” single outcome of PBLH.

The DECIS algorithm, in general,‘slices’ the probed atmosphere for a vertically pointing lidar system (favorably) or slant one, at steps of 200 m (user dependable) and measures the t (Eqs. 8, 9) from the low-height lidar measurements (typically up to ~l200m at mid-latitudes, due to usually intense aerosol backscattering) and for typical PBLH values from 2500 to 4000m. The DECIS algorithm can be used at any height or distance (vertical, slant), in order to provide the blurriness of the atmosphere and to assist the MUSTI algorithms to avoid any miscalculations. Then, it calculates t at every atmospheric‘slice’ and if t exceeds a barrier of 2*l0 ³ at 355nm, where with C(7,R) = 10 sr, is the beginning of ‘Moderate Fog’ for a_3bG(l, R) according to “WEATHER PHENOMENA” algorithm, at 355nm, then it ‘counts’ this atmospheric slice as a blurry one (any other arrangement can be made by the user, according to the site location, the range resolution of the measurements, etc.).

“DECIS” counts the number‘slices’ which are declared as blurry and at which altitude (lower than l200m or in the limits of the area height at possible PBLH, as denoted by user), from the whole of the 10-15 data files (time recording of 15’ or more - user depending). Then, it addresses MUSTI-PBLH to the use of MUSTI-PBLH LOWER ALT-l, LOWER ALT-2, LOWER ALT-3, LOWER ALT-4 or MUSTI-PBLH HIGHER ALT. These last algorithms are all the same MUSTI-PBLH algorithm, but with different settings, which are already inputted for processing, accordingly.

These different settings have to do, firstly with the search height area and with the agreement of 3, 4, 5, 6, 7 or even 8 of continuously‘space filters’ and its produced matrices from MUSTI algorithms. For example, M2DR<(-0.2l2) and M3DR <(-0.212) and M5DR<(-0.2l2) must occur in order to be able to continue with“MUSTI-PBLH” and finally, produce possible PBLHs at PPBLH matrix, at a‘Clear Atmosphere’. For this case, the setting chosen will be MUSTI-PBLH LOWER ALT-l .

In general, different combinations of conditions (‘space filter’ products agreement) must be fulfilled in order for the preset “MUSTI-PBLH” algorithms to run accordingly (user adjustable). The idea of selection and the use of the above mentioned matrices and the required conditions for the production of the preset of the algorithms have to do with the automatic selection of the most effective“MUSTI-PBLH” algorithm. The“correct” choice of the agreed conditions (space filter and produced matrices), depends on the level of blurriness and the altitude of occurrence thereof, that weather data files present and lidar data shows. Then, MUSTI-PBLH LOWER ALT-l, LOWER ALT-2, LOWER ALT-3, LOWER ALT-4 or MUSTI-PBLH HIGHER ALT are being used to produce accurately the PBLH.

DECIS triggers the way that PBLH will be calculated out of the PPBLH matrix. So, firstly, it “defines” the blurriness of the atmosphere and the heights where it occurs (with a resolution of 200m) and secondly, it“drives” the way that“MUSTI-PBLH” algorithms will work. For example, in a blurry atmosphere till up to 3000m height, as denoted by DECIS, MUSTI- PBLH will try to find PBLH most probably around the height of 3000m, if recordings like these, are elements of the PPBLH matrix and the MUSTI-PBLH rules of processing, are met (user dependable).

The idea is that, if a blurry atmosphere continuously exists from ground up to a certain height (during a small period of time of recording, typically l5-30min) and above this height level and up to the height limit set as common PBLH areas,“DECIS” is not denoting blurry ‘slices’, .then it becomes most probable that a molecular atmosphere is being reached. So, when molecular conditions are met (p_aer^,R) and/or a_aer(7,R) become equal to zero), the PBLH is found at lower heights. Thus, with the“MUSTI-PBLH” algorithm we will try, in a “smart” way, to look, with DECISs’ assistance, for the presence of PBLH, at these heights, starting above the compact and blurry atmosphere up to the upper height of commonly PBLH area (user adjustable).

Another example that can be mentioned is the one of a layered atmosphere, where“DECIS” finds layers of blurriness below the height at which we have set as the upper limit of possible PBLH. Then, taking into account the way that MUSTI algorithms works, with‘time filters’, we can set the recording time close to 30 min or even longer (but less than 60 min). Then, the remaining records of common‘endings’ of MUSTI are taking into account the fact that PBLH values can vary within l50m (for 15’ recording time), 300m (for 30’ recording time) or about 500m (for 60’ recording time), according to the PBL theory of the dynamic atmosphere by Stull (R. B. Stull,“An Introduction to Boundary Layer Meteorology”, Kluwer Academic Publishers, Dordrecht, (2013)).

On the other hand, in the case of a‘Clear Atmosphere’ or‘Crystal Clear Atmosphere’ where DECIS finds minimum (or no) blurriness of the atmosphere, at all available heights, then MUSTI mainly depends on the agreement of certain‘space filters’ and their products, which lead to the remaining of common‘endings’, through the whole time records of 15’ or more. So, as we can see that the combination of DECIS and MUSTI algorithms enables to study the PBLH variations at any case of atmospheric conditions.

To summarize, the“DECIS” algorithm is being used after the DENOISING algorithm has been applied and it decides whether to use MUSTI-PBLH LOWER ALT-l, LOWER ALT- 2, LOWER ALT-3, LOWER ALT-4 or MUSTI-PBLH HIGHER ALT (or other form of settings upon the user’s decision) for the retrieval of the PBLH (cf. Fig. 11).

DECIS IN PROCESS

In this section we present in an analytic way how“DECIS” algorithm works. In particular:

(i) We start from the processed lidar data providing the values of p_aeitA,R) (and/or a_aerOCR), RCS) with range R. Then, we divide R into‘distance slices’ of 200m (recommended). These ‘distance slices’ might be in a vertical or slant pointing direction. In case of vertical‘distance slices’ we can limit our search, typically, below 1500 to 2600m according to PBLH values over the site under study (user adjustable) and also from 2l0m to l200m (in case of low PBLHs where intense backscattering is observed). Although this presentation deals with vertical measurements, it is estimated that the following steps are also valid for slant measurements using the FASTPLAN technique that is described below.

(ii) Then, we calculate the optical depth (t) (Eqs. 8, 9) of every‘distance slice’, after the application of the “WEATHER PHENOMENA” algorithm and classification of the atmosphere, where lidar is pointing, with range. The intention is, at the next step, to compare this (t) of every‘distance slice’ of 200m with the one for‘Moderate Fog’, according to “WEATHER PHENOMENA” algorithm for 355nm, or, according to new ‘atmospheric classification limits’ at other wavelengths.

(iii) If t>2*10 ³, which, according to“WEATHER PHENOMENA” algorithm and according to [R.M. Measures,“Laser remote sensing. Fundamentals and Applications”, Krieger, Sys No 9247, MEA 621.3678 (1992], is the limit of the beginning of ‘Moderate Fog’ at 355nm (recommended for aeronautical use), we assume this ‘distance slice’ as a blurry one. Otherwise we could assume that this‘distance slice’ could have a value limit of x > 2*l0 ², for extremely blurry atmospheres, which according to “WEATHER PHENOMENA” algorithm and according to [R.M. Measures, “Laser remote sensing. Fundamentals and Applications”, Krieger, Sys No 9247, MEA 621.3678 (1992] is the limit of‘Cumulus Cloud’ at 355 nm, with C(k,R)=l 0sr as well. By these means, the user could give a different meaning to a ‘Blurry Atmosphere’ within the “DECIS” algorithm and, thus, follow a different approach. In our case we keep the limit value of x>2* 10 ³.

(iv) Then, we check which of the‘distance slices’ show a t > 2*10 ³ (cf. step (iii)) and mark the distance points (R_y) where this phenomenon occurs. According to the results of that procedure we can have a clear understanding of the blurriness of the atmosphere where lidar is pointing on, in terms of‘Blurry Slices’ with height. (v) If the lidar is pointing in the vertical direction, then the“DECIS” algorithm immediately addresses the continuity of this process to a certain MUSTI algorithm, with the appropriate settings (by means of the agreement of certain‘space filters’ and their product matrices). For example, for an‘Almost Blurry Atmosphere’ from 600 to 1200 m of height distance, the agreement of the 3^rd, 4^th and 5^th ‘space filters’ and their product matrices (first order derivatives, second order derivatives etc.), have to be in agreement, in order to surpass the obstacle of miscalculations of a faulty PBLH to be produced.

(vi) It is now time for us, to define the blurriness of the atmosphere at the wavelength of 355 nm. For example for the height limit of commonly occurring PBLHs between 210 to 2600 m, and in accordance with the blurriness-cloudiness of the sky, we derive blurriness for each recording file as“Crystal Clear Atmosphere”, if in all or 11/12 of the‘distance slices’ are not denoted as blurry;“Clear Atmosphere” if the 2/12 to 3/12 of the‘distance slices’ are denoted as blurry ones;“Almost Blurry Atmosphere” if this occurs at the 4/12 to 5/12 of the‘distance slices’;“Blurry Atmosphere” if this occurs at the 6/12 to 9/12 of the‘distance slices’ and “Extremely Blurry Atmosphere” if t>2*10 ³ in 10/12 to 12/12 of the‘distance slices’.

The upper limit of the PBLH can be set according to the available common PBLH data of the area under study, as already described, by setting blurriness of the atmosphere and limitation for the“DECIS” algorithm. This kind, of categorization of blurriness is recommended and the user can choose to define its own metric system of blurriness to be used within“DECIS”. In any case, height-distance of‘blurry slices’ and their height formation, is the key for MUST! algorithms to be used, in terms of agreement of conditions to be met among matrices and their products, for PBLH to be calculated later at the exact heights. In Figs. 6 and 7, we show two random cases for“DECIS” and“MUSTI-PBLH” algorithms in a smooth collaboration.

(vii) In the following step an average value of the blurriness of the atmosphere for the whole recording files (15 or more record data files) is being produced for specific PBLH settings to be used in that period of time.

(viii) Each record data file can drive the correct use of the PBLH settings, separately, according to their“DECIS” output results and/or the PBLH settings can be used from the “DECIS” results for the total atmospheric volume studied, in terms of‘averaged’ values of the atmosphere. It is, thus, recommended that“DECIS” should be used, separately, for each record data file, to drive its individual MUSTI PBLH algorithm.

The hereinabove DECIS Algorithm is beneficial in that it provides categorization- classification of atmospheric conditions in 3D with accurate Lidar measurements and in that it can be used for providing a precise layout of the blurriness of the atmosphere in 3D and/or to advantageously drive the use of other algorithms, such as MUSTI-L/D/PBLH).

3D PBLH RETRIEVAL - FASTPLAN TECHNIQUE

Up to now we used only lidar data from a vertically pointing lidar. In this work, in parallel to the MUSTI and DECIS algorithms, the slant range measurements are also proposed, to provide PBLH retrievals. One of the main advantages of a 3D lidar is that it can be deployed at any location and it can perform slant range and thus, conical 3D measurements. In this case a new technique called FASTPLAN, is being applied for slant range measurements to derive the PBLHs.

A 3D scanning lidar starts performing slant measurements at a zenith angle of 30° (or differently, user dependable) in order to acquire, at first, the molecular atmospheric backscattering signals at the vertical or near the vertical direction. We, then, set the rotation direction of the 3D scanning lidar to be equal to the rotation direction of the Earth (from the west to the east direction), especially in the case of low air mass speeds above the lidar’s location at heights where the commonly occurring PBLHs appear.

The reason why we recommended the 3D scanning lidar to have the same rotation direction as that of the Earths’, is because of the inertia of the air masses when they meet the topography of the Earth’s surface. In general, it is preferred to scan at the opposite side of that of the moving air direction, for common PBLH heights according to different seasons of the year.

The idea is that, by scanning a wide area above the location of a 3D scanning lidar (the same or even wider areas than those scanned with vertical pointing lidars), we can obtain the same results at even smaller timeframes (for example, less than l5min), as those produced by steady vertically pointing systems at much longer timeframes. In this way, the 3Dscanning lidar is able to scan the same sky area (or even wider ones) as the one during longer time periods, when pointing vertically, with less lidar data recording files, at much less time periods (half of the time or even less).

The idea of the greater time periods required by the MUSTI and DECIS algorithms (higher than 15min), had to do with the conception of producing as many lidar backscattering data recording files, over a wide sky area, in order to derive common valuable information about the PBLH.

In Fig. 9a a 3D scanning lidar is shown performing the scans for the FASTPLAN technique, by scanning uniformly, at steps of 15° or 7.5° or even 5°, in order to acquire a scan of a wide area of the sky for PBLH retrieval. It is noted that, if meteorological data denotes high air mass speeds over a commonly occurring PBLH height area, at the direction of the earth’s rotation (similar or higher air speeds of that of the linear speed of our 3D scanning lidar), then, an opposite scanning direction should be chosen by the user. The same stands for common great air mass speeds at commonly occurring PBLHs and at locations where the air mass movements have to do with the steady flow of air masses globally.

This is done in order to avoid recordings of the same area of the sky, by means of acquiring the characteristics of the same air masses that move along with the lidar’s direction scanning speed. In that case, the user could calculate a faulty PBLH for that time period for a dynapiic atmosphere or a correct one (e.g. in the case of a steady layered atmosphere). Of course, this situation can be corrected by the use of the MUSTI and DECIS algorithms at the next lidar scan applying the FASTPLAN technique, where a new air mass is scanned.

For slant range measurements of b_3a·(l,II) (and/or a_aer(k,R), RCS) the 3D Stepping technique could be of valuable use. This technique for PBLH retrieval could produce the above values safely, but in order to do so and the time needed to be spent, we could lose our main objective, that is the PBLH retrieval. We have to note that the 3D Stepping technique is very effective when applied to slant range measurements, in order to provide us with p_aer(k,R) (and/or a_aeiA,R), RCS) by going from vertical or near vertical, to slant - horizontal measurements. In this case though, the multi-angle method as described by V. A. Kovalev, W. E. Eichinger,“Elastic Lidar”, Wiley Interscience (2004), pages 295-296 could be used.

In the latter method, slant range measurements of b_3eG(l,B) (and/or a_aer(k,R), RCS) are performed with the hypothesis of horizontal atmospheric homogeneity (presuming homogeneity at“thin horizontal atmospheric slices”) of a wider area, at all heights and the values of b_3eG(l,B) (and/or a_aer(k,R), RCS) are presumed to be the same as the ones in vertical measurements, above the location of the lidar system.

In the case of the PBLH retrieval, we do not use the hypothesis of a homogeneous atmosphere and the operational restrictions it arises, as our goal is not to perform slant and horizontal concise measurements and check for their‘agreement’ in the vertical over our lidar’s location, but to make use of those values (p_aer(k,R)) in the MUSTI and DECIS algorithms with their presented settings, for PBLH retrieval purposes only.

At vertical or near vertical pointing angles (usually 20° to 30° off zenith), the FASTPLAN technique has to do with the ability of the lidar backscattering signal to reach molecular heights. With this technique and in order to use MUSTI and DECIS algorithms with the settings as presented here, we make no use of the locality of the lidar system, as long as the above requirements occur; but we make use of the aerosol layering as retrieved, if vertical measurements were taking place and PBL was positioned hypothetically at higher heights: PBLH=cos ) *R (50)

where, f is the angle with respect to the vertical (near vertical).

“DECIS” and“MUSTI-L/D” algorithms, as they are presented here, can be used also for slant lidar measurements. There are only two steps that the user should consider to change and these are the steps (iv) and (vi) of MUSTI algorithms concerning the differences between the _aer(k,R) (and/or a_aer(k,R), RCS)) values between consecutive‘space filters’. The original settings of these differences were produced statistically during the development of these algorithms, for vertical measurements and should be considered to be set at lower values.

So, for both steps, the cosine rule should be applied, according to angle (f) at each slant - near vertical measurement:

R = h/cos( ) =>h = R*cos< ) (51)

h_C} = G2 * cos (f) = 3 * 7.5 (m) * cos( ) (52)

where, h, is the corresponded height (in vertical measurements) when applying the multi angle method and 1IG3, is an example of vertical appearance, of the original slant one, of the slant‘space filter’ of 3 bins (G3), at slant measurements. Eqs. 38 and 39 should be taken into consideration at step (iv), accordingly.

At step (vi) Eqs. 40 and 41 should be applied in terms of specific set processing, of first order derivative of (M) matrices (of each‘space filter’ at each (M) matrix), the M3DR etc., matrices, and thus:

dh

— = cos (f) -R * sin( ) (53), and

dR

[cos( )-/?*sin( )] dR

where, the first order derivative of b_3bG(l,1i) (and/or a_aer(k,h), RCS) is being calculated versus height and a change in the variable from height (h) to (slant) range (R), depending on angle (f), is being taken into account for M3DR etc. to be produced in slant-near vertical measurements, so that MUSTI or DECIS algorithms do not change their remaining settings or notion of processing.

So, in the case of slant measurements, “DECIS” will drive the appropriate settings for MUSTI-L/D to provide the user with layering and distribution of p_aer(k,R) (and/or a_aer(k,R), RCS) of an atmospheric parcel over the lidar location. When the MUSTI-PBLH algorithms are used, a conversion of these values using the multi-angle method into hypothetical ones, along the vertical direction, without any restriction of the locality of the lidar system. Thus, we use the multi-angle method, in order to produce, hypothetical values of b _aer (l,II) (and/or a_aer(^R), RCS) and not real ones, along the vertical direction. This aims to help the processing in our algorithms (“MUSTI-PBLH”), because these columns carry the same PBLH information, in terms of height, to be extracted. The layering of these vertical columns should differ between them and the layering and the input values of b_3bi·(l,II) (and/or a_aer( ,R), RCS) into the“MUSTI-L/D” algorithms will not be easy to produce common ‘endings’ (Fig. 9b). Even if, some results could be produced in a layered atmosphere, it will not be of any valuable scientific value, because the produced vertical columns are hypothetical and cannot be fully trusted for their relevance between slant range and vertical pointing, in a real atmosphere.

Under operational use at this layered atmosphere case, the hypothetical vertical columns may present the same layering after the use of“MUSTI-L/D” algorithms; but thanks to the fact that the MUSTI algorithms are using vertically pointing measurements, the PBLH can also be produced safely. It has to be noted that the deficiency of operational usage of the multi-angle method is in advantage of the FASTPLAN technique regarding the PBLH retrieval.

In general, at all slant range measurements to derive the PBLH using the“MUSTI-PBLH” algorithm, we must keep in mind that Eqs. 36 to 49, have to be converted from slant range into vertical ones and, then, continue with MUSTI processing and the application of the multi-angle method. As the slant values of b_3bi·(l^) (and/or a_aer(k,R), RCS) will be converted into vertical ones, they will be, still, keeping their original values (b_3a(l^)=b_3a-(l,1i) and/or a_3ei·(l^)=a_3ei·(l,]i), RCS = hCS (height Corrected Signal)) but now, the spatial resolution of 1 bin will become smaller (cfi Eq. 51 and Fig. 9b). On the other hand, one could convert the vertical values of b_3bi·(l,E) (and/or a_aer(k,R), RCS) into slant ones (in any case, if other settings are to be used, then, new‘space filters’, other than the ones presented here, could be originally used).

For example, Eq. 26 could be used for a larger number of bins produced at each‘space filter’ (e.g. for 1 bin = 3.75m, f=45° and M3=3 (*7.5 m) = 3 (*“old bins”), then the“new bins”=6 (*3.75 m)) according to their vertical angle and slant range. In case of other settings, new space filters, other than the ones presented here, could be originally used.

To summarize, in a steady and layered atmosphere, homogeneity could be normally applied in the area of search and then,“MUSTI-L/D/PBLH” algorithms could work normally and produce layering and distribution of the atmosphere in slant and vertical use, as well as PBLH, according to above mentioned process (Eq. 35 in mind for comparing common layer ‘endings’, for PBLH retrieval). In the case of a disturbed atmosphere or generally in an undisputed, non-homogeneous air environment, which is the most common one in a real atmosphere, the hypothetical vertical columns that would be produced with b^l,E) (and/or a_aer(UR), RCS) will not respond to the reality of the atmosphere above lidar’s location. This situation helps the“MUSTI-PBLH” algorithm the most, with its processing, in order to work with the calculated hypothetical vertical value element columns. The values of b_3bG(l,B) (and/or a_aCrU,R), RCS) probably will not“match” at all, when appeared in vertical columns in terms of height, because they have been produced from different vertical angles and the only common layer‘ending’ altitude that would be produced, will be the PBLHs, as it will appear at the same height, at all hypothetical vertical columns, at the area of search (cf. Eq. 48).

In this way, the multi-angle method becomes operational and applicable under any weather conditions, but only for PBLH retrieval purposes. The continuously slant lidar measurements at different angles (f), produce different backscattering signals for processing, but with one common information, that, of same PBLH. The wider the area of the sky scanned, an averaged PBLH produced. It is very difficult in slant measurements, to exceed the limitation of 300 to 500m of common PBLH, for the time frame of 30’ or an hour, accordingly, of vertical measurements. In general, the smaller the area of the sky searched with slant range measurements, the more accurately is the PBLH produced in the local area, but in steady and layered weather conditions, a wider area probably needs to be searched. In disturbed weather conditions, the wide area will not be required (user dependable), because of the seriously unbalanced values of _aer(L,R) (and/or a_aer(L,R), RCS) in terms of common layer height ‘endings’ to be processed, with“MUSTI-PBLH”.

In either case, if the atmosphere is very disturbed or holds heavy clouds at an almost overcast formation, the laser beam of the lidar would probably need all of its power transmittance to reach molecular heights at the receiver and vertical measurements should be preferred. FASTPLAN technique as well as MUSTI and DECIS algorithms, will be imported in a 3Dlidar system to be tested, as they have already developed through real atmospheric backscattering lidar data, but only in vertical pointing, at the LRSU of the National Technical University of Athens (NTUA).

In the following we will briefly present a comparison between the DECIS plus the MUSTI algorithms with the well-established EKF technique as presented by D. Lange et al. [“Atmospheric boundary layer height monitoring using a Kalman filter and backscatter lidar returns”, IEEE Transactions on Geoscience and remote sensing, 52, 4717-4728 (2014)] and [“Atmospheric boundary layer height estimation using a Kalman filter and a frequency- modulated continuous- wave radar”, IEEE transactions on geoscience and remote sensing, 53, 3338-3349 (2015)] with a scope of being able to validate the above algorithms of the invention under different atmospheric conditions.

The Extended Kalman Filter (EKF) as presented in the above references is an algorithm according to which the PBLH transition, modeled by an over-simplified Erf-like curve, is parameterized by four time-adaptive coefficients and noise covariance information estimated by the filter.

In the case of the PBLH retrieval, the EKF is based on the combination of present and past estimates together with a priory model, in order to provide continuous estimations of the PBLH. Alexiou et al. [“Planetary boundary layer height variability over Athens, Greece, based on the synergy of Raman lidar and radiosondes data: Application of the Kalman filter and other techniques (2011-2016)”, Proc. of the 28^th International Laser radar Conference, 25-30 June 2017, Bucharest, Romania. EPJ Web of Conferences 176, 06007 (2018)] showed that the EKF is able to provide realistic values of the PBLH over a specific time, in most of the studied cases, under different meteorological conditions, when compared to local radiosounding data and other techniques, as presented hereinabove: gradient, inflection point, variance and wavelet covariance techniques.

In this work we compare the results obtained by the application of the EKF and“MUSTI- PBLH” for different cases of atmospheric loadings, as presented in the above reference to Alexiou et al., on the PBLH retrievals.

In Figs. 8a, 8b the left column presents the RCS plots produced by the LRSU NTUA lidar system of the present invention, where the curved dashed line (left side) represents the PBLH retrieved when applying the EKF method; the arrows represent the time periods used to derive the PBLH by the VEDRE plus MUSTI algorithms. In the same figure (right side) a zoomed part of the RCS plots, including the same time periods defined by the arrows, showing by the white dashed lines the respective PBLH retrievals based on the VEDRE plus MUSTI algorithms. As a first example (cf. Fig. 8a), for the“etesian case”, the outcome of the algorithms of the invention gave an average PBLH of the order of l350m (07:56:50-08:22:00 UTC), while the PBLH retrieved by the EKF is totally different and close to 2100m. For the“sea breeze” case (cf. Fig. 8b), our algorithms gave an average PBLH of the order of l478m (09:04:30 to 09:28:00 UTC), while the PBLH retrieved by the EKF was l450m, very similar to ours. For the“clear sky” case, the PBLH retrieved by MUSTI plus DECIS was H25m (06:19:20- 06:42:50 UTC) very close to the one retrieved by EKF (1150m). In the“dust” case the algorithms of the invention gave an average PBLH of the order of 743m (05:52:50-06:16: 10 UTC), while the PBLH retrieved by the EKF was of the order of 720m. In the“clouds” case, our algorithms gave an average PBLH of the order of l 873m (07:44:30-08:07:50 UTC), while the PBLH retrieved by the EKF was of the order of 1550m.

It is obvious that the MUSTI and DECIS algorithms are able to retrieve the correct PBLH values, in most cases, as we take also into account realistic values of the PBLH variations, following Stull’s PBLH observations referred to hereinabove.

We can further improve the accuracy of our algorithms regarding the PBLH retrievals, by adjusting the initial settings according to the meteorological conditions prevailing over the studied site. The algorithms of the invention do not need a meta-PBLH estimation, as in the case of the EKF technique, where the user provides the initial’guessed’ value of the PBLH. Moreover, DECIS plus MUSTI are able to provide PBLH estimation, automatically, without user’s intervention, in near real time.

The algorithms presented were tested for a period of 20 to 30’ in every case, as this time window seems to be the most suitable one to produce the common PPBLH matrix recordings and avoid cloud and other layer‘endings’, mistaken, as PBLH. The comparison presented here, showed that DECIS and MUSTI algorithms are on the right track to produce the aerosol layering and its spatial distribution, the _aer(k,R) (and/or a_aer(k,R), RCS), as well as the PBLH variation. In the case where the PBLH retrieved by the EKF was faulty (cf. the“Etesian” and “Clear sky” cases), MUSTI and DECIS algorithms were able to keep the real PBLH values, thanks to their novel working principle without the necessity of a previous PBLHs knowledge and an a priori model (cf. Fig. 8b,“Sea Breese” case, after 11 :30 UTC (Up) and“Cloud case” (Down) in the timeframe studied). The FASTPLAN Technique makes use of the aerosol layering as retrieved, if vertical measurements are taken and the PBL is positioned, hypothetically, at higher heights. Further the FASTPLAN Technique provides the capacity of acquiring PBLH with less Lidar measurements in substantially less time of measuring.

THE WEATHER PHENOMENA-PBLH Algorithm

The “WEATHER PHENOMENA-PBLH” algorithm basically uses “WEATHER PHENOMENA” algorithm in a totally new form for a totally new usage. In this work there is a transformation of the outcomes of “WEATHER PHENOMENA” algorithm through a combination with Stull’s directions as disclosed in [“An Introduction to Boundary Layer Meteorology”, Kluwer Academic Publishers, Dordrecht, (2013), p. 3] for the meteorological existence of PBL, where Stull describes the existence of fog, as a stratocumulus cloud that touches the ground, as an existent boundary-layer phenomenon.

So, by using“WEATHER PHENOMENA” algorithm outcomes we can always set type of fog heights as possible heights of PBLH. Specifically, the“WEATHER PHENOMENA” is able to provide meteorological conditions like clouds, haze, type of fog, etc., versus distance (R), for any vertical or slant range lidar measurements. The same 3D scanning lidar and the above-mentioned techniques and data, can be used to characterize the atmospheric conditions (type of fog, haze, etc.), plus the visibility. The“WEATHER PHENOMENA-PBLH” algorithm follows the steps (i)-(ix) as described in the“WEATHER PHENOMENA” algorithm of the present invention.

Following the above values of a_acr(k,R) obtained in accordance to the above, we conclude the corresponding requirements for an airport visibility up to 10 km. A color discrimination of the‘present’ weather conditions is been acquired and a visualization is being made. The values of C(1,R) that has been used through this process were according to Figs. 2a-2c as explained hereinabove.

(x) If b _aer (}„R) < 2* 1 O ⁷ m ¹ sr ¹ (a_aer(k,R) < 1 *10 ⁵ rrf¹), no color discrimination is being made and the visibility goes higher than 10 km.

(xi) A check of the AOD total value in comparison with a_aer(A,R) is calculated and visibility versus R can also be constructed.

After following the above steps and the classification of the atmospheric conditions made by the“WEATHER PHENOMENA” algorithm, the usage for the PBLH finding is being introduced. So, as mentioned before, in conformity with Stull’s directions for fog existence as a boundary-layer phenomenon, extra steps are taken where the“WEATHER PHENOMENA- PBLH” works from the maximum valuable range-height where useful signal (Signal to Noise Ratio - SNR over 3 to 5, adjustable) of a 3D lidar pointing to higher molecular heights, towards lower range-heights where types of fog (thin fog, light fog etc) exist and gives these heights as possible PBLHs (Fig. 10). Specifically:

(xii) Starting from maximum valuable range, where useful signal (Signal to Noise Ratio - SNR over 3 to 5, adjustable) is observed, to lower ranges - heights.

(xiii) If (for 355nm lidar wavelength) p_aerU,R) > 2* l0 ⁷ (with C(k,R) = 50sr and a_aer^,R) > 1 *10 ⁵) which correspond to the lower limit of‘Sky Clear’ and p_aer(k,R) > 4*l0 ⁵ (with mean C(k,R) = 6.875 sr (6.25 < C(k,R) < 7.5 sr) and a_aer(k,R) > 3* l0 ⁴) which correspond to the upper limit of‘Thin Fog’ identification and/or _ae,(k,R) > 8*1 O ⁵ (with mean C(k,R) = 8.125 sr (6.25 < C(k,R) < 10 sr) and a_aer(k,R) > 5* 1 O ⁴) which corresponds to the upper limit of ‘Light Fog’ identification, and/or (for colder environments) _aer(k,R) > 2* 1 O ⁴ with C(L,R) = 10 sr and a_aer(k,R) > 2* 10 ³ which corresponds to the upper limit of‘Moderate Fog’, and/or make use of‘Haze’ spectrum (check steps (vi) and (vii) above on this algorithm) for hotter climates, but basically trying to find a type of ‘Fog’ according to “WEATHER PHENOMENA” algorithm, then the outcome is possible PBLH.

(xiv) The highest of the above possible PBLHs is being given as PBLH with preferably limits first for‘Light fog’, then‘Thin Fog’ then‘Haze’ for hooter and dryer climates or‘Moderate fog’ for more humid climates .

If used in combination with“MUSTI-PBLH” and“DECIS” algorithms then follow one of the below mentioned paths (xv, xvi, xvii) and/or (xviii):

(xv) Perform the above mentioned steps (i to xiii) for as many signal files as selected for “MUSTI-PBLH” and create a medium value from each of the possible PBLHs found from step (xiii) of this algorithm, in terms of being inside the spectrum of 10 to 15 range bins (100 to l50m) of most of the outcome heights of this specific step (xv).

(xvi) Proceed to the step after step (xxiv) of“MUSTI-PBLH” and give the outcomes of step (xv) of this algorithm to be taken as most probable PBLHs in terms of 10 to 15 range bins (100 to l50m) and for the rest possible PBLHs of “MUSTI-PBLH”, at this point, to be excluded. (xvii) Proceed to step (xxv) of“MUSTI-PBLH” and continue with the procedure of that algorithm.

(xviii) Go to the input part of data for b_3bi-(l,II) and/or a_aeif ,R) from“MUSTI - L/D” algorithm (which is part of“MUSTI - PBLH”) at step (i) and if p_aer^,R) > 2* 10-7 (with C(L,R) = 50sr and a_aer^,R) > 1*10 ⁵) which correspond to the lower limit of‘Sky Clear’ and p_aer(W > 4*l0 ⁵ (with mean C(/.,R) = 6.875 sr (6.25 < C(/ R) < 7.5 sr) and a_aer(L,R) > 3* 10- ⁴) which correspond to the upper limit of‘Thin Fog’ identification and/or p_aer(7,R) > 8*10 ⁵ (with mean C( ,R) = 8.125 sr (6.25 < C(7,R) < 10 sr) and a_aer(^R) > 5* 1 O ⁴) which corresponds to the upper limit of ‘Light Fog’ identification, and/or (for more cold environments) p_aerU,R) > 2*l0 ⁴ with C(/^R) = 10 sr and a_aer(7,R) > 2*1 O ³ which corresponds to the upper limit of‘Moderate Fog’, and/or make use of‘Haze’ spectrum for hooter climates, but basically input data for type of fog according to “WEATHER PHENOMENA” algorithm and then the final outcome is most probable PBLH (which is step (xiii) of this (“WEATHER PHENOMENA - PBLH”) algorithm).

MUSTI, DECIS and WEATHER PHENOMENA-PBLH algorithms make a most useful combination for the PBLH retrieval and can be used from aviation, meteorologist as well as from the scientific community and commerce, with remarkable results, over 93% on cases studied, for a 3D or vertical pointing lidar or signal, strongly reacting with the atmosphere.

Of course, the limits where“WEATHER PHENOMENA-PBLH” determines can also be altered for a total Sky clear. These limits will be given once“DECIS” characterizes the above atmosphere as totally clear atmosphere.“WEATHER PHENOMENA-PBLH” can also be used in combination with FASTPLAN technique in order to acquire PBLH in 3D and can work independently or in combination with the above mentioned algorithms in vertical and in 3D pointing. Of course in 3D pointing and“DECIS” algorithm’s help (user dependable) “WEATHER PHENOMENA-PBLH” works in most effective way with remarkable results (Fig. 9c).

The WEATHER .PHENOMENA - PBLH Algorithm provides the beneficial advantage of acquiring PBLH with just a single Lidar measurement in a single beam even with a single laser shoot. It is herein noted that the combination of WEATHER PHENOMENA-PBLH and MUSTI-L/D/PBLH Algorithms, such combinatory application being depicted in the Flowchart of Fig. 11, provides a maximally accurate retrieval of PBLH.

THE SIBESMEA Algorithm

Another novel algorithm called“SIBESMEA” (Single BEam Speed MEAsurement) has been created for the use of a single beam, of a 3D or vertical pointing lidar, to be able to be used as a wind or observed object speed measurement equipment through time. In particular a single beam is being pointing to a certain direction creating recorded signal data files to be used by the above mentioned algorithms, and maybe in other algorithms chosen by the user.

Then taking over the recordings backwards in time it is possible to track the speed of certain atmospheric phenomena - conditions described in WEATHER PHENOMENA and MUSTI algorithms and/or of an object, that is detected by the lidar or signal beam and if it is still tracked in that beam. In other words, identification of the atmospheric layer and its configuration in signal analysis is being made and tracking, though time but backwards, is being made, for the same or almost the same atmospheric and signal condition, for the extraction of speed of that atmospheric condition or object travelling and tracked through time. Specifically:

(i) The 3D Stepping technique is being used in every measurement for 3D lidar measurements. (ii) Recording of 15 or more (adjustable by the user) lidar or signal data files is being kept.

(iii) DENOISE 1/2 or other noise subtraction algorithm is been applied on the received signals.

(iv) Eqs. (1) to (7) are been used to provide the appropriate coefficients if Klett and the 3D Stepping techniques are being used for a 3D lidar and in general any kind of technique used to provide safely (less errors) valuable signals (any coefficient versus Range in 2 or 3D) for processing after denoising.

(v) Continuously recorded signals are being placed in one plane for 2D signals (or for each of the 2D signals being formed out of the 3D or multi-D signals).

(vi) The acreage of every 3 or 5 or more range bins (user adjustable) is being calculated for each of the continuously recorded signals, in order to have square bins or other unit of measurement (m² etc) calculations. A ratio (similarity) of continuously acreages or other shape similarity technique used, between two of more continuously recorded data files- signals, greater than 0.9 (user adjustable), from one signal - recorded data file to the next, could denote tracking of the same atmospheric condition-object, inside the spatial resolution of the ±3 or ±5 or more range bins the user decides to use between the continuously recorded data files - signals. The point in this step is to create variable range filters like in MUSTI space filters creation (‘frames’), in order to find matching - similar, small parts of the total signal figures, for each one of the signal figures being recorded, in order to identify them travelling through time from one signal - recorded data file to the next, with a speed given by the range where the identification has been made and the speed - time of the continuously signals recorded for processing.

(vii) By having time and range in a plane for 2D processing, one is able to have the speed of the object or atmospheric condition as a medium speed (depending on the speed of data recorded and accuracy of the whole equipment and error made by the previous processing techniques until then) or, for further recording, of the acceleration or degradation of speed of the atmospheric condition or object still captured by the same beam, without the need of expensive Doppler equipment. In other words a capture of‘frames’ of the situation recorded by the data files is being made and a rearward analysis in time is being made, in order to find the speed of the similar atmospheric condition or same object in close range travelling.

(viii) The number of‘frames’ - data signals or recording files, has to do with the history we want to observe and the expected minimum and maximum speed of the object we want to observe in our full range. In general, the highest the expected speed of the atmospheric conditions or object, directs to more recorded data files and/or higher repetition frequency of the lidar or signal beam, for pulsed signal. The lower the expected speed of the atmospheric conditions directs to less recorded data files and/or lower or medium repetition frequency of the lidar or signal beam, for pulsed signal. In the latter case and in order to measure very low speeds, higher repetition frequency of the lidar or other signal beam might be used in order to capture very low speeds.

(ix) Wherever the lidar or other signal beam is pointing in a slant direction, angle (cp), h and R from Eq. (51) are being used for one plane presentation and 2D.

(x) If the acreage of step (v) of this algorithm seems to vain through time from a signal - recorded data file to the next, by keeping the above mentioned similarity (ratio greater than 0.9 - user adjustable) and for high repetition frequency of the lidar or other signal beam for pulsed signal, could be safely assumed that the rest of the matter - air mass is being travelling to the 3^rd dimension and an estimation of the air speed in the 3^rd dimension could be given. The air mass missing between continuously measurements, applies for 3^rd dimension wind movement, by presuming high repetition frequency of beam in pulsed signals. The speed of the 3^rd dimension movement of the air mass could be calculated through the rate of air mass loss in that (+ or -) direction which direction is based on step (x) below.

(xi) For the 3^rd dimension estimation, the multi-angle technique could be of use, similar to the one discussed hereinabove. A measurement at a different azimuth angle is being made with the use of steps (i) to (viii) and a history through data recorded files is being created for a period of a few minutes.

(xii) A small angle change is made (recommended of 1° or less) and the same steps are followed, but now the similarity occurs also between the azimuth change angle in order for a beginning of + or - (3^rd dimension) of the wind to be created and recorded - tracked thereafter. Those last two steps (x) and (xi), are followed only once, at the beginning of the 3^rd dimension wind speed measurement or whenever the user needs to recalibrate his system.

(xiii) A 3D wind direction and speed is being created and recorded each time.

This novel algorithm could replace expensive Doppler equipment and make 3D wind speed retrieval possible from a single beam of a 3D lidar or other 3D signal pointing equipment, in a (f) angle from the ground (see Figs. l2a-c).

THE WEATHER PHENOMENA-WIND ALGORITHM

This algorithm (“WEATHER PHENOMENA-WIND”) has to do with the tracking of specific atmospheric layers (e.g. clouds, type of fog etc) and their speed measurement in 2D or 3D formation. Basically the“WEATHER PHENOMENA” algorithm and steps (i) to (xi) thereof is being used in order to create a visualization of the layering of the atmosphere at the 3D pointing direction each time.

Then the use of a second data recording file is being made followed by the same sequence as - above and after that a third one or even fourth one etc (user dependable) in order for tracking of certain atmospheric layers (e.g. clouds, type of fog etc) to be made possible. Afterwards at step:

(xii) Tracking of certain layers is being identified from sequential data recording files inside a range spectrum of the centerline of the specific layer (independently of the width of the specific layer),

(xiii) The speed of the specific layer is being calculated depending on the range the centerline of the layer travelled through time and the repetition frequency the equipment (lidar) is pulsing (or the capture of each“frame” repetition frequency is being made in a Continuous Waved lidar).

(xiv) Once the layer is not yet tracked through sequential data recording files, the procedure stops the tracking and speed calculation of that layer (e.g. clouds, type of fog etc) and the tracking of the rest of the layers continues as long as new layers are produced from new sequential data recording files and new layer tracking followed by their speed is being conducted.

(xv) Sequential layer tracking speed is being given to the user in 2D or 3D when “SIBESMEA” algorithms is followed from its step (vi) and after, until the end of

“SIBESMEA” algorithm for tracking same atmospheric layers in different, but close (sequential), polar angles f and Q for vertical and azimuth setting - pointing Of the equipment.

(xvi) 2D and 3D layer tracking speed measurement is being made available to the user. With this algorithm the tracking of certain layers with certain characteristics, can be made available with the ability of having the speed of every layer tracked in the pointing direction of a 3D or vertical lidar. So, a good use of this algorithm can be made in the airports for safety and urban areas or burned forests - other airborne agents, accidents, in order the locate, identify and track specific important agent carrier layers and their speed in their headings. This novel algorithm could replace expensive Doppler equipment and make 3D layer tracking speed retrieval possible from a single beam of a 3D lidar or other 3D signal pointing equipment, in a (f) angle from the ground (Figs. l3a-l3b).

In summarizing the above, novel algorithms and techniques have been presented applied to 3D scanning lidar slant range measurements to evaluate horizontal, slant and vertical visibility for tower aircraft controllers, meteorologists, but also from pilot’s point of view, as well as for the detection of atmospheric layering of the‘present’ weather and the distribution of p_aer ,R) and/or a_aer(k,R) in three dimensions.

Some typical examples of applications of these techniques with“real” lidar data acquired by the LRSU-NTUA lidar system have been provided. More than 50 case studies were evaluated (it is estimated that a number of cases studies around 50 is ample for providing reliable results) with the above mentioned algorithms and found to be in very good agreement (statistically more than 90%).

In the case studies presented herein, the fully satisfactory performance of the algorithms has been demonstrated as applied on“real” lidar data obtained at 355 nm, versus“real” lidar data 1064 nm, both obtained by the LRSU-NTUA lidar system.

A procedure of classification of atmospheric layers was disclosed with scaled values of p_aer(k,R) and/or a_aer^,R), using Koschmieder’s law and Wright’s diagrams produced by MAPP (SRI), in combination with the herein disclosed 3D stepping technique.

Thus, a good parameterization of the vertical / slant / horizontal range visibility is being accomplished, according to ICAO and WMO rules during daytime conditions. Additionally, conditions for eye safety must be met, according to Measures and ANSI (American National Standards Institute), at 355 nm, therefore the Maximum Permissible Exposure (MPE) should be 10 ³ J/cm ². So, the device used in testing, i.e. a pulsed 3D scanning lidar fulfilling the above mentioned requirements or the EU standard on laser safety EN 60825-1 :2007, proved to be a good choice.

Therefore, the application of the proposed 3D Stepping technique in 3 dimensional pulsed multi-wavelength lidar systems operated at airports, in conjunction with the“WEATHER PHENOMENA” and“VISIBILITY” algorithms of the invention have much to contribute to the estimation of the horizontal, slant and vertical visibility for tower aircraft controllers, meteorologists, but also from the pilot’s point of view.

Further, in order to view layering and PBLH of the atmosphere, we took the same denoised RCS and passed it through DECIS and MUSTI algorithms. These algorithms are able to successfully produce atmospheric layering and PBLH at the vicinity of the lidar pointing at vertical or slant range. They are also able to be incorporated into a 3D lidar and a full processing of the produced RCS can be made, with the above mentioned substantial results. Meteorological reports (METARs - METeorological Aerodrome Reports) can find these algorithms extremely helpful and together with 3D Stepping and FASTPLAN techniques, “WEATHER PHENOMENA” and“VISIBILITY” algorithms, a very attractive and useful combination can be made, especially to airport tower controllers, meteorologists and aircraft pilots as well as for scientific purposes. The most important is that the above algorithms were tested and evaluated for its relevance, mainly with LRSU NTUA elastic-Raman lidar system in many cases, at the same or different days, in different meteorological conditions, different times of the year and through a number of years. They were found to be 100% successful in providing layering and meteorological conditions change retrieval (from type of fog to haze and then to cloud etc) and more than

85% successful on the retrieval of PBLH. This outcome is also in relevance with the results produced by EKF for PBLH retrieval as disclosed hereinabove, plus the additional advantages that MUSTI and DECIS algorithms provide.

These advantages are: ( 1 ) the small time that needs to be spent in processing (very few seconds), (2). the need of minimum time recording, (3) There is no, necessarily, previous knowledge or demand of former PBLHs, like with the EKF method, (4) The near real time produced results, that is totally different from EKF, (5) That there is no necessity or dependence of the user’s ability, to point at initial conditions of processing, to follow PBLH retrieval and then, if conditions change (like the appearance of clouds and parody, change the PBLH, for a matter of time) to lose track of the PBLH, producing for a long period of time PBLH in areas of great dispute. The MUSTI and DECIS algorithms are fully operational, under any circumstances, at any geographical location and at any time of the year, standalone, without the need of any partner to participate for the outcomes to be produced in near real time. A flow chart of all algorithms and their collaboration is presented in Fig. 14. Also, with the assistance of FASTPLAN technique, slant range PBLH retrieval measurements can take place, widening the notion of usage of the latter algorithms also in slant range measurements application.

Claims

1. Method being performed by a lidar device located on ground, maritime or space enviromnent and a lidar data processing unit operatively associated with said lidar device to provide real-time monitoring of meteorological parameters through detection of atmospheric layering including detection of the Planetary Boundary Layer Height (PBLH), said lidar device being configured to scan the aerosol layers of the atmosphere by emitting a plurality of pulses at a wavelength (l) and receiving a plurality of lidar return signals, each of said signals providing return signal parameters (p) comprising a Range Squared Corrected lidar Signal (RCS) that is the received power P’( R) after atmospheric and electronic noise background (BG) correction, a range (R) dependent variable extinction coefficient a_aer(L,R) and a range (R) dependent variable backscattering coefficient p_aer(L,R), characterized in that:

said lidar device being selectively operational in a vertical, slant or horizontal direction, in a two-dimensional or 3-dimensional operating mode;

said lidar device being calibrated through setting said extinction coefficient a_aer(^R) and said backscattering coefficient p_aer(L,R) to a zero value for non-aerosol presence and detection of molecular layer whilst emitting in a vertical direction and is continuously re-calibrated during sequential movement at predetermined slant range steps through setting said extinction coefficient a_aer(L,R) and said backscattering coefficient p_aer(L,R) to new calibrating values that retain as reference the values of a_aer(L,R) and/or or _aer(L,R) of the immediately precedent step;

said lidar device operating with a lidar ratio C(L,R)= a_aer(L,R)/p_aer(L,R), said lidar ratio (C) being set at a predetermined value that corresponds to the type of atmospheric layer being detected, said predetermined value being deducted through the following steps:

obtaining average values of a_aer(L,R) and corresponding range limits from an established diagram of a_aer(L,R) versus visibility depicting measurements employing a wavelength of 550 nm;

calculating a ratio for each pair of adjacent atmospheric conditions, including a ratio for sky crystal clear / sky clear, a ratio of sky clear / light haze, a ratio of light haze / haze, a ratio for haze / thin fog, a ratio for thin fog / light fog, a ratio of light fog / moderate fog and cloud;

applying said ratios for pairs of adjacent atmospheric conditions, said average values of a_aer(A,R) and corresponding range limits in an established diagram of a_aer( R) and p_aer( R) versus wavelength to obtain values of a_aer(L,R) and a_er(L,R) at the wavelength (l) of operation of said lidar device;

said lidar data processing unit operatively associated with said lidar device being configured to process said lidar return signals to generate a quantitative and a qualitative report of atmospheric layering, detection of the (PBLH), wind measurement and estimation of visibility.

2. Method being performed by the lidar device and operatively associated lidar data processing unit thereof according to claim 1, wherein a “WEATHER PHENOMENA” algorithm is adapted to run with lidar wavelength set at 355 nm for generating a qualitative report of atmospheric layering, being adapted to check the retrieved values, through the respective C(L,R), of p_aer^,R) and a_aer(L,R) providing a qualitative report of atmospheric layers as follows:

Values of p_aer(k,R) > 2* 1 O ³ m ¹ sr¹ with C(L,R) = 10 sr and a_aer(L,R) > 2*1 O ² nr¹ indicate Cl =“Cloud and no visibility”;

5 m ¹, indicate SC=“Sky Clear, visibility 6 to 10 km“.

3. Method being performed by the lidar device and operatively associated lidar data processing unit thereof according to claim 1, wherein an algorithm MUSTI L/D is adapted to run for generating a quantitative report of atmospheric layering, said algorithm MUSTI L/D comprising the steps of:

Setting said lidar ratio (C) at a predetermined value of and setting the lidar device into vertical or near vertical direction;

Selecting said variable backscattering coefficient p_aer(L,R) or a_aer(L,R) as the return signal parameter (p) to be monitored and creating a matrix (B) with elements of values of said backscattering coefficient p_aer^,R) or a_aer(L,R) and a matrix (C) with elements of values of the range (R) wherein the values of said backscattering coefficient p_aer^,R) or a_aer(L,R) have been recorded;

Creating a plurality (n) of space filters Gl, G2, G3, ....Gn, each space filter representing a predetermined integer multiple of bins, each bin being identical to the resolution of said lidar device;

Creating a plurality (n) of matrices Ml, M2, M3, ...Mn, each with records of the median values of the elements of matrix (B) and corresponding values of the range (R) of said space filters Gl, G2, G3, ....Gn respectively;

Checking sequential records in each one of said matrices Ml, M2, M3, ...Mn to locate sequential records with a value of

having a difference that exceeds 10% of the immediately precedent value of _aer(L,R) or a_aer(L,R), such difference being associated of a possible change in meteorological conditions that might indicate a“beginning” or an“ending” of an atmospheric layer; Creating a plurality (n) of matrices Hl, H2, H3, ...Hn, each with records of the median values of b_3eG(l,11) or a_aer( ,R) that exceed the immediately precedent value of b_3ei-(l,II) or a_aer(7,R) by more than 10% of said space filters Gl, G2, G3, ....Gn respectively;

Excluding from said matrices HI, H2, H3, Hn those records with values of b_36G(l^) or a_aer( ,R) greater than 8*10-7 m ¹ sr¹ and C( ,R) = 50 sr (a_aer(X,R) > 4*10 ⁵ m->);

Creating matrices Ml DR, M2DR, M3 DR, ..., MnDR of first order derivatives for every space filter from said matrices Ml, M2, M3, ..., Mn;

Creating matrices M1DRP, M2DRP, M3DRP, ..., MnDRP of positive first order derivatives and matrices M1DRN, M2DRN, M3DRN, ..., MnDRN of negative first order derivatives contained in said matrices Ml DR, M2DR, M3DR, ..., MnDR;

Creating matrices DP1, DP2, DP3, ..., DPn for positive first order derivatives of respective matrices M1DRP, M2DRP, M3DRP, ..., MnDRP which are in agreement with respective said matrices HI, H2, H3, ..., Hn, by means of values other than zero, said matrices DPI, DP2, DP3, ..., DPn indicating“beginnings” of an atmospheric layer and matrices DN1, DN2, DN3, ..., DNn for negative first order derivatives of respective matrices M1DRN, M2DRN, M3DRN, ..., MnDRN which are in agreement with respective said matrices HI, H2, H3, ..., Hn, by means of values other than zero, said matrices DN1, DN2, DN3, ..., DNn indicating“endings” of an atmospheric layer;

Creating matrices MIDDR, M2DDR, M3 DDR, ..., MnDDR of second order derivatives for every space filter from said matrices Ml, M2, M3, ..., Mn;

Creating matrices M1DDRP, M2DDRP, M3DDRP, ..., MnDDRP of positive second order derivatives, contained in said matrices MIDDR, M2DDR, M3DDR, ..., MnDDR, which are in agreement with respective said matrices Hl, H2, H3, ..., Hn, and matrices M1DDRN, M2DDRN, M3DDRN, ..., MnDDRN of negative second order derivatives, contained in said matrices Ml DDR, M2DDR, M3 DDR, ..., MnDDR, which are in agreement with respective said matrices HI, H2, H3, ..., Hn;

Creating matrices MIS, M2S, M3S,..., MnS of “beginnings” that include corresponding data of said matrices DP1, DP2, DP3, ·..., DPn and of said matrices M1DDRP, M2DDRP, M3DDRP, ..., MnDDRP in agreement and matrices M1F, M2F, M3F,..., MnF of “endings” that include corresponding data of said matrices DN1, DN2, DN3, .. DNn and of said matrices M1DDRN, M2DDRN, M3DDRN, . MnDDRN in agreement;

Excluding each matrix M_mS of said matrices MIS, M2S, M3S,..., MnS with characteristics that, at the same corresponding range, are not shared by an immediately precedent matrix M_ra-iS and by an immediately following matrix M_m+iS;

Creating atmospheric layering matrices LAYl, LAY2, LAY3 ... LAYn that include the agreement, in terms of range, of said matrices DPI, DP2, DP3, ..., DPn and of said matrices M1DDRN, M2DDRN, M3DDRN, ..., MnDDRN and atmospheric layering matrices NLAY1, NLAY2, NLAY3, ..., NLAYn that include the agreement, in terms of range, of said matrices DN1, DN2, DN3, ..., DNn and of said matrices M1DDRP, M2DDRP, M3DDRP, ..., MnDDRP; Averaging non-zero values of the elements in said atmospheric layering matrices LAY1, LAY2, LAY3 ... LAYn and in said atmospheric layering matrices NLAY1, NLAY2, NLAY3, ..., NLAYn and creating an eventual LAYER matrix that includes elements of LAY1, LAY2, LAY3 ... LAYn and elements of NLAY1, NLAY2, NLAY3, ..., NLAYn respectively identifying atmospheric layers.

4. Method being performed by the lidar device and operatively associated lidar , data processing unit thereof according to claims 2 and 3, wherein data from said

“WEATHER PHENOMENA” algorithm are used by said“MUSTI-L/D” algorithm, at first upon initiation of return signal processing to facilitate processing and upon termination thereof to provide verification of results.

5. Method being performed by the lidar device and operatively associated lidar data processing unit thereof according to claim 3, wherein an algorithm named “MUSTI-PBLH” is adapted to run for the detection of PBLH and comprises the steps of:

Running “MUSTI-L/D” algorithm (n) times and creating an MLAYER matrix containing records being sequentially obtained during a predetermined period of time of said matrices LA YER1, LA YER2, ..., LA YERn;

introducing a plurality of classified settings MUSTI-PBLH LOWER ALT-1, LOWER ALT-2, ..., LOWER ALT-n, each of said settings defining a low threshold value of the altitude at which PBLH is being searched, said low threshold value being dependent on the type of the atmosphere in a specific location and season, each one of said classified settings MUSTI-PBLH LOWER ALT-1, LOWER ALT-2, ..., LOWER ALT-n comprising a plurality of operational modes, each operational mode with predetermined varying values of said plurality (n) of space filters Gl, G2, G3, ....Gn and,

Creating n matrices ALAYER1, ALAYER2, ..., ALA YERn with elements of the most probable PBLH values in the corresponding sequentially obtained records of said matrices LAYER1, LAYER2, ..., LA YERn, wherein said MLAYER matrix is being transformed through substitution of the elements of said matrices LAYER1, LAYER2, ..., LA YERn by the elements of said matrices ALAYERl, ALAYER2, ...,

ALAYERn;

Creating n matrices NLAYER1, NLAYER2, ..., NLA YERn containing the endings of the values of heights contained in said matrices ALAYERl, ALAYER2, ..., ALAYERn, wherein said MLAYER matrix is being transformed again through substitution of the elements of said matrices ALAYERl, ALAYER2, ..., ALAYERn by the elements of said matrices NLAYER1, NLAYER2, ..., NLA YERn;

Excluding PBLH records in each matrix of said matrices NLAYER1, NLAYER2, ..., NLA YERn, that fail to agree with possible PBLHs in the location being monitored;

Comparing atmospheric layer beginnings and endings contained in said ALAYERl, ALAYER2, ..., ALAYERn matrices with endings contained in said MLAYER matrix to discard those endings, that designate atmospheric layer, non-probable PBLH endings and especially those in the range of 3 bins and creating a BLAYER matrix to contain the remaining atmospheric layer endings; Comparing PBLH records in said BLAYER matrix with previously detected PBLH values to discard records acquired above a threshold PBLH rate of change of the order of 500 m/hr and creating an XL AYER matrix to contain the remaining atmospheric layer endings;

Searching within (n) records of said XLAYER matrix to locate common endings within a range of the order of 10-15 bins and calculating an average of detected common endings;

Creating a PPBLH matrix with averaged detected common endings indicating possible values of PBLH;

Searching within said PPBLH matrix to locate most common PBLH values and check the frequency of occurrence thereof;

If a single most common PBLH appears with a substantially increased frequency of occurrences within said most common PBLH values, provide said single most common PBLH as a final PBLH value, or else calculate an average of said most common PBLH values within a range of the order of 10-15 bins and repeat searching within the averaged values and provide a single most common PBLH appearing with a substantially increased frequency of occurrences within said most common PBLH values as a final PBLH value, or else provide as a final PBLH value this averaged most common PBLH value that most proximally matches customary PBLH values of the specific location and time of the year.

6. Method being performed by the lidar device and operatively associated lidar data processing unit thereof according to claim 5, wherein an auxiliary“DECIS” algorithm is executed to provide said“MUSTI-PBLH” algorithm with infonnation appropriate for selecting one of said classified settings MUSTI-PBLH LOWER ALT- 1, LOWER ALT-2, ..., LOWER ALT-n and one of said plurality of operational modes therein, said“DECIS” algorithm comprising the steps of:

slicing the atmosphere in a plurality of slices of a predetermined thickness, said predetermined thickness being set at 200 m for said lidar device operating at a wavelength of 355 nm;

R

calculating an Aerosol Optical Depth (t) given by ^ (o,^) = j ^'a_c (A,R')dR'

0

and a total Aerosol Optical Depth (x_totai) given by _T,_{M l (}o, R) = _TCI + _TC + ... + r_c

for each slice of the plurality of slices, where, R=200 m and R^' is the distance at which a_aer,cn(U,R') has been measured with a range resolution of 1 bin and corresponds to each different layer within the range of R=200 m, and Ttotai is the sum of xcn, which is the individual optical depth (n=T,2,3,...) of the (n) areas of measuring;

characterizing each slice of the plurality of slices of predetermined thickness of 200 m as a“blurry slice” if, in accordance with the“WEATHER PHENOMENA” algorithm, the Aerosol Optical Depth (x) > 2 TO ³ and as a“non-burry slice” if (x) < 2 T0 ³;

providing said“MUSTI-PBLH” algorithm with information on“blurry” and“non- blurry” slices of the atmosphere thereby guiding selection of one most appropriate of said classified settings MUSTI-PBLH LOWER ALT-1, LOWER ALT-2, ..., LOWER ALT-n and one of said plurality of operational modes therein;

if, assuming an expected range of 210-2600 m for commonly occurring PBLH, thereby having obtained measurements of 12 slices of predetermined thickness of 200 m, providing an estimate of the overall blurriness of the atmosphere in accordance with categories as follows:

a“Crystal Clear Atmosphere” if at least 11 slices have been classified as non- blurry,

a“Clear Atmosphere” if at least 10 slices have been classified as non-blurry, an“Almost Blurry Atmosphere” if at least 5 slices have been classified as blurry;

a“Blurry Atmosphere” if 6-9 slices have been classified as blurry, and an“Extremely Blurry Atmosphere” if t>2*10 ³ in 10-12 of the slices have been classified as blurry.

7. Method being perfonned by the lidar device and operatively associated lidar data processing unit thereof according to claims 1, 5 and 6, wherein a technique designated as FASTPLAN technique is being used to perform slant range measurements for obtaining lidar response signal data, said data being processed by said“MUSTI-L/D” algorithm to acquire atmospheric layering and by said“MUSTI- PBLH” and“DECIS” algorithms to acquire PBLH within a substantially reduced processing time, said FASTPLAN technique comprising the steps of:

Setting said lidar device in a vertical or near vertical direction to detect molecular layer at a range Rref, said range Rref being thereafter retained to indicate a maximum value;

Stepping said lidar device downwardly from said vertical or near vertical direction in angular steps of a predetermined value in a direction, same or opposite to the direction of rotation of the Earth, depending on a prevailing direction of movement of aerial masses and on related meteorological conditions, to acquire sequential packages of values of said return signal parameters (p) in a series of sequential steps, said sequential packages of values of said return signal parameters (p) being thereafter processed by said“MUSTI-L/D” algorithm to acquire atmospheric layering and by said“MUSTI-PBLH” and“DECIS” algorithms to acquire PBLH, wherein, prior to their provision for processing, said values are being brought in conformity with

PBLH=cos )*R, where f is the vertical angle (near vertical)

R = h/cos ) =>h = R*cos )

h_c 3 = G3 * cos( ) = 3 * 7.5(m) * cos(<p)

where, h, is the corresponded height (in vertical measurements) when applying the multi-angle method and h_G3, is an example of vertical appearance, of the original slant one, of the slant‘space filter’ of 3 bins (G3), at slant measurements,

ά(b„{lM _MZDR 1 M_eJ ,R))

dh [cos(_<^) - R * sin(^)] dR

where, the first order derivative of p_aer(L,h) (and/or a_aei<Y,h), RCS) is being calculated versus height and a change in the variable from height (h) to (slant) range (R), depending on angle (f), is being taken into account for M3DR etc. to be produced in slant-near vertical measurements.

8. Method being performed by the lidar device and operatively associated lidar data processing unit thereof according to claim 2, wherein an algorithm designated as “WEATHER PHENOMENA-PBLH” algorithm is being used to acquire PBLH, wherein said“WEATHER PHENOMENA-PBLH” algorithm takes account of the qualitative report of atmospheric layering being produced by said “WEATHER PHENOMENA” algorithm and concentrates on data indicating type of fog existence as a boundary-layer phenomenon and starting from a maximum valuable range-height, wherein said boundary-layer phenomenon occurs through indication of a‘Fog’ and proximal qualitative layers, and moving towards lower ranges-heights where layers present varying fog designations, to acquire values of possible PBLHs, eventually providing the highest of said values of possible PBLHs as the PBLH obtained by said “WEATHER PHENOMENA-PBLH” algorithm.

9. Method being performed by the lidar device and operatively associated lidar data processing unit thereof according to claim 8 and claims 5-6, wherein, said “WEATHER PHENOMENA-PBLH” algorithm is adapted to be used in combination with said“MUSTI-PBLH” and said“DECIS” algorithms being executed for as many signal files as selected for said“MUSTI-PBLH” algorithm,

wherein processing of“MUSTI-PBLH” is initiated with said matrix (B) of“VASPA- L/D” algorithm containing solely those elements of values of said backscattering coefficient p_aer( ,R) or a_aer(^R) that are indicative of type of fog and proximal qualitative layers and/or a medium value from each of the possible PBLHs found from execution of said“WEATHER PHENOMENA-PBLH” algorithm is calculated, said possible PBLHs being within a spectrum of 10 to 15 range bins or within a range of 100- 150m, wherein other possible PBLHs being outside said spectrum of 10 to 15 range bins or range of 100- 150m is discarded from said XLAYER matrix of“MUSTI- PBLH” algorithm and thereafter execution of said“MUSTI-PBLH” algorithm is continued to provide the PBLH.

10. Method being performed by the lidar device and operatively associated lidar data processing unit thereof according to claim 1 , wherein an algorithm designated as “SIBESMEA” algorithm is being used to acquire wind speed data, said“SIBESMEA” algorithm comprising the steps of:

Setting said lidar device pointing in a single vertical or slant direction and selecting a higher or lower signal emission frequency for assessment of estimated higher or lower wind speeds respectively;

Obtaining a plurality of records of sequential lidar response signals from said lidar device emitting in said single vertical or slant direction, the number of records obtained being higher or lower for assessment of estimated higher or lower wind speeds respectively; Placing said records of sequential lidar response signals onto a single plane and calculating the acreage of each of said sequential records of sequential lidar response signals at a predetermined range of a selected number of range bins;

Identifying identically configured parts in said records of sequential lidar response signals and establishing a similarity factor being employed to provide a resemblance condition of consecutive records of said sequential lidar response signals on the basis of comparison of said identically configured parts, wherein said resemblance condition indicates a substantially constant movement of aerial mass within the range being covered in the time elapsing between a pair of sequential lidar response signals;

Initiating comparison between consecutive signals moving backwards from the last record of said plurality of records of sequential lidar response signals and derive an outcome of an average value of wind speed by the ratio of the total range within which consecutive signals satisfy said resemblance condition and the total time that has elapsed during movement of aerial mass within said total range being covered.

1 1. Method being performed by the lidar device and operatively associated lidar data processing unit thereof according to claim 2, wherein an algorithm designated as “WEATHER PHENOMENA- WIND” algorithm is being executed to acquire wind speed data, wherein said “WEATHER PHENOMENA- WIND” algorithm takes account of the qualitative report of atmospheric layering being produced by said “WEATHER PHENOMENA” algorithm and thereafter comprises the steps of:

setting said lidar device pointing at horizontal, vertical or slant from lidar’ s device location, selecting a signal emission frequency thereof and initiating monitoring of atmospheric layers being provided by said“WEATHER PHENOMENA” algorithm by tracking movement of a centerline of each atmospheric layer;

Obtaining a plurality of records of sequential lidar response signals from said lidar device emitting in said single vertical, horizontal or slant direction, said records constituting ‘frames’ of movement of the layers, wherein the number of records required increases as wind speeds to be measured are estimated to have higher values;

Deriving an outcome of an average value of wind speed of each atmospheric layer, by a ratio of the total range being covered and of the time that has elapsed during the period of receipt of sequential lidar response signals, said total range being derived from the displacement of the centerline of each said atmospheric layer as measured between the first and the last of said plurality of records of sequential lidar response signals and said time elapsed being derived by the number of records obtained and the signal emission frequency that has been selected.

Providing an outcome of zero wind speed value for an atmospheric layer being found stationary through said total range being covered and continue monitoring this atmospheric layer, as well as each atmospheric layer being provided by said “WEATHER PHENOMENA” algorithm including monitoring of newly arising atmospheric layers if any.

12. Method being performed by the lidar device and operatively associated lidar data processing unit thereof according to claim 10 or 11, wherein said lidar device is adapted to provide a 2-dimensional and/or 3-dimensional wind speed measurement, wherein a 2-dimensional wind speed measurement provided with said lidar device being arranged to point at a slant direction is analyzed in a wind speed component lying along a first vertical direction and a wind speed component lying along a second horizontal direction, said second direction perpendicularly oriented to said first vertical direction, and

a 3 -dimensional wind speed measurement further provides wind speed measurement along a third direction perpendicular to a plane of said first direction and said second direction of wind speed measurement;

wherein an air mass loss in an atmospheric layer, detected during processing, by means of deterioration of the lidar response signal received therefrom, along said first and said second direction of consecutive measurements of said lidar device indicates movement of said air mass being lost, along said third direction, whereby wind speed along said third direction is calculated through monitoring a rate of said air mass loss and processing said lidar response signal along said third direction;

wherein lidar device is being sequentially oriented in a direction differing from said first vertical, horizontal or slant direction by an incremental change of angle (f) and/or in a direction differing from said third direction by an incremental change of angle (Q), wherein said lidar device is recalibrated during sequential movement at predetermined slant range steps of angle (f) or polar range steps of angle (Q), through setting said extinction coefficient a_aer(T,R) and said backscattering coefficient _aer( ,R) to new calibrating values that retain as reference the immediately precedent values of a_aer( R) and/or or p_aer(T,R), whereby processing of consecutively obtained records is performed to provide a 2-dimensional and/or 3 -dimensional wind speed profile.

13. Method being performed by the lidar device and operatively associated lidar data processing unit thereof according to claim 2, wherein a“VISIBILITY” algorithm is being executed to acquire visibility data along a vertical, horizontal or slant direction from the position of installation of said lidar device, comprising the steps of: setting said lidar device pointing at horizontal, vertical or slant direction and selecting a signal emission frequency thereof;

taking account of the qualitative report of atmospheric layering being produced by said“WEATHER PHENOMENA” algorithm as derived from processing data of said lidar device operating at a wavelength of 355 nm;

R

calculating an Aerosol Optical Depth (t) given by ^ (o , R) = a_c (X,R')dR

⁰

and a total Aerosol Optical Depth (x_totai) given by (o , R) = r_Ci + r_C; + ... + r_c

for each aerosol layer contained in said report of atmospheric layering being produced by said“WEATHER PHENOMENA” algorithm, where, R^' is the distance at which oi_aer,cn(^,R') has been measured with a range resolution of 1 bin and corresponds to each different layer within the maximum range R of the lidar signal, and r_totai is the sum of Te_n, which is the individual optical depth (n=l,2,3,...) of the (n) areas of measuring; calculating the visibility Vis (km) versus the atmospheric extinction coefficient a (knr ^!) being given by:

Vis- 3/a_aer( ,R) (for e=0.05 at 550 nm) or preferably by

Vis= 3.9\2_i ^lcc_aer( ,R) (for e=0.02 at 550 nm)

where, e is a pure number, showing a contrast threshold, as a difference of the self- luminance of any object and the general luminance of the area viewed from a standing position

examining whether the total Aerosol Optical Depth (xtotai)> 2 TO ², said value of (xtotai) corresponding to said extinction coefficient

and to said backscattering coefficient b_36ΐ-(l,II)> 2T0 ³m ¹sr¹ with a lidar ratio C( ,R)=10 sr indicating cumulus cloud according to“WEATHER PHENOMENA” algorithm at 355 nm and if this condition is satisfied provide said range R^' of the aerosol layer being investigated as the value of visibility range from the position of installation of said lidar device or else continue investigation in other aerosol layers within said maximum range R of the lidar signal until the above condition of the total Aerosol Optical Depth (xtotai)> 2 TO ² is being fulfilled and then provide the resulting value of visibility range, whereas following completion of investigation of all aerosol layers within said maximum range R of the lidar signal without said condition of the total Aerosol Optical Depth (xtotai)> 2 TO ² being fulfilled, provide an outcome of visibility being said maximum range R of the lidar signal.

14. Method being performed by the lidar device and operatively associated lidar data processing unit thereof according to claim 13, wherein said “VISIBILITY” algorithm is being executed to acquire visibility data at an airport, said lidar device being oriented at a vertical, horizontal or slant direction and said operatively associated lidar data processing unit providing data of human eye visibility of the airport tower controller in the direction of an aircraft approaching for landing and of visibility of the pilot of the approaching aircraft in the direction of the line of approach of said aircraft to the runway.