An Improved Calibration Method for Photonic Mixer Device Solid-State Array Lidars Based on Electrical Analog Delay

As a typical application of indirect-time-of-flight (ToF) technology, photonic mixer device (PMD) solid-state array Lidar has gained rapid development in recent years. With the advantages of high resolution, frame rate and accuracy, the equipment is widely used in target recognition, simultaneous localization and mapping (SLAM), industrial inspection, etc. The PMD Lidar is vulnerable to several factors such as ambient light, temperature and the target feature. To eliminate the impact of such factors, a proper calibration is needed. However, the conventional calibration methods need to change several distances in large areas, which result in low efficiency and low accuracy. To address the problems, this paper presents an improved calibration method based on electrical analog delay. The method firstly eliminates the lens distortion using a self-adaptive interpolation algorithm, meanwhile it calibrates the grayscale image using an integral time simulating based method. Then, the grayscale image is used to estimate the parameters of ambient light compensation in depth calibration. Finally, by combining four types of compensation, the method effectively improves the performance of depth calibration. Through several experiments, the proposed method is more adaptive to multiscenes with targets of different reflectivities, which significantly improves the ranging accuracy and adaptability of PMD Lidar.


Introduction
Three-dimensional information acquisition has gained extensive attention in the field of computer vision, robot navigation, human-computer interaction, automatic driving, etc. [1]. Generally, the methods to obtain three-dimensional information mainly include stereo vision [2,3], structural light [4], single-pixel 3D imaging [5] and time-of-flight (ToF) [6]. Stereo vision needs advanced matching algorithm to obtain accurate depth information, which is vulnerable to ambient light. Structural light needs projection optimization compensation, which requires high performance of the processing system. Compared with the two methods above, the ToF sensor system utilizes active infrared laser to achieve depth information acquisition, which has the advantages of low cost, high frame frequency and high reliability [7,8].
Photonic mixer device (PMD) solid-state array Lidar, as one typical kind of the ToF sensor system, is widely used in computer vision [9,10]. However, there inevitably exist several errors sources (such as ambient light, integration time, temperature drift and reflectivity), which reduce the performance of the ToF sensor significantly. Hence, the equipment needs to be properly calibrated to achieve reliable depth information acquisition [11].
Signal demodulation is the key step during the working process of the PMD solid-state array Lidar, which is shown in Figure 2. Two different capacitors (CA and CB) with two phase windows (0° and 180°) are set under each pixel of the ToF chip. The differential correlation sampling (DCS) method was used to demodulate the received signal. In general, the sampling number determined the accuracy of the demodulation, while the efficiency could be accordingly influenced. In this paper, the four-step phase-shift method was adopted for sampling. In other words, the process of demodulation was to sample the received signal at four different phases (0°, 90°, 180° and 270°) respectively by using the capacitors of two phase windows, and then suppresses the noise by obtaining the difference between these capacitors. The phase shift of the modulated light was calculated according to the sampling amplitude. At last the target distance was calculated from the phase shift. The specific process of the four-step phase-shift method was performed as follows. The emitted light signal can be presented as ( ) Acos( )  Signal demodulation is the key step during the working process of the PMD solid-state array Lidar, which is shown in Figure 2. Two different capacitors (C A and C B ) with two phase windows (0 • and 180 • ) are set under each pixel of the ToF chip. The differential correlation sampling (DCS) method was used to demodulate the received signal. In general, the sampling number determined the accuracy of the demodulation, while the efficiency could be accordingly influenced. In this paper, the four-step phase-shift method was adopted for sampling. In other words, the process of demodulation was to sample the received signal at four different phases (0 • , 90 • , 180 • and 270 • ) respectively by using the capacitors of two phase windows, and then suppresses the noise by obtaining the difference between these capacitors. The phase shift of the modulated light was calculated according to the sampling amplitude. At last the target distance was calculated from the phase shift.  Signal demodulation is the key step during the working process of the PMD solid-state array Lidar, which is shown in Figure 2. Two different capacitors (CA and CB) with two phase windows (0° and 180°) are set under each pixel of the ToF chip. The differential correlation sampling (DCS) method was used to demodulate the received signal. In general, the sampling number determined the accuracy of the demodulation, while the efficiency could be accordingly influenced. In this paper, the four-step phase-shift method was adopted for sampling. In other words, the process of demodulation was to sample the received signal at four different phases (0°, 90°, 180° and 270°) respectively by using the capacitors of two phase windows, and then suppresses the noise by obtaining the difference between these capacitors. The phase shift of the modulated light was calculated according to the sampling amplitude. At last the target distance was calculated from the phase shift. The specific process of the four-step phase-shift method was performed as follows. The emitted light signal can be presented as ( ) Acos( ) , while the received signal can be presented as Where ω is the modulation frequency, A is the amplitude of the emitted signal, ϕ Δ means the phase shift between the emitted signal and received signal, B is the noise The specific process of the four-step phase-shift method was performed as follows. The emitted light signal can be presented as E(t) = kAcos(ωt), while the received signal can be presented as R(t) = B + kAcos(ωt + ∆ϕ). Where ω is the modulation frequency, A is the amplitude of the emitted signal, ∆ϕ means the phase shift between the emitted signal and received signal, B is the noise signal Sensors 2020, 20, 7329 5 of 22 generated during the transmission of light and k means the signal attenuation coefficient. The sampling process can be expressed in Equation (1): where Q DC i 1 and Q DC i 2 , are the integral values of capacitors C A and C B at sampling point i, respectively.
The distance is calculated by Equation (3):

