This manuscript aims to reconstruct global satellite Chlorophyll a (Chl) fields over the last 25 years, combining satellite observations with BGC-Argo float data and spatial kriging techniques. Gap-free Chl fields, including polar regions during wintertime, are vital for applications that require complete inputs, such as ocean CO2 sink assessment. This paper is well-motivated and enjoyable to read, especially the discussion section, which was strong and well-argued. The associated dataset on Zenodo is very well-referenced and is likely to be useful for the community.

That said, I have several concerns and suggestions, notably on some of the methods employed that I believe should be addressed before publication.

General comments :

My main comment about this manuscript is that it foregoes an important limitation of satellite-based Chl reconstructions: photoacclimation under clouds. By computing monthly averages from only clear-sky days, one neglects the increased intracellular Chlorophyll content of phytoplankton under reduced light under clouds. This increase can achieve a tenfold factor. Consequently, reconstructions based solely on clear-sky spatial patterns (whether kriging, DINEOF, or others) underestimate Chl in regions with significant cloud coverage. While I recognize this limitation may be beyond the scope of the present study, I strongly suggest including a discussion of photoacclimation as a key caveat for Chl gap-filling.

I would also moderate some of your conclusions about the Wintertime gap-filling. Although there is a relatively good fit between the Wintertime gap filled and the Valente HPLC data, all other comparisons, by any metrics employed, are not what I would qualify as “good”. I believe the results against fluorometric data warrant a discussion on sources of disagreement and what further steps can be used to reduce this gap. Integrating over the whole penetration depth rather than 20meters and changing the technique used for the fluorometric correction would likely improve the performance of the reconstructions, see details below.

From personal preference, throughout this manuscript (and in the figures), I would consider using linear scale units for the reader to relate to the Chl values present. In the figures, I would suggest putting the Chl values in mg m-3 and using a logarithmic axis for increased readability.  

Specific comments:

BGC-Argo usage: Although BGC-Argo floats are a formidable tool to validate and complement remote sensing data, some slight methodological changes could result in an improved product.

• BGC-Argo data estimate Chl from fluorescence, which is not directly relatable to Chl from Satellite. In Roesler (2017), which you cite, they identify an average factor of 2 difference, which has a large variability across regions (up to a factor of 6 in the Southern Ocean !). In Section 2.4 you mention that you do account for this, but I don’t believe that the method employed is accurate. The Roesler paper does not suggest using a single value bias but rather using a “Slope factor”, which is much more accurate than a single bias value (that affects very differently small and large Chl concentrations). Accurately applying a Slope factor should significantly improve the relationship between BGC-Argo and OC-CCI (See Xing (2011) on a method to compute it from float radiometry). It is unclear to me if Figure S1 is prior or after the fluorescence-bias correction, but I would expect it to be much closer to a 1:1 line after a slope factor correction and reduce your “high intercept values” that you mention in paragraph 115.
• You explain computing the mean Chl value from the first 20meters of each profile. This likely underestimates your Chl compared to what the Satellite measures, as increased chlorophyll concentration (such as in the DCM) can be found at depths deeper than 20meters, yet within the layer visible from the satellite. When comparing Satellite data with in-situ profiles, a commonly accepted technique is, for a given profile, to integrate/average over the penetration depth (Zpd) as this is considered a good approximation of what the satellite sees. It is computed as Zpd = 1/Kd(490). I would suggest retrieving Kd(490) either from the Satellite pixel or, for more accuracy, to compute it from a BGC-Argo float Ed(490) profile, see Xing (2020). 
• Some of the information on the correction that is in 2.4 would probably be more appropriate in Section 2.2, so the reader knows at once how the BGC-Argo data were processed.

On the Spatial Kriging:

• Please quantify the fraction of ocean pixels filled by kriging versus BGC-Argo. Supplementary Fig. S3 shows temporal coverage, but a spatial map distinguishing contributions of each method would be more informative. Additionally, it is important to mention that this Kriging method is effectively not filling specifically cloudy values but rather any pixels that have been permanently obscured for a given month. This should be emphasized, as persistent coccolithophore blooms have also resulted in pixels being flagged. I understand that most empty pixels are caused by clouds, but the technique here is not specific to clouds.  
• There are also numerous BGC-Argo floats and in-situ datapoints from the Valente dataset in the area filled by the spatial Kriging. Although the paper’s main point is not on this already published method, it would be strengthened by the evaluation of the performance of the Spatial Kriging.

On the wintertime reconstruction:

In general, I thought the wintertime BGC-Argo reconstruction method could benefit from additional details, some reorganizing, and perhaps a schematic? I had to reread this section several times, and I am still convinced I have not understood this section correctly.

