remote Constructing a New Inter-Calibration Method for DMSP-OLS and NPP-VIIRS Nighttime Light

: The anthropogenic nighttime light (NTL) data that are acquired by satellites can characterize the intensity of human activities on the ground. It has been widely used in urban development assessment, socioeconomic estimate, and other applications. However, currently, the two main sensors, Defense Meteorological Satellite Program’s Operational Linescan System (DMSP-OLS) and Suomi National Polar-orbiting Partnership Satellite’s Visible Infrared Imaging Radiometer Suite (NPP-VIIRS), provide inconsistent data. Hence, the application of NTL for long-term analysis is hampered. This study constructed a new inter-calibration method for DMSP-OLS and NPP-VIIRS nighttime light to solve this problem. First, NTL data were processed to obtain vicarious site across China. By comparing di ﬀ erent candidate models, it is discovered the Biphasic Dose Response (BiDoseResp) model, which is a weighted combination of sigmoid functions, can best perform the regression between DMSP-OLS and logarithmically transformed NPP-VIIRS. The coe ﬃ cient of determination of BiDoseResp model reaches 0.967. It’s residual sum of squares is 6.136 × 10 5 , which is less than 6.199 × 10 5 of Logistic function. After obtaining the BiDoseResp-calibrated VIIRS (BDRVIIRS), we smoothed it by a ﬁlter with optimal parameters to maximize the consistency. The result shows that the consistency of NTL data is greatly enhanced after calibration. In 2013, the correlation coe ﬃ cient between DMSP-OLS and original NPP-VIIRS data in the China region is only 0.621, while that reaches to 0.949 after calibration. Finally, a consistent NTL dataset of China from 1992 to 2018 was produced. When compared with the existing methods, our method is applicable to the full dynamic range of DMSP-OLS. Besides, it is more suitable for country or larger scale areas. It is expected that this method can greatly facilitate the development of research that is based on the historical NTL archive. in the dynamic range is consistent. All of these indicate that the consistency between BDRVIIRS and DMSP-OLS is strong.


Introduction
Human activities are gradually increasing, which have significantly affected the human environment. The nighttime light (NTL) retrieved by satellite provides a great opportunity to detect human activities on a wide scale and for a long period. Currently, the widely used NTL data is mainly obtained by Defense Meteorological Satellite Program's Operational Linescan System (DMSP-OLS) and the National Polar-orbiting Partnership Satellite's Visible Infrared Imaging Radiometer Suite (NPP-VIIRS). National Oceanic and Atmospheric Administration (NOAA) released the DMSP-OLS annual NTL mean product 2 of 15 from 1992 to 2013, and the NPP-VIIRS monthly NTL mean product from April 2012 to the present. The NTL data that were retrieved through these two sensors have been widely used in multiple research areas, which characterize human activities from different perspectives. The Evaluation of socioeconomic indices based on NTL can cross administrative boundary and provide reference for the areas without census data [1][2][3]. Urban mapping and impervious surface identifying are boosted by NTL due to its nature of portrayal of urban structure [4][5][6][7]. In addition, NTL provides a new opportunity for estimating energy consumption and greenhouse gas emissions [8][9][10][11]. It also promotes research on the distribution of light pollution and its impact on human health and protected areas [12][13][14].
Despite the wide applications of NTL, only the data of single sensor are used in most existing studies. Few works connect the two sensors. The main barriers include 1) the lack of on-board calibration for DMSP-OLS, 2) the difference in overpass time, and 3) the significant overglow effect of DMSP-OLS [15][16][17]. However, the longer-term analysis and application of NTL can be deeply facilitated by connecting the two sensors. Accordingly, the inter-calibration between these two sensors is now one of the most eagerly awaited researches [18][19][20].
Until now, just a few studies have tried to make the inter-calibration between DMSP-OLS and VPP-VIIRS. In the earliest study, taking some NPP-VIIRS radiance values as reference, a power function was used to calibrate the DMSP-OLS data [21]. However, they used daily data and their method was not suitable for temporally average product. In the following researches, the calibration of temporally average NTL data has been focused on. Jeswani used a logarithm and linear function to perform the inter-calibration [22], yet the overglow effect of DMSP-OLS [23,24] was not considered in their work. In the study of the Syrian Civil War, a power function combined with Gaussian low pass filter was used to make an inter-calibration between monthly average NTL data [25], while the dynamic range after calibration is incomplete (0 to 50). Wu et al. [26] adopted a similar method to correct the annual NTL data of Beijing region, while their result has similar problems with that of Li et al. [25]. The city-scale NPP-VIIRS data were converted into DMSP-like data from 1996 to 2017 by a residuals-corrected geographically weighted regression model [27], while the time range of NTL was not completely expanded [20]. All of these works have promoted the applications of NTL, but these functions still have some limitations in different aspects. A more robust function is expected to achieve better calibration.
The objective of this study is to construct a better calibration function to improve the performances in terms of spatial scale, dynamic range, time range, and simplicity of application, etc. The study area (the whole region of China) and NTL data are introduced in Section 2. Subsequently, the methods of calibration, including a flexible function-BiDoseResp model and the way to find the optimal parameters for smoothing NPP-VIIRS are given in Section 3. Subsequently, Section 4 shows the results of calibration. Finally, in Section 5, calibration performance and the differences between existing methods are evaluated and discussed in detail.