Analysis of Error Sources of PMD Solid-State Array Lidar
The PMD Lidar is vulnerable to several factors such as internal non-uniformity of the ToF sensor, demodulation process, temperature drift, ambient light, etc. Some of the factors have been discussed in our previous work [42], which will not be mentioned in this paper. However, factors like integration time, temperature drift and DRNU still need to be considered. The analysis of these factors and the qualitative study are described in detail as follows.

Integration Time
Integration time is the time span to output individual data. In general, the integration time can be set from tens to thousands of microseconds. Too short integration time brings the loss of local information, as shown in Figure 3a. While too long integration time will exceed the ToF sensor's receiving range, leading to a local saturation, as shown in Figure 3c. A proper value is needed to capture a sufficient number of photoelectrons without saturation, as shown in Figure 3b.

Analysis of Error Sources of PMD Solid-State Array Lidar
The PMD Lidar is vulnerable to several factors such as internal non-uniformity of the ToF sensor, demodulation process, temperature drift, ambient light, etc. Some of the factors have been discussed in our previous work [42], which will not be mentioned in this paper. However, factors like integration time, temperature drift and DRNU still need to be considered. The analysis of these factors and the qualitative study are described in detail as follows.

Integration Time
Integration time is the time span to output individual data. In general, the integration time can be set from tens to thousands of microseconds. Too short integration time brings the loss of local information, as shown in Figure 3a. While too long integration time will exceed the ToF sensor's receiving range, leading to a local saturation, as shown in Figure 3c. A proper value is needed to capture a sufficient number of photoelectrons without saturation, as shown in Figure 3b.

Temperature Drift
The ToF sensor is susceptible to environment temperature and the heat generated by itself during its working, leading to an uneven temperature distribution. However, the electron mobility in the sensor is temperature dependent. The higher the temperature, the lower the electron mobility, leading to a non-uniformity of measurement, as shown in Figure 4. The ToF sensor is susceptible to environment temperature and the heat generated by itself during its working, leading to an uneven temperature distribution. However, the electron mobility in the sensor is temperature dependent. The higher the temperature, the lower the electron mobility, leading to a non-uniformity of measurement, as shown in Figure 4. In addition, the illumination driver and the external circuit also have temperature dependent demodulation delay, which affect the distance measurement.

Distance Response Non-Uniformity (DRNU)
The ToF sensor typically have many analog to digital converters (ADCs) arranged along the columns of the pixel-field. The different ADCs have slightly different behaviors and result in a non-uniformity between the columns. In addition, this type of error exists in row and due to the non-uniformity of the row addressing signals. These two types of non-uniformities lead to differences of demodulation from pixel to pixel, which is called distance response non-uniformity (DRNU). This type of error also needs to be compensated.
For instance, the phase shift is calculated with Equation (4): While the real phase shift without DRNU compensation is calculated with Equation (5): Figure 5 demonstrates the depth image without DRNU compensation, where the non-uniformity was obviously unneglectable.

Methodology
The conventional calibration methods are time-consuming and complex due to the existence of multiple error sources. Based on the previous work, this paper put forward an improved calibration In addition, the illumination driver and the external circuit also have temperature dependent demodulation delay, which affect the distance measurement.

Distance Response Non-Uniformity (DRNU)
The ToF sensor typically have many analog to digital converters (ADCs) arranged along the columns of the pixel-field. The different ADCs have slightly different behaviors and result in a non-uniformity between the columns. In addition, this type of error exists in row and due to the non-uniformity of the row addressing signals. These two types of non-uniformities lead to differences of demodulation from pixel to pixel, which is called distance response non-uniformity (DRNU). This type of error also needs to be compensated.
For instance, the phase shift is calculated with Equation (4): While the real phase shift without DRNU compensation is calculated with Equation (5): Figure 5 demonstrates the depth image without DRNU compensation, where the non-uniformity was obviously unneglectable. The ToF sensor is susceptible to environment temperature and the heat generated by itself during its working, leading to an uneven temperature distribution. However, the electron mobility in the sensor is temperature dependent. The higher the temperature, the lower the electron mobility, leading to a non-uniformity of measurement, as shown in Figure 4. In addition, the illumination driver and the external circuit also have temperature dependent demodulation delay, which affect the distance measurement.