•  The first BGC-Argo floats used in this study were deployed in 2008, and yet Figure S2 and your text mention that the maximum time lag between OC-CCI and BGC-Argo used to fill a gap was 9 months. I am therefore unclear on how the Wintertime reconstruction was performed for those 9 years before the first Argo profile? From looking at Figure S2, I am hypothesizing that you used the BGC-Argo float to create some kind of monthly climatology in pixels, but I was unable to find explicit mention of this in Section 2.4.  You mention “For each time lag, the median percentage difference was calculated between the OC-CCI and bias-corrected BGC-Argo chl-a in mg m-3 […] on a pixel-by-pixel basis”. I am again assuming you mean across all 25 years, within a given pixel, you find all BGC-Argo profiles that occurred within a month of the last OC-CCI measurement and compute the median, before repeating this for the 2-month lag and so on, but I believe the reader would really benefit from a rewriting of this section.
•  I would put emphasis in the discussion that the wintertime reconstruction is based on data from 2008 on and acknowledge the fact that this reconstruction is weighted around the time period in which there are more floats. This information is presented in Figure 1, but the limitations associated with this technique and uneven sampling frequency is only quickly mentioned in paragraph 365, and would benefit from a more thorough discussion, linking it to areas that have experienced rapid change in productivity and ice coverage over the last 10 years, notably in the Arctic.

In Figure 2, the boxplots show a very large spread in the percent difference. It would be interesting to see if there are spatial patterns around this spread, notably if some areas of the Northern Hemisphere decrease in Chl more rapidly than others.

For the in-situ Valente data of Chl fluorescence, has a conversion been applied similarly to the fluorescence by Argo (the Slope factor from Roesler, (2017)) ?  This might significantly help improve the comparison with Chl reconstructed from BGC-Argo.

Figure 3(b)-(f) It is impossible to distinguish between the two timeseries on the plots due to space constrain (the x-axis is squished). I would suggest either removing 1or 2 timeseries graph or splitting them over two rows, as currently I cannot draw any conclusion from those.

Figure 3 and 4: I would put the time series values back in linear scale unit.

Roesler, C., Uitz, J., Claustre, H., Boss, E., Xing, X., Organelli, E., Briggs, N., Bricaud, A., Schmechtig, C., Poteau, A., D’Ortenzio, F., Ras, J., Drapeau, S., Haëntjens, N., & Barbieux, M. (2017). Recommendations for obtaining unbiased chlorophyll estimates from in situ chlorophyll fluorometers: A global analysis of WET Labs ECO sensors. Limnology and Oceanography: Methods, 15(6), 572–585. https://doi.org/10.1002/lom3.10185

Xing, X., Morel, A., Claustre, H., Antoine, D., D’Ortenzio, F., Poteau, A., & Mignot, A. (2011). Combined processing and mutual interpretation of radiometry and fluorimetry from autonomous profiling Bio-Argo floats: Chlorophyll a retrieval. Journal of Geophysical Research, 116(C6), C06020. https://doi.org/10.1029/2010JC006899

Xing, X., Boss, E., Zhang, J., & Chai, F. (2020). Evaluation of Ocean Color Remote Sensing Algorithms for Diffuse Attenuation Coefficients and Optical Depths with Data Collected on BGC-Argo Floats. Remote Sensing, 12(15), 2367. https://doi.org/10.3390/rs12152367

This paper proposes a procedure to generate a gap-free time series of global monthly chlorophyll maps based on satellite data. The novelty is that this product covers also the polar oceans thanks to the availability of bio-argo profiles during the polar night. This effort is motivated by the common need in climate studies for spatially complete datasets spanning the longest possible time interval. The time series begins in 1998 with the onset of ocean colour missions and is openly available on Zenodo, representing a valuable resource for numerous climate research applications.

Considering the usefulness and the importance of the proposed data set I am convinced that the paper merits to be published but, at the same time, I also have some major concern about the reconstruction methods applied, that need to be better explained and qualified.  In addition, I am not convicted on the claimed better performance of the Kriging method respect other methods commonly used in literature or routinely applied by operational centres such as the CMEMS Ocean Colour Thematic Assembly Center (surprisingly is never cited in the manuscript).

The method used to reconstruct the chlorophyll field during the polar night, while rather convoluted and extremely crude, it is better or, at least, less wrong than use a constant value representing a small step forward in the production of global gap-free satellite images.

Specific Comments:

2.2 Satellite observational data:

1 - Daily maps used to produce monthly means are the result of a composition of several passage acquired by several satellite missions. Can you report here how the daily aggregated maps are produced or, at least include a citation where how the daily composite maps are produced is described?

2 – line 100: “(1 sigma; given as the root mean square difference; RMSD)”. RMSD of what?

3 – line 100-102. Also assuming that “spatial uncertainties within adjacent cells are dependent and spatially correlated” I am not totally sure that the mean of the 4 km uncertainties is the mean of the single standard deviation. Can you please indicate where in Taylor 1997 this specific point is discussed and proved?

4 – line 113-115: The authors use a type II regression: are errors on BGC-Argo and OC-CCI data that enter in the regression comparable?