Study Areas and Datasets
The whole region of China is taken as the research area. Since the 1980s, China has experienced rapid economic development and urban expansion every year. Studies have shown that NTL in China have significantly changed since 1992 [28][29][30]. Therefore, it is particularly important to establish a better calibration function for DMSP-OLS and NPP-VIIRS nighttime light of China region.
DMSP-OLS was originally designed to detect the presence of night clouds, so it is not equipped with on-board calibration. When the sky is cloudless, city illumination can be captured and recorded by the sensor. However, it keeps changing gain coefficients according to the dynamic scene luminosity. So, the gain coefficients are not recorded systematically [31]. Furthermore, due to its large footprint, there is an obvious overglow effect which causes the lit area in the data larger than actual light source area. The annual DMSP-OLS stable light product used in this study was downloaded from the website of Earth Observation Group (EOG) affiliated with NOAA (https://ngdc.noaa.gov/eog/dmsp/ downloadV4composites.html). It is available from 1992 to 2013, with a spatial resolution of 30 arc On the other hand, the NPP-VIIRS Day/Night Band (DNB) is dedicated to acquiring NTL data on the ground. It adopts a solar diffuser for on-board calibration. Three different gains stages are used to record the ground illumination under different bright conditions. In the processing phase, according to known coefficients, raw data are converted into radiance. In addition, the overglow effect in the NPP-VIIRS data is effectively reduced by aggregating and deleting stretched pixels [32]. The NPP-VIIRS monthly product from April 2012 to the present is also available on the website of EOG (https://ngdc.noaa.gov/eog/viirs/download_dnb_composites.html). It records the radiance of ground illumination. Its unit is nanoWatts/cm 2 /Sr. The spatial resolution is 15 arc seconds, the quantization ability is 14 bits for high gain stage (HGS) and 13 bits for medium/low gain stage (MGS/LGS). The overpass time is 1:30 [33]. In this study, the NTL data of these two sensors in the overlapping year of 2013 were used to construct the calibration function. Table 1 shows the properties of DMSP-OLS and NPP-VIIRS.

Methodology
It is assumed that NPP-VIIRS data are the record of actual ground radiation. The variation of gain coefficient and the overglow effect are regarded as DMSP satellite's non-linear transformation to actual ground radiation. Reconstructing the non-linear transformation can be divided into two parts: (1) constructing a robust calibration function from NPP-VIIRS to DMSP-OLS, and (2) smoothing the calibrated VIIRS to adapt to the overglow effect of DMSP-OLS. Figure 1 shows the flowchart of this study.
Number (DN). Its quantization ability is six bits. The overpass time of the satellite is between 19:30 and 21:30.
On the other hand, the NPP-VIIRS Day/Night Band (DNB) is dedicated to acquiring NTL data on the ground. It adopts a solar diffuser for on-board calibration. Three different gains stages are used to record the ground illumination under different bright conditions. In the processing phase, according to known coefficients, raw data are converted into radiance. In addition, the overglow effect in the NPP-VIIRS data is effectively reduced by aggregating and deleting stretched pixels [32]. The NPP-VIIRS monthly product from April 2012 to the present is also available on the website of EOG (https://ngdc.noaa.gov/eog/viirs/download_dnb_composites.html). It records the radiance of ground illumination. Its unit is nanoWatts/ cm 2 /Sr. The spatial resolution is 15 arc seconds, the quantization ability is 14 bits for high gain stage (HGS) and 13 bits for medium/low gain stage (MGS/LGS). The overpass time is 1:30 [33]. In this study, the NTL data of these two sensors in the overlapping year of 2013 were used to construct the calibration function. Table 1 shows the properties of DMSP-OLS and NPP-VIIRS.

Preprocessing
The NTL data of both sensors were re-projected to the Albers Equal Area Conic Projection. An annual synthesis of NPP-VIIRS data was carried out since NPP-VIIRS data are monthly mean product, while DMSP-OLS data are annual mean product. Subsequently, it was resampled to 1 km × 1 km to match the DMSP-OLS data. In addition, there are temporary illumination and background noise in the NPP-VIIRS data; instead, only long-term and stable light exists in DMSP-OLS data. A simple method that was proposed by Li et al. [34] was adopted to eliminate the temporary illumination and background noise. After that, a small number of negative pixels were set to 0. Finally, the method that was proposed by Elvidge et al. [35] was used to improve the continuity of the DMSP-OLS data from 1992 to 2013.
A vicarious site is required in inter-calibration to provide samples [36][37][38]. However, there is no special calibration site designed for NTL in the world. This study used the coefficient of variation (CV) to select the calibration site to deal with this problem, referring to Jeswani [22]. It is assumed that the small CV stands for the stable NTL radiation. First, the CVs of 3 × 3 windows of all pixels in DMSP-OLS data and NPP-VIIRS data were computed, respectively. Subsequently, the pixels with values less than 20 were extracted, which represented the radiation-stable areas in the two datasets. Since there are several hours' difference between the overpass time of the two sensors, some of NTL on the ground is changed during this period. An intersection operation was made between the radiation-stable areas of DMSP-OLS and NPP-VIIRS to obtain the stable pixels both spatially and temporally. Finally, the pixels in the intersection constituted the vicarious site. Figure 2 shows the distribution of the vicarious site across China. gain coefficient and the overglow effect are regarded as DMSP satellite's non-linear transformation to actual ground radiation. Reconstructing the non-linear transformation can be divided into two parts: (1) constructing a robust calibration function from NPP-VIIRS to DMSP-OLS, and (2) smoothing the calibrated VIIRS to adapt to the overglow effect of DMSP-OLS. Figure 1 shows the flowchart of this study.

Preprocessing
The NTL data of both sensors were re-projected to the Albers Equal Area Conic Projection. An annual synthesis of NPP-VIIRS data was carried out since NPP-VIIRS data are monthly mean product, while DMSP-OLS data are annual mean product. Subsequently, it was resampled to 1 × 1 to match the DMSP-OLS data. In addition, there are temporary illumination and background noise in the NPP-VIIRS data; instead, only long-term and stable light exists in DMSP-OLS data. A simple method that was proposed by Li et al. [34] was adopted to eliminate the temporary illumination and background noise. After that, a small number of negative pixels were set to 0. Finally, the method that was proposed by Elvidge et al. [35] was used to improve the continuity of the DMSP-OLS data from 1992 to 2013.
A vicarious site is required in inter-calibration to provide samples [36][37][38]. However, there is no special calibration site designed for NTL in the world. This study used the coefficient of variation (CV) to select the calibration site to deal with this problem, referring to Jeswani [22]. It is assumed that the small CV stands for the stable NTL radiation. First, the CVs of 3 × 3 windows of all pixels in DMSP-OLS data and NPP-VIIRS data were computed, respectively. Subsequently, the pixels with values less than 20 were extracted, which represented the radiation-stable areas in the two datasets. Since there are several hours' difference between the overpass time of the two sensors, some of NTL on the ground is changed during this period. An intersection operation was made between the radiation-stable areas of DMSP-OLS and NPP-VIIRS to obtain the stable pixels both spatially and temporally. Finally, the pixels in the intersection constituted the vicarious site. Figure 2 shows the distribution of the vicarious site across China.

Constructing Inter-calibration Function
As the dynamic range and the luminance distribution of DMSP-OLS and NPP-VIIRS are different, it is necessary to narrow those differences before constructing a calibration function. Accordingly, the NPP-VIIRS data underwent a logarithmic transformation, which cannot only

Constructing Inter-calibration Function
As the dynamic range and the luminance distribution of DMSP-OLS and NPP-VIIRS are different, it is necessary to narrow those differences before constructing a calibration function. Accordingly, the NPP-VIIRS data underwent a logarithmic transformation, which cannot only compress its dynamic range, but also made its luminance distribution more similar to that of the DMSP-OLS data [39].
A calibration function is the key to improve the consistency between DMSP-OLS and NPP-VIIRS. First, the logarithmically transformed VIIRS and DMSP-OLS pixels within the vicarious site were extracted. Afterwards, according to the characteristics of these pixel pairs in scatter plot, Sigmoid Remote Sens. 2020, 12, 937 5 of 15 functions were used to build the relationship. A Logistic model and a BiDoseResp model were used as the candidate functions in order to find a mathematical expression that can achieve the best regression. Although both can show S-shaped growth trends, the two are different in details. The Logistic model is subject to the characteristic of strict central symmetry, which might lead to its lack of flexibility in regression. However, the BiDoseResp model is a weighted combination of two Logistic models (c1 and c2). The combination can make it free of central symmetry [40]. Therefore, the BiDoseRsep model has greater flexibility. The Logistic model is described as in Equation (1): where Y denotes the calibrated value from NPP-VIIRS to DMSP-OLS; the Bottom denotes value of Y for minimal curve asymptote; the Top denotes value of Y for maximal curve asymptote; independent variable Log[V] is the logarithm of the NPP-VIIRS value; LogMean denotes the mean of independent variable; h is the slop of the curve. The BiDoseResp model is defined as in Equation (2): where w is the weight of one of the Logistic curves; LogMean1 and LogMean2 denote the mean of independent variables of c1 and c2, respectively; h1 and h2 are the slops of c1 and c2 respectively. The coefficient of determination (R 2 ) and the residual sum of square (RSS) were taken as the evaluation indicators. After comparing the regression results, the model with better performance was used to process the NPP-VIIRS across China. In this way, the calibrated NPP-VIIRS was generated.

Mitigating Spatial Variability
There is still a gap between DMSP-OLS and the calibrated NPP-VIIRS due to the obvious overglow effect of DMSP-OLS. It is necessary to smooth the calibrated NPP-VIIRS to adapt to the overglow effect in order to maximize the consistency between the two. Some studies have adopted Gaussian low pass filter to smooth NPP-VIIRS, so that it can get similar to the DMSP-OLS [25]. However, how the different parameters (σ, w) of the filter affect the consistency has not been clearly revealed. σ was changed by iteration from 0.20 to 5.00 in steps of 0.01 and w was changed from 3 to 29 in steps of 2 to generate 6734 parameter pairs to generate a better understanding of this. Each parameter pair was taken as the input of Gaussian low pass filter to smooth the calibrated NPP-VIIRS. RSS was taken to characterize the gap between DMSP-OLS and the calibrated NPP-VIIRS. The larger the RSS, the greater the gap, and vice versa. With σ, w, and RSS as the axes of coordinates, the computing results were plotted in a 3D space to clearly reveal the effect of different parameters on the consistency. After that, the parameters corresponding to the lowest point were used as the input of the filter to smooth the calibrated NPP-VIIRS, so that it can be consistent with the DMSP-OLS data to the greatest extent. Figure 3a shows the original NPP-VIIRS and DMSP-OLS pixel pairs in the vicarious site. All of the pixel pairs are concentrated along the y-axis; therefore, it is impossible to directly make a good regression. However, after the logarithmic transformation of NPP-VIIRS, it is discovered that these pixel pairs present a significant S-shaped distribution in the scatter plot of Figure 3b. Different from the previous one, the S-shaped distribution shows the potential relationship between DMSP-OLS and NPP-VIIRS, and allows for constructing a mathematical expression to achieve the inter-calibration. Both the BiDoseResp model and the Logistic model achieve good regressions, with R 2 reaching 0.967, as shown in Table 2. However, the RSS (6.136 × 10 5 ) of the BiDoseResp model is significantly less than that (6.199 × 10 5 ) of the Logistic model. It indicates that the BiDoseResp model can bring less residuals and higher calibration accuracy. This is due to the weighted combination, which gives the BiDoseResp model more flexibility. In addition, the power and linear functions that are proposed in the previous methods were also used to perform regression. It turns out that the power function cannot converge, and the fitting of the linear function is not as good as the BiDoseResp model. Hence, when compared with other models, the BiDoseResp model can achieve the best regression. The yellow curve in Figure 3b stands for the regression curve of the BiDoseRsep model, and its parameters are shown in Table 3. Finally, the BiDoseResp-calibrated VIIRS (BDRVIIRS) was obtained by applying the model on the NPP-VIIRS data across China.  Figure 3a shows the original NPP-VIIRS and DMSP-OLS pixel pairs in the vicarious site. All of the pixel pairs are concentrated along the y-axis; therefore, it is impossible to directly make a good regression. However, after the logarithmic transformation of NPP-VIIRS, it is discovered that these pixel pairs present a significant S-shaped distribution in the scatter plot of Figure 3b. Different from the previous one, the S-shaped distribution shows the potential relationship between DMSP-OLS and NPP-VIIRS, and allows for constructing a mathematical expression to achieve the inter-calibration. Both the BiDoseResp model and the Logistic model achieve good regressions, with 2 reaching 0.967, as shown in Table 2. However, the RSS (6.136 × 10 5 ) of the BiDoseResp model is significantly less than that (6.199 × 10 5 ) of the Logistic model. It indicates that the BiDoseResp model can bring less residuals and higher calibration accuracy. This is due to the weighted combination, which gives the BiDoseResp model more flexibility. In addition, the power and linear functions that are proposed in the previous methods were also used to perform regression. It turns out that the power function cannot converge, and the fitting of the linear function is not as good as the BiDoseResp model. Hence, when compared with other models, the BiDoseResp model can achieve the best regression. The yellow curve in Figure 3b stands for the regression curve of the BiDoseRsep model, and its parameters are shown in Table 3. Finally, the BiDoseResp-calibrated VIIRS (BDRVIIRS) was obtained by applying the model on the NPP-VIIRS data across China.

Mitigating Spatial Variability
After the different Gaussian low pass filters processed BDRVIIRS, the RSS values between BDRVIIRS and DMSP-OLS were obtained. Cast into a 3D space, these 6734 RSS values construct such a surface, as shown in the Figure 4. For the first time, the surface makes it possible for us to clearly understand how the different parameters of the filter affect the consistency. In the beginning, the RSS is 9.238 × 10 7 (as indicated by Point a). After the filter is performed, the RSS begins to change. As σ increases, RSS dramatically decreases. When σ gets to 1.00, the filter becomes normally distributional at this time, which significantly improves the consistency (as indicated by the blue line). However, the current consistency has not reached its best. When σ becomes 1.51 and w becomes 15, the lowest point of the surface is eventually reached (as indicated by Point b). The minimum RSS is 2.932 × 10 7 and the corresponding root mean squared error (RMSE) is 3.976. As parameters continue to change, RSS turns to increase, and the consistency starts to decrease. This is because BDRVIIRS is over-smoothed at this time. Finally, the parameters corresponding to the lowest Point b was selected as the input of the filter to process BDRVIIRS. In this way, it was smoothed to be consistent with DMSP-OLS to the greatest extent. After the different Gaussian low pass filters processed BDRVIIRS, the RSS values between BDRVIIRS and DMSP-OLS were obtained. Cast into a 3D space, these 6734 RSS values construct such a surface, as shown in the Figure 4. For the first time, the surface makes it possible for us to clearly understand how the different parameters of the filter affect the consistency. In the beginning, the RSS is 9.238 × 10 7 (as indicated by Point a). After the filter is performed, the RSS begins to change. As σ increases, RSS dramatically decreases. When σ gets to 1.00, the filter becomes normally distributional at this time, which significantly improves the consistency (as indicated by the blue line). However, the current consistency has not reached its best. When σ becomes 1.51 and w becomes 15, the lowest point of the surface is eventually reached (as indicated by Point b). The minimum RSS is 2.932 × 10 7 and the corresponding root mean squared error (RMSE) is 3.976. As parameters continue to change, RSS turns to increase, and the consistency starts to decrease. This is because BDRVIIRS is over-smoothed at this time. Finally, the parameters corresponding to the lowest Point b was selected as the input of the filter to process BDRVIIRS. In this way, it was smoothed to be consistent with DMSP-OLS to the greatest extent.   After the different Gaussian low pass filters processed BDRVIIRS, the RSS values between BDRVIIRS and DMSP-OLS were obtained. Cast into a 3D space, these 6734 RSS values construct such a surface, as shown in the Figure 4. For the first time, the surface makes it possible for us to clearly understand how the different parameters of the filter affect the consistency. In the beginning, the RSS is 9.238 × 10 7 (as indicated by Point a). After the filter is performed, the RSS begins to change. As σ increases, RSS dramatically decreases. When σ gets to 1.00, the filter becomes normally distributional at this time, which significantly improves the consistency (as indicated by the blue line). However, the current consistency has not reached its best. When σ becomes 1.51 and w becomes 15, the lowest point of the surface is eventually reached (as indicated by Point b). The minimum RSS is 2.932 × 10 7 and the corresponding root mean squared error (RMSE) is 3.976. As parameters continue to change, RSS turns to increase, and the consistency starts to decrease. This is because BDRVIIRS is over-smoothed at this time. Finally, the parameters corresponding to the lowest Point b was selected as the input of the filter to process BDRVIIRS. In this way, it was smoothed to be consistent with DMSP-OLS to the greatest extent.

Trend of Calibrated NTL of China
Due to the lack of on-board calibration for DMSP-OLS sensors, there is the lack of continuity between the original DMSP-OLS data from 1992 to 2013, so the total light value (TLV) between different years and different satellites fluctuates severely, as shown in Figure 6a. A quadratic polynomial method proposed by Elvidge et al. [35] was used to make an inter-calibration between the DMSP-OLS data of different years and different satellites to solve these problems. As a result, the continuity of DMSP-OLS data was enhanced after inter-calibration, and the fluctuation of TLV between different years and different satellites was significantly alleviated. On the other hand, since the non-calibrated NPP-VIIRS data is inconsistent with DMSP-OLS data, there is a big gap between the TLV of the two kinds of sensors. Therefore, the annual NPP-VIIRS data from 2012 to 2018 was calibrated by the method that was proposed in this paper. Eventually, as presented in Figure 6b, the TLV of the DMSP-OLS data and that of the BDRVIIRS data show a good continuity. The TLV of China has been increasing from 1992 to 2018, which is consistent with the economic development of China during this period.  Figure 4. Figure 5a shows the scatter density plot of BDRVIIRS and DMSP-OLS of the whole China region without a Gaussian low pass filter, which corresponds to Point a in Figure 4. The scatter distribution is very discrete. This is due to the conflict between the overglow effect of DMSP-OLS and the high-precision of unsmoothed BDRVIIRS. The Gaussian low pass filter (σ = 1.51, w = 15) can make BDRVIIRS adapt to DMSP-OLS as much as possible. The smoothed BDRVIIRS reduces the scatter dispersion significantly, which corresponds to Point b in Figure 4, as Figure 5b shows. At this time, it shows a significant positive correlation with DMSP-OLS where the correlation coefficient (r) reaches 0.949. Due to the lack of on-board calibration for DMSP-OLS sensors, there is the lack of continuity between the original DMSP-OLS data from 1992 to 2013, so the total light value (TLV) between different years and different satellites fluctuates severely, as shown in Figure 6a. A quadratic polynomial method proposed by Elvidge et al. [35] was used to make an inter-calibration between the DMSP-OLS data of different years and different satellites to solve these problems. As a result, the continuity of DMSP-OLS data was enhanced after inter-calibration, and the fluctuation of TLV between different years and different satellites was significantly alleviated. On the other hand, since the non-calibrated NPP-VIIRS data is inconsistent with DMSP-OLS data, there is a big gap between the TLV of the two kinds of sensors. Therefore, the annual NPP-VIIRS data from 2012 to 2018 was calibrated by the method that was proposed in this paper. Eventually, as presented in Figure 6b, the TLV of the DMSP-OLS data and that of the BDRVIIRS data show a good continuity. The TLV of China has been increasing from 1992 to 2018, which is consistent with the economic development of China during this period.

Evaluating Image Texture
After the nationwide calibration, Beijing, Shanghai, and Hangzhou were selected to conduct a detailed evaluation. Figure 8 presents images of DMSP-OLS, NPP-VIIRS, and BDRVIIRS in the three cities. DMSP-OLS images show urban centers, intercity roads, and suburbs simultaneously. However, the original NPP-VIIRS images only show urban centers, different from DMSP-OLS. Through the inter-calibration, BDRVIIRS images can show urban centers, intercity roads and suburbs simultaneously just like DMSP-OLS images. Therefore, it is concluded that the inter-calibration makes image textures of BDRVIIRS data more similar to that of DMSP-OLS data. After the nationwide calibration, Beijing, Shanghai, and Hangzhou were selected to conduct a detailed evaluation. Figure 8 presents images of DMSP-OLS, NPP-VIIRS, and BDRVIIRS in the three cities. DMSP-OLS images show urban centers, intercity roads, and suburbs simultaneously. However, the original NPP-VIIRS images only show urban centers, different from DMSP-OLS. Through the inter-calibration, BDRVIIRS images can show urban centers, intercity roads and suburbs simultaneously just like DMSP-OLS images. Therefore, it is concluded that the inter-calibration makes image textures of BDRVIIRS data more similar to that of DMSP-OLS data.  a, b, and c); NPP-VIIRS images in the three cities (d, e, and f); BDRVIIRS images of the three cities (g, h, and i). Figure 9 shows the transects of NTL data of the three cities. It is obvious that all of the transects in the areas far away from urban centers are low and elevated toward urban centers. When the transects pass through water body, the brightness suddenly decreases. By comparing NPP-VIIRS transects with DMSP-OLS transects, it can be found that there is a significant difference between the two, although both are based on the same spatial resolution. There are many dramatic variations in NPP-VIIRS transects, parts of which are towering peaks of brightness, while some other parts are narrow valleys. Furthermore, the variations from peaks to valleys are almost vertical. However, in DMSP-OLS data, the transects are much smoother due to the overglow effect. Similarly, the BDRVIIRS transects are also smoother than the original NPP-VIIRS transects due to the Gaussian low pass filter. Overall, BDRVIIRS transects and the DMSP-OLS transects are similar to each other.  Figure 9 shows the transects of NTL data of the three cities. It is obvious that all of the transects in the areas far away from urban centers are low and elevated toward urban centers. When the transects pass through water body, the brightness suddenly decreases. By comparing NPP-VIIRS transects with DMSP-OLS transects, it can be found that there is a significant difference between the two, although both are based on the same spatial resolution. There are many dramatic variations in NPP-VIIRS transects, parts of which are towering peaks of brightness, while some other parts are narrow valleys. Furthermore, the variations from peaks to valleys are almost vertical. However, in DMSP-OLS data, the transects are much smoother due to the overglow effect. Similarly, the BDRVIIRS transects are also smoother than the original NPP-VIIRS transects due to the Gaussian low pass filter. Overall, BDRVIIRS transects and the DMSP-OLS transects are similar to each other.

Evaluating Transect
Besides, there are some other details worth noting. The changes in human activities during this period lead to the variations in brightness in some areas because of the difference in overpass time. Therefore, the transects of DMSP-OLS and BDRVIIRS are not exactly coincident. In addition, a saturation effect exists in DMSP-OLS data, while there is no saturation effect in the BDRVIIRS data. Accordingly, BDRVIIRS transects show more details than DMSP-OLS transects in the urban centers. Meanwhile, this also indicates that the underestimation of BDRVIIRS exists in areas where DMSP-OLS data are saturated. Besides, there are some other details worth noting. The changes in human activities during this period lead to the variations in brightness in some areas because of the difference in overpass time. Therefore, the transects of DMSP-OLS and BDRVIIRS are not exactly coincident. In addition, a saturation effect exists in DMSP-OLS data, while there is no saturation effect in the BDRVIIRS data. Accordingly, BDRVIIRS transects show more details than DMSP-OLS transects in the urban centers. Meanwhile, this also indicates that the underestimation of BDRVIIRS exists in areas where DMSP-OLS data are saturated. Figure 10 shows the scatter plots and the histograms of DMSP-OLS and BDRVIIRS of the three cities. These scatter plots indicate a strong correlation between DMSP-OLS and BDRVIIRS. The regression analysis illustrates that the correlation coefficients of these scatter plots are greater than 0.97, RMSE are less than 4.30, and the slopes of the regression lines are close to 1. In addition, the histograms of BDRVIIRS and DMSP-OLS are consistent. This means that the distribution of the brightness of BDRVIIRS and DMSP-OLS in the dynamic range is consistent. All of these indicate that the consistency between BDRVIIRS and DMSP-OLS is strong.  Figure 10 shows the scatter plots and the histograms of DMSP-OLS and BDRVIIRS of the three cities. These scatter plots indicate a strong correlation between DMSP-OLS and BDRVIIRS. The regression analysis illustrates that the correlation coefficients of these scatter plots are greater than 0.97, RMSE are less than 4.30, and the slopes of the regression lines are close to 1. In addition, the histograms of BDRVIIRS and DMSP-OLS are consistent. This means that the distribution of the brightness of BDRVIIRS and DMSP-OLS in the dynamic range is consistent. All of these indicate that the consistency between BDRVIIRS and DMSP-OLS is strong.

Comparison with Related Works
Several studies have attempted to make the inter-calibration between DMSP-OLS and NPP-VIIRS. Table 4 shows the comparison with previous works. The differences of these studies are mainly reflected in the data characteristics and calibration functions. In terms of data characteristics, early researches used daily and monthly data [21,25], but those methods cannot be applied on annual data, so the subsequent studies focused on annual data. In terms of calibration function, according to available information, the BiDoseResp function in this study performs best, whose 2 is 0.967 and RMSE is 1.75. In addition, the method that is proposed in this study has significant advantages in applicable scope, accuracy, and efficiency.

Comparison with Related Works
Several studies have attempted to make the inter-calibration between DMSP-OLS and NPP-VIIRS. Table 4 shows the comparison with previous works. The differences of these studies are mainly reflected in the data characteristics and calibration functions. In terms of data characteristics, early researches used daily and monthly data [21,25], but those methods cannot be applied on annual data, so the subsequent studies focused on annual data. In terms of calibration function, according to available information, the BiDoseResp function in this study performs best, whose R 2 is 0.967 and RMSE is 1.75. In addition, the method that is proposed in this study has significant advantages in applicable scope, accuracy, and efficiency.
One of the major advantages is applicable scope, including spatial area, dynamic range, and time range. First, in spatial area, most of the studies are based on administrative units as research areas. Li et al. [25] and Jeswani [22] selected Syria and India as research areas, respectively, while recent studies narrow the research areas and focus on cities for better performance. However, the research area of this study is larger than that of previous researches. Accordingly, this method is more suitable for promotion in larger areas. Secondly, the dynamic range of these studies is incomplete [21,25,26]. In contrast, this study's dynamic range is 0-63. Third, previous studies have not fully exploited the potential of historical NTL archiving, as described in the recent review [20]. Instead, our method, combined with Elvidge's one [35], can maximize the time range of the annual NTL sequence. Finally, the annual NTL sequence of China was produced from 1992 to 2018. Hence, this study is excellent in spatial area, dynamic range, and time range.
Another advantage lies in accuracy. To make a comparison with Wu et al. [26] and Zheng et al. [27], the calibrated NTL in the same areas was extracted. With the help of available information, it is found that the r and RMSE in this study are both superior to the studies of Li et al. [25] and Wu et al. [26]. In addition, r values of Zheng et al. [27] in Hangzhou and Beijing (0.993 and 0.995) are slightly higher than the values (0.983 and 0.975) of this study in the corresponding regions. While their r value in Shanghai (0.966) is lower than 0.986 of this study. Therefore, the accuracy of this study is superior to Li et al. [25] and Wu et al. [26], and close to Zheng et al. [27].
This study is better not only in applicable scope and accuracy, but also in efficiency. First, the method that was proposed by Zheng et al. [27] requires input data to undergo a complex preprocessing, while the preprocessing of this study is much simpler. Secondly, the functions used in Zheng's [27] and Wu's [26] methods are just suitable for city-level applications. What is more, these functions need to be reconstructed when they are shifted to other areas. On the contrary, this method can achieve one-time calibration at the country level or even global, and then, the comparison between multiple cities can be carried out directly. Therefore, the method of this study is simpler and more efficient. When compared with the existing studies, this one has produced the calibrated NTL sequence with the vastest research area, longest time range, and complete dynamic range. In addition, the evaluation indicators show that the calibration has high accuracy. Therefore, this method is efficient and worth promoting. The products can be used as extended DMSP-OLS data.

Conclusions
In this study, a new inter-calibration method was designed to improve the consistency between NPP-VIIRS and DMSP-OLS. The core of this method is the BiDoseResp model. Its coefficient of determination reaches 0.967, and residual sum of squares is 6.136 × 10 5 , less than 6.199 × 10 5 of Logistic function. Subsequently, the Gaussian low pass filter with the optimal parameters promoted the consistency further. After that, the correlation coefficient between DMSP-OLS and BDRVIIRS of the whole China region reaches 0.949. Taking the results of Beijing, Hangzhou, and Shanghai for a detailed analysis, it is found that the results perform well in image feature, transect, and pixel relationship. Finally, the consistent NTL sequence of China from 1992 to 2018 is generated. The result shows that the NTL in China have increased continuously and significantly during this period.
Although the method that is proposed in this study is better in many ways, there are also some limitations. First, the NPP-VIIRS data were aggregated temporally and spatially according to the DMSP-OLS data, which caused the NPP-VIIRS data to lose their finer temporal and spatial resolution. Secondly, the overglow effect of NTL data makes the areas around the light sources also recorded as bright pixels, which makes it larger than the actual light source area. However, as compared with the NTL data of single sensor, the integrated NTL sequence can promote the mining of richer and deeper information about human settlement.
It is expected that this method can facilitate the development of researches on NTL and that more novel information can be mined out. For example, the BDRVIIRS data would not be completely consistent with the DMSP-OLS data due to the difference in overpass time. The inconsistency indicates the changes in NTL on the ground (enhancement or dimming). Therefore, the way forward could be analyzing the changes in NTL during this period to reveal human activities further. In addition, future directions also include the development of a calibration method from DMSP-OLS to NPP-VIIRS, which is expected to generate a consistent NTL dataset that is characterized by higher quality.