Distance Response Non-Uniformity (DRNU)
The ToF sensor typically have many analog to digital converters (ADCs) arranged along the columns of the pixel-field. The different ADCs have slightly different behaviors and result in a non-uniformity between the columns. In addition, this type of error exists in row and due to the non-uniformity of the row addressing signals. These two types of non-uniformities lead to differences of demodulation from pixel to pixel, which is called distance response non-uniformity (DRNU). This type of error also needs to be compensated.
For instance, the phase shift is calculated with Equation (4): While the real phase shift without DRNU compensation is calculated with Equation (5): Figure 5 demonstrates the depth image without DRNU compensation, where the non-uniformity was obviously unneglectable.

Methodology
The conventional calibration methods are time-consuming and complex due to the existence of multiple error sources. Based on the previous work, this paper put forward an improved calibration

Methodology
The conventional calibration methods are time-consuming and complex due to the existence of multiple error sources. Based on the previous work, this paper put forward an improved calibration Sensors 2020, 20, 7329 7 of 22 method based on the electrical analog delay. The method fused various error compensations into a comprehensive calibration model, where the grayscale image and the depth image were calibrated jointly, as shown in Figure 6. The lens distortion was corrected using a self-adaptive interpolation algorithm based on Zhang's [49] calibration method. For grayscale image correction, DSNU and PRNU were compensated based on an integration time simulating based method. For depth information correction, the grayscale image was used to estimate the parameters of ambient light compensation. After calculating the raw depth data, the pixel fix pattern noise was eliminated by DRNU error compensation. The temperature drift error was compensated at last.
Sensors 2020, 20, x FOR PEER REVIEW 7 of 22 method based on the electrical analog delay. The method fused various error compensations into a comprehensive calibration model, where the grayscale image and the depth image were calibrated jointly, as shown in Figure 6. The lens distortion was corrected using a self-adaptive interpolation algorithm based on Zhang's [49] calibration method. For grayscale image correction, DSNU and PRNU were compensated based on an integration time simulating based method. For depth information correction, the grayscale image was used to estimate the parameters of ambient light compensation. After calculating the raw depth data, the pixel fix pattern noise was eliminated by DRNU error compensation. The temperature drift error was compensated at last.

Self-adaptive Grayscale Correlation based Depth Calibration Method
Depth Image Calibration  Figure 6. The comprehensive calibration model process.

Lens Distortion Correction
The internal and external parameters of PMD solid-state array Lidar were obtained through Zhang's [49] calibration method, which is not detailed in this paper. Different from the grayscale image, there may exist the holes on the depth image (the received signal is too low to demodulate a valid signal), resulting in the inapplicability of the traditional correction algorithm. A pixel adaptive interpolation strategy we proposed in [42] was utilized in this paper to solve the problem, which is presented in Table 1.

Case
Pixel Adaptive Interpolation Strategy Figure 6. The comprehensive calibration model process.

Lens Distortion Correction
The internal and external parameters of PMD solid-state array Lidar were obtained through Zhang's [49] calibration method, which is not detailed in this paper. Different from the grayscale image, there may exist the holes on the depth image (the received signal is too low to demodulate a valid signal), resulting in the inapplicability of the traditional correction algorithm. A pixel adaptive interpolation strategy we proposed in [42] was utilized in this paper to solve the problem, which is presented in Table 1.

Case Pixel Adaptive Interpolation Strategy
Sensors 2020, 20, x FOR PEER REVIEW 8 of 23

Case Pixel Adaptive Interpolation Strategy
Green points represent the projection of the corrected pixels on the raw image. Purple points Sensors 2020, 20, x FOR PEER REVIEW 8 of 23

Case Pixel Adaptive Interpolation Strategy
Green points represent the projection of the corrected pixels on the raw image. Purple points are the surrounding pixels of the green ones, in which the solid purples mean the real distance pixel Sensors 2020, 20, x FOR PEER REVIEW 8 of 23

Case Pixel Adaptive Interpolation Strategy
Green points represent the projection of the corrected pixels on the raw image. Purple points are the surrounding pixels of the green ones, in which the solid purples mean the real distance pixel Sensors 2020, 20, x FOR PEER REVIEW 8 of 23

Case Pixel Adaptive Interpolation Strategy
Green points represent the projection of the corrected pixels on the raw image. Purple points are the surrounding pixels of the green ones, in which the solid purples mean the real distance pixel while the hollow ones represent the holes. D0, D1, D2, D3 and D4 denote the pixels to be interpolated

Case Pixel Adaptive Interpolation Strategy
Green points represent the projection of the corrected pixels on the raw image. Purple points are the surrounding pixels of the green ones, in which the solid purples mean the real distance pixel while the hollow ones represent the holes. D0, D1, D2, D3 and D4 denote the pixels to be interpolated and its lower-left pixel, top-left pixel, lower-right pixel and top-right pixel, respectively. (xp0,yp0), Green points represent the projection of the corrected pixels on the raw image. Purple points are the surrounding pixels of the green ones, in which the solid purples mean the real distance pixel while the hollow ones represent the holes. D 0 , D 1 , D 2 , D 3 and D 4 denote the pixels to be interpolated and its lower-left pixel, top-left pixel, lower-right pixel and top-right pixel, respectively. (x p0 ,y p0 ), (x p1 ,y p1 ), (x p2 ,y p2 ), (x p3 ,y p3 ) and (x p4 ,y p4 ) are the coordinates of the pixels to be interpolated and its lower-left pixel, top-left pixel, lower-right pixel and top-right pixel, respectively. α x = x p0 − x p1 , α y = y p0 − y p1 and (u,v) are the coordinates of the pixels to be interpolated under a barycentric coordinate system.