Section 2.3: Spatial Kriging for cloud gap filling

1 – line 130: “…. used to estimate the chl-a concentration”,…  add citation.

2 - The choice to reconstruct the field over data voids using Ordinary Kriging is primarily based on the work of Stock et al. (2020) which compares Ordinary Kriging, DINEOF, and several widely used AI methods. Optimal interpolation is not considered and other advanced methods based on Singular Spectra Analysis (Kondrashov, D., & Ghil, M., 2006 Spatio-temporal filling of missing points in geophysical data sets. Nonlinear Processes in Geophysics, 13(2), 151-159) are not even mentioned.  I understand that the Ordinary Kriging method is significantly less computationally demanding compared to some of more sophisticated methods (e.g. SSA or DINEOF or Optimal Interpolation); however, if this is the case, it should be clearly stated in the text, rather than simply claiming that Kriging and DINEOF perform better.

it should be discussed and what is the advantage of using ordinary kriging respect to other kriging methods such as the method adopted by CMEMS L4 Global chlorophyll product (Saulquin et al, 2018).

In addition, it will be important to compare the proposed product with other monthly L4 global chlorophyll products such as those distributed by CMEMS (OCEANCOLOUR_GLO_BGC_L4_NRT_009_102, OCEANCOLOUR_GLO_BGC_L4_MY_009_104) and discuss the results.  

3 -  How is polar night defined? Below what solar elevation value is it defined as “polar night”?

4 – “…. The semi-variogram was fit to a ~5% subset of the OC-CCI observations that were equally distanced in space, for a monthly varying latitude band where at least 20% of the OC-CCI observations are available….”  Does this mean that the parameters of the exponential function are calculated for each latitude band and, consequently, the fitted function depends on latitude? Please clarify.

5 – Finally, considering that the cloud gap filling kriging approach uses observations in the vicinity of the empty pixel (see first line of the Uncertainty propagation section) is your kriging including a definition of an influential distance that limits the search radius?

Section 2.4: BGC-Argo Wintertime filling

1 –  In section 2.2, the authors correctly note that substantial differences can occur between OC-CCI (based on empirical algorithms that use HPLC data) and BGC-Argo chlorophyll measurements. However, in section 2.4, they address this issue by simply applying two constant bias corrections to the BGC-Argo data, one for each hemisphere, justifying this approach using the results shown in Figure S1.

Since in Figure S1 it is evident that the difference between the two dataset is not limited to bias, the question that arises here is: why not also account for the slope of the relationship, which would likely allow for a more accurate correction?

2 – Figure 1b shows that before 2010 BGC-Argo profiles are not available. If I have correctly understood, the filling procedure adopted by authors requires to have satellite data in the next spring and/or in the previous autumn and BGC-Argo profiles in between. In the absence of BGC-Argo profiles, it is unclear how the filling procedure can be applied and how the data gaps are filled. Please clarify the methodology adopted.

3 – line 179-181: “…..For each wintertime pixel, the time lag between the autumntime and springtime OC-CCI observations was calculated, and the relationship with the lowest time lag was used”.

 Since both values have been calculated, what prevents the use of a weighted average of the two chlorophyll values, in analogy with the approach used in the case of Kriging?

4 – The authors wrote: “…..Any remaining pixels that were not gap filled by any of the previous procedures are filled with a final kriging pass……”. In this regard, it would be useful to quantify the percentage of sea pixels that remain unfilled after applying the wintertime filling procedure.

Section 2.5: Uncertainty propagation

1 – lines 215-219: The description of the method used to estimate uncertainty in percentage differences is somewhat convoluted and difficult to follow. How is absolute deviation converted into an equivalent standard deviation? How are the two sources of uncertainty propagated through the analysis? And how is the Monte Carlo approach applied?

Section 3. Results

1 – Figure 3 shows a gap-free global map in which chlorophyll data are available at all latitude, Figure 1b shows that BGC-Argo profiles are available only below 75° N. It will be important to highlight in the figure caption that area covered by 90% of ice where set to 0.1 mg/m3.

Section 4. Discussion

1 - In the discussion line 317-319 it is recognized that produce higher temporal resolution datasets is possible but the results of the reconstruction over areas of large and persistent cloud cover could be questionable. While this is certainly correct, it raises the question of how much of the reconstruction difficulty is due to the use of purely spatial interpolation, as opposed to spatio-temporal interpolation employed in other approaches, such as Optimal Interpolation.

In addition, in this section will be important to compare the proposed product with other L4 available chlorophyll products used by the user community and discuss the difference and similarity.

References:

Saulquin, B., Gohin, F., and Fanton d’Andon, O.: Interpolated fields of satellite-derived multialgorithm chlorophyll-a estimates at global and European scales in the framework of the European Copernicus-Marine Environment Monitoring Service, J. Oper. Oceanogr., 12, 47–57, 2018, www.tandfonline.com/doi/full/10.1080/1755876X.2018.1552358.