Grayscale Image Calibration
A grayscale image is acquired by TOF chip when PMD solid-state Array Lidar works in the passive mode. The acquisition process is consistent with the intensity image obtained by traditional complementary metal oxide semiconductor (CMOS) chip, which can be used to characterize the intensity of ambient light.
In general, the grayscale image is vulnerable to DSNU and PRNU, where DSNU represents the differences of gray values between pixels captured under the dark condition and PRNU represents the differences of gray values between pixels captured under the common condition, respectively.
The most commonly used approach to correct the influence of DSNU and PRNU has a standardized process, which can be found in european machine vision association (EMVA) standard 1288 [50]. In this approach, one dark image and one bright image are captured under the same exposure condition. DSNU and PRNU are then calculated from the images. However, there exist some limitations of the approach. For instance, the images are captured under a specific condition, leading to bad applicability. The ambient light is set artificially, which introduces extra error. Based on this approach, an integration time simulating based method was proposed in this paper. The main contributions of the method mainly include: (I) Instead of setting the ambient light artificially, the levels of ambient light were simulated by setting the integration times and (II) the non-uniform of the exposure was eliminated by calculating the mean value of multiple ambient light levels. The process of the method is given as follows.
(1). Set the integration time to 0 µs to simulate the dark condition. Collect N = 100 frames of grayscale images and calculate the mean value.
(2). Change the integration time to simulate different ambient light levels. Collect N = 100 frames of grayscale images under amplitudes of 10%, 30%, 50% and 80% respectively. Similarly, the mean values with different amplitudes are obtained.
(3). Calculate the spatial variances under different ambient levels. Spatial variance is simply an overall measure of the spatial nonuniformity, which is helpful to estimate DSNU and PRNU.
where Y dark−avg is the mean value of grayscale images under the dark condition, Y AL−avg is the mean value of grayscale images under ambient light, N means the number of frames, S 2 dark and S 2 AL are spatial variances under dark and ambient light conditions, respectively, b DSNU is the offset of DSNU, k PRNU is the gain of PRNU and I corr (x, y) is the compensation value of grayscale image after calibration.

Ambient Light Compensation
Due to its special structure, the PMD solid-state array Lidar has the ability to obtain a grayscale image and depth image simultaneously. Meanwhile, the amplitude of a grayscale image has a close relationship with ambient light. Based on this, the grayscale image was used to estimate the parameters of ambient light compensation in depth calibration, which is the self-adaptive grayscale correlation based depth calibration method (SA-GCDCM). The basic idea of the method is to eliminate the influence of ambient light in the sampling stage by introducing an ambient light correction factor K AL . The factor K AL is calculated from several DCs sampled under different ambient light levels. The ambient light is controlled accurately by adjusting the integration time. K AL is then utilized to correct the DCs in the real sampling process. There exists internal noise and external error during the sampling. The errors are corrected in the calculation. The process of the method is concluded as (all the following measurements use the spatial average of the region of interest (ROI) within the coordinates (100,70) and (220,165) and the temporal average of 100 frames as a default): (1). Turn the ambient light on and record the amplitude of the grayscale image as Q gray . (2). Change the amplitude of the grayscale image to 0.5 times that of Q gray by adjusting the integration time. Measure the DC0/2 and record as DC0 setting1 and DC2 setting1 , respectively. (3). Change the amplitude of the grayscale image to 1.5 times that of Q gray by adjusting the integration time. Measure the DC0/2 and record as DC0 setting2 and DC2 setting2 , respectively. (4). Turn the ambient light off and measure the DC0/2, which are recorded as DC0 no and DC2 no , respectively. (5). Calculate four measurements, Q0 1 , Q0 2 , Q2 1 and Q2 2 . (6). Correct the errors generated in the sampling. There inevitably exists internal noise and external error. The internal noise mainly comes from the internal circuit and can be eliminated by subtracting two samples at the same phase. The external error mainly comes from the instability of the environment. It can be suppressed by calculating the mean value of the samples.
As introduced in Section 2.1, distance measurements are taken by acquiring the four DCs and calculated pixel-by-pixel during runtime as Equation (15): The equation is revised after ambient light compensation as Equation (16): where DC1 corr (x, y) and DC0 corr (x, y) are sample values corrected by K AL . Through SA-GCDCM, the disturbance of ambient light could be effectively eliminated. Compared with the traditional joint method using a common RGB camera with a ToF camera, this method has no requirement of coordinate transformation and feature matching, leading to a better data consistency and self-adaptability.

Demodulation Error Correction
In general, the sinusoidal wave is adopted as the modulated continuous wave signal, which can be represented as E(t) = kAcos(ωt). Similarly, the received signal is deemed as a sinusoidal wave in demodulation as well. However, because of the limitations of the generator bandwidth, the actual received signal is similar to a rectangular wave [42], as shown in Figure 7. Therefore, the rectangular wave was used for demodulation analysis in this paper.
Sensors 2020, 20, x FOR PEER REVIEW 10 of 22 (7). The ambient light correction factor KAL is calculated by Equation (13): where DC0/1/2/3 are the sample values acquired at 0°, 90°, 180° and 270° respectively. KAL is used to compensate the ambient light error during the demodulation process as Equation (14) As introduced in Section 2.1, distance measurements are taken by acquiring the four DCs and calculated pixel-by-pixel during runtime as Equation (15): a t a n 2 2 2 2 , 0 , The equation is revised after ambient light compensation as Equation (16): a t a n 2 2 2 2 , 0 , where ( ) Through SA-GCDCM, the disturbance of ambient light could be effectively eliminated. Compared with the traditional joint method using a common RGB camera with a ToF camera, this method has no requirement of coordinate transformation and feature matching, leading to a better data consistency and self-adaptability.

Demodulation Error Correction
In general, the sinusoidal wave is adopted as the modulated continuous wave signal, which can be represented as ( ) Acos( ) Similarly, the received signal is deemed as a sinusoidal wave in demodulation as well. However, because of the limitations of the generator bandwidth, the actual received signal is similar to a rectangular wave [42], as shown in Figure 7. Therefore, the rectangular wave was used for demodulation analysis in this paper. Different from Section 2.1, the sampling process is modified as Equation (17): Different from Section 2.1, the sampling process is modified as Equation (17): Sensors 2020, 20, 7329

of 22
Thus, an extra error, which is called the demodulation error is generated in the revising of the sampling process. The method to correct demodulation error has been discussed extensively in [42], which will not be introduced in this paper.

DRNU Error Compensation
The calibration was based on the electrical analog delay method. Instead of changing the real distance, delay-locked loop was used to simulate the distance in this method. The simulated distance is composed of two parts, as shown in Equation (18). The first part is the simulated distance of DLLs. The system contains several DLLs and each DLL represents a specific simulated distance, e.g., 0.3 m. The second part represents the real distance between the calibration plate and the PMD Lidar. Through combining the two parts, multiple distances could be simulated without moving the calibration plate. However, DLL is susceptible to temperature changing and circuit delay, leading to a deviation between simulated distance and the set distance. Thus, a DRNU error compensation was conducted as Equation (19): where d DLL represents the simulated distance of a single DLL, n is the number of DLLs, o zero is the distance between the PMD solid-state array Lidar and the reflecting plate, D sim (x,y) represents the overall simulated distance, D cal (x,y) represents the corrected distance after compensation and DRNU(x,y) is the compensation value.
Since the DRNU error is related with distance, and the limited number of compensate values cannot completely cover the whole distance. A linear interpolation was carried out to obtain a continuous offset curve. The interpolate method is quite basic and will not be illustrated here.

Temperature Compensation
Several research have reported the influence of the temperature for the ToF camera [7,15,35]. In this paper, the main components related to temperature error in the PMD Lidar were further classified into three parts, which are the ToF sensor, the illumination driver and the external circuit, as discussed in Section 2.2.2. It was found from experiments that the error showed a linear relation with temperature, from which the higher the temperature, the higher the error was. Since the error arose from temperature drift was compensated with a joint equation, as shown in Equation (20): where D final (x,y) is the corrected distance after temperature compensation, T act represents the acting temperature, T cal means the temperature during the calibration, TC pix is the temperature coefficient of the pixel, TC laser means the temperature coefficient of the illumination unit and TC DLL represents the temperature coefficient of DLL stage. For the device in this paper, the TC pix was 11.3 mm/K, the TC laser was 2.7 mm/K and the TC DLL was 0.7 mm/K. It is worth mentioning that the parameters were obtained by the specific device, which means the parameters are not applicable for each device. This is due to the diversity of the circuit board, chip and other components derived from fabrication.
After utilizing the temperature compensation, the multiscene adaptive multifactor calibration model was established. Then the feasibility of the model was verified by experiments in Section 4.

Experimental Settings
The PMD solid-state array Lidar is shown in Figure 8, which is mainly composed of four parts: the emitting unit, the receiving unit, the processing unit and the transmission unit. The emitting unit generates and emits the NIR light with the VCSEL. The receiving unit receives the returned light with a CMOS sensor and converts the optical signal into an electrical signal. The processing unit calculates the distance data by demodulating the phase delay between the emitted and the detected signal. The transmission unit transmits the distance data to the computer. was 0.7 mm/K. It is worth mentioning that the parameters were obtained by the specific device, which means the parameters are not applicable for each device. This is due to the diversity of the circuit board, chip and other components derived from fabrication.
After utilizing the temperature compensation, the multiscene adaptive multifactor calibration model was established. Then the feasibility of the model was verified by experiments in Section 4.

Experimental Settings
The PMD solid-state array Lidar is shown in Figure 8, which is mainly composed of four parts: the emitting unit, the receiving unit, the processing unit and the transmission unit. The emitting unit generates and emits the NIR light with the VCSEL. The receiving unit receives the returned light with a CMOS sensor and converts the optical signal into an electrical signal. The processing unit calculates the distance data by demodulating the phase delay between the emitted and the detected signal. The transmission unit transmits the distance data to the computer. As illustrated in Figure 9, the experimental settings were established. The grayscale image calibration system, as shown in Figure 9a, was used to calibrate the lens distortion and eliminate the effect of DSNU and PRNU. The depth image calibration system, as shown in Figure 9b, was used to compensate multierror sources in distance measurements. As illustrated in Figure 9, the experimental settings were established. The grayscale image calibration system, as shown in Figure 9a, was used to calibrate the lens distortion and eliminate the effect of DSNU and PRNU. The depth image calibration system, as shown in Figure 9b, was used to compensate multierror sources in distance measurements. was 0.7 mm/K. It is worth mentioning that the parameters were obtained by the specific device, which means the parameters are not applicable for each device. This is due to the diversity of the circuit board, chip and other components derived from fabrication.
After utilizing the temperature compensation, the multiscene adaptive multifactor calibration model was established. Then the feasibility of the model was verified by experiments in Section 4.

Experimental Settings
The PMD solid-state array Lidar is shown in Figure 8, which is mainly composed of four parts: the emitting unit, the receiving unit, the processing unit and the transmission unit. The emitting unit generates and emits the NIR light with the VCSEL. The receiving unit receives the returned light with a CMOS sensor and converts the optical signal into an electrical signal. The processing unit calculates the distance data by demodulating the phase delay between the emitted and the detected signal. The transmission unit transmits the distance data to the computer. As illustrated in Figure 9, the experimental settings were established. The grayscale image calibration system, as shown in Figure 9a, was used to calibrate the lens distortion and eliminate the effect of DSNU and PRNU. The depth image calibration system, as shown in Figure 9b, was used to compensate multierror sources in distance measurements.  The grayscale image calibration system mainly includes the checkerboard, the PMD solid-state array Lidar, the clamping device and the rail. Based on the system, the grayscale image calibration was conducted as follows.
(1). Clamp the PMD Lidar on the clamping device. (2). Adjust the clamping device to a proper location where the checkerboard is suitable in size and position in the field of view of the PMD Lidar. (3). Calibrate the grayscale image based on integration time simulating. (4). Obtain several grayscale images with checkerboard in different directions to calibrate the lens distortion. (5). Utilize the pixel adaptive interpolation strategy to fill the holes. The depth image calibration system mainly includes the reflecting plate, the PMD solid-state array Lidar, the cylindrical tube, the ambient light source, the clamping device and the rail. The cylindrical tube was used to protect the ToF sensor from affecting the stray light. The ambient light source was used to provide the assistant lighting. Based on the system, the depth image calibration was conducted as follows.
(1). Install the cylinder on the PMD Lidar and clamp the PMD Lidar on the clamping device. (2). Adjust the clamping device to a proper location where the quality of the light spot projected on the reflecting plate is optimized. (3). Change the distance with the electrical analog delay method to perform the depth calibration.
Multiple error compensation is included in this step. (4). Change the reflecting board to adjust the method with objects of different reflectivities. (5). Conduct the interpolation on the data to obtain the continuous offset curves.

Lens Distortion Correction
Several grayscale images with checkerboard in different directions were obtained to calibrate the lens distortion. The pixel adaptive interpolation strategy was used to fill the holes. The results are shown Figure 10.
(1). Clamp the PMD Lidar on the clamping device. (2). Adjust the clamping device to a proper location where the checkerboard is suitable in size and position in the field of view of the PMD Lidar. (3). Calibrate the grayscale image based on integration time simulating. (4). Obtain several grayscale images with checkerboard in different directions to calibrate the lens distortion. (5). Utilize the pixel adaptive interpolation strategy to fill the holes. The depth image calibration system mainly includes the reflecting plate, the PMD solid-state array Lidar, the cylindrical tube, the ambient light source, the clamping device and the rail. The cylindrical tube was used to protect the ToF sensor from affecting the stray light. The ambient light source was used to provide the assistant lighting. Based on the system, the depth image calibration was conducted as follows.
(1). Install the cylinder on the PMD Lidar and clamp the PMD Lidar on the clamping device. (2). Adjust the clamping device to a proper location where the quality of the light spot projected on the reflecting plate is optimized. (3). Change the distance with the electrical analog delay method to perform the depth calibration.
Multiple error compensation is included in this step. (4). Change the reflecting board to adjust the method with objects of different reflectivities. (5). Conduct the interpolation on the data to obtain the continuous offset curves.

Lens Distortion Correction
Several grayscale images with checkerboard in different directions were obtained to calibrate the lens distortion. The pixel adaptive interpolation strategy was used to fill the holes. The results are shown Figure 10. The internal parameters and the distortion coefficients are shown in Table 2.  The internal parameters and the distortion coefficients are shown in Table 2. The influences of DSNU and PRNU were eliminated based on the integration time simulating method introduced in Section 3.2. Firstly, the raw grayscale images were acquired under multiple ambient light levels through setting the integration time in several values. Then the spatial variances were calculated from the images to estimate the dark signal non-uniformity (DSNU) and photo response non-uniformity (PRNU). At last the influence of DSNU and PRNU were eliminated by a correction algorithm.

Internal Parameters
To evaluate the effectiveness of the method, several experiments were conducted in different scenes. A checkerboard and a flat white board were used as the test scenes to conduct a qualitative analysis and a quantitative analysis. Then gray images were captured under two real scenes to verify the feasibility of the method. The results are shown in Figure 11.
The images in the left column were captured before calibration, while the images in the right column were captured after calibration. Figure 11a,b were captured with the checkerboard. Compared with the image before calibration, two types of non-uniformity were compensated. In the raw image, the central area was brighter while the surroundings were darker because of the uneven exposure. Meanwhile, there existed distinct light and dark stripes in vertical. In the images after calibration, these two phenomena were suppressed obviously.
To better prove the effectiveness of the method, a quantitative analysis was then conducted. A flat white board was suitable to conduct the analysis because of its good flatness and smoothness. The images are shown in Figure 11c,d. The improvement in visual was consistent with results of the checkerboard. Two types of non-uniformities were compensated obviously. To better verify the improvement, mean value, root mean square error (RMSE) and peak signal to noise ratio (PSNR) [51] were chosen to evaluate the quality of the images, and the results are shown in Table 3. The mean value shows little difference before and after the calibration, which means the method did not change the overall sampling of the grayscale signal. However, the RMSE got a significant reduction after calibration, which indicates the uniformity of the grayscale signal was improved distinctly. In addition, the PSNR improved after calibration, which means the noise derived from PRNU and DSNU was reduced.
Experiments were then conducted in two real scenes to verify the applicability of the method in reality, as shown in Figure 11e-h. It can be seen that the vertical stripes were effectively suppressed after calibration. Meanwhile, the uneven exposure, which leads to uneven brightness of the image, was obviously suppressed. Similarly, quantitative analysis was conducted and the results are shown in Table 3. The mean values showed no distinct differences before and after the calibration, while the reduction of the RMSEs was obvious. It means that the non-uniformity of the grayscale images was suppressed. The PSNRs were higher after calibration, which means the noise was effectively reduced. The results in real scenes were in accordance with results in test scenes, which mean the calibration method is feasible in reality. correction algorithm.
To evaluate the effectiveness of the method, several experiments were conducted in different scenes. A checkerboard and a flat white board were used as the test scenes to conduct a qualitative analysis and a quantitative analysis. Then gray images were captured under two real scenes to verify the feasibility of the method. The results are shown in Figure 11. The images in the left column were captured before calibration, while the images in the right column were captured after calibration. Figure 11a,b were captured with the checkerboard. Compared with the image before calibration, two types of non-uniformity were compensated. In the raw image, the central area was brighter while the surroundings were darker because of the uneven exposure. Meanwhile, there existed distinct light and dark stripes in vertical. In the images after calibration, these two phenomena were suppressed obviously.
To better prove the effectiveness of the method, a quantitative analysis was then conducted. A Figure 11. Results of grayscale image calibration. Images (a, c, e, g) are captured before calibration, while images (b, d, f, h) are captured after calibration.

Result with Depth Image Calibration
Depth image calibration was carried out after grayscale image calibration. The results are shown in Figure 12. In the left image, which was obtained before calibration, there existed many incorrect even invalid data points, while the bright and dark stripes can be observed distinctly. The non-uniformity of the whole image indicates the depth data was untrustworthy before calibration. After utilizing ambient light compensation, demodulation error correction, DRNU error compensation and temperature compensation, the quality of the depth image improved significantly. The number of noise points reduced obviously. The confidence of depth information was greatly improved.
Sensors 2020, 20, x FOR PEER REVIEW 16 of 22 non-uniformity of the whole image indicates the depth data was untrustworthy before calibration. After utilizing ambient light compensation, demodulation error correction, DRNU error compensation and temperature compensation, the quality of the depth image improved significantly. The number of noise points reduced obviously. The confidence of depth information was greatly improved.

Ranging Accuracy Verification under Real Environment
To verify the ranging accuracy and the adaptability of the calibration method, several tests under the real environment were conducted. The test system is shown in Figure 13. Firstly, the test system was placed in the dark environment (the ambient light was about 0 Lux), indoor environment (about 500 Lux) and outdoor environment (about 1200 Lux), respectively. Then the reflecting plate was set as 80% reflectivity, 50% reflectivity and 20% reflectivity, respectively. In each test, the distance between the PMD solid-state array Lidar and the reflecting plate changed from 0.5 to 5 m in a gradient. The mean value in ROI (1000 pixels in the central region) was recorded as the measured distance. Then the distance error was calculated. The test results are illustrated in Figure 14.

Ranging Accuracy Verification under Real Environment
To verify the ranging accuracy and the adaptability of the calibration method, several tests under the real environment were conducted. The test system is shown in Figure 13. Firstly, the test system was placed in the dark environment (the ambient light was about 0 Lux), indoor environment (about 500 Lux) and outdoor environment (about 1200 Lux), respectively. Then the reflecting plate was set as 80% reflectivity, 50% reflectivity and 20% reflectivity, respectively.

Ranging Accuracy Verification under Real Environment
To verify the ranging accuracy and the adaptability of the calibration method, several tests under the real environment were conducted. The test system is shown in Figure 13. Firstly, the test system was placed in the dark environment (the ambient light was about 0 Lux), indoor environment (about 500 Lux) and outdoor environment (about 1200 Lux), respectively. Then the reflecting plate was set as 80% reflectivity, 50% reflectivity and 20% reflectivity, respectively. In each test, the distance between the PMD solid-state array Lidar and the reflecting plate changed from 0.5 to 5 m in a gradient. The mean value in ROI (1000 pixels in the central region) was recorded as the measured distance. Then the distance error was calculated. The test results are illustrated in Figure 14. In each test, the distance between the PMD solid-state array Lidar and the reflecting plate changed from 0.5 to 5 m in a gradient. The mean value in ROI (1000 pixels in the central region) was recorded as the measured distance. Then the distance error was calculated. The test results are illustrated in Figure 14. Compared with the method in [42], the method proposed in this paper effectively reduced the error in the distance range of 0.5-5 m. Meanwhile, the adaptivity under different environments improved a lot. Table 4 gives the detailed performance comparison results of the two methods. The maximal error, average error and RMSE were chosen to compare the detailed performance. From Table 4, the maximum error was reduced distinctly in the proposed method, while the non-uniformity reduced a lot, too. Though the difference of average error was not distinct as other two indicators, the proposed method had better adaptability. In other words, the proposed was more adaptive to multiscene and different reflectivities.
The proposed method was compared with several traditional methods as well, and the result is shown in Table 5. It can be obviously figured out that the proposed method has better performance on range accuracy compared with methods in [13,17,26]. Although the mean distance error shows no distinct improvement compared with method in [19], the proposed method has better Compared with the method in [42], the method proposed in this paper effectively reduced the error in the distance range of 0.5-5 m. Meanwhile, the adaptivity under different environments improved a lot. Table 4 gives the detailed performance comparison results of the two methods. The maximal error, average error and RMSE were chosen to compare the detailed performance. From Table 4, the maximum error was reduced distinctly in the proposed method, while the non-uniformity reduced a lot, too. Though the difference of average error was not distinct as other two indicators, the proposed method had better adaptability. In other words, the proposed was more adaptive to multiscene and different reflectivities.
The proposed method was compared with several traditional methods as well, and the result is shown in Table 5. It can be obviously figured out that the proposed method has better performance on range accuracy compared with methods in [13,17,26]. Although the mean distance error shows no distinct improvement compared with method in [19], the proposed method has better performance on calibration time and scene scope. In addition, the results of experiments expressed that the proposed method had an outstanding performance on adaptability.

Conclusions
To improve the range accuracy of the PMD solid-state array Lidar, this paper presents a self-adaptive grayscale correlation based depth calibration method (SA-GCDCM) based on electrical analog delay. Based on the characteristic of the PMD solid-state array Lidar, the grayscale image was used to estimate the parameters of ambient light compensation in depth calibration. To obtain uniform and stable grayscale image, an integration time simulating based method was proposed for eliminating the influence of DSNU and PRNU. Combining SA-GCDCM and demodulation error correction, DRNU error compensation and temperature compensation, a comprehensive, multiscene adaptive multifactor calibration model was established. A series of experiments were conducted to verify the ranging accuracy and the adaptability of the method. Compared with the prior work, the maximum error has reduced distinctly, meanwhile the RMSE was reduced as well, indicating the proposed method had better accuracy and adaptability, respectively. Compared with the traditional methods, the proposed method had better performance on range accuracy and calibration time and scene scope. The proposed method was more adaptive to multiscenes with targets of different reflectivities, which significantly improved the ranging accuracy and adaptability of PMD Lidar.
Author Contributions: X.W. and P.S. developed the improved calibration method and procedure for PMD solid-state array Lidar. X.W. and W.Z. established the calibration system and verification test platform, meanwhile, analyzed the experiment results. All authors have read and agreed to the published version of the manuscript.

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