Assimilation of river discharge in a land surface model to improve estimates of the continental water cycles

River discharge plays an important role in earth's water cycle, but it is difficult to estimate due to un-gauged rivers, human activities and measurement errors. One approach is based on the observed flux and a simple annual water balance model (ignoring human processes) for un-gauged rivers, but it only provides annual mean values which is insufficient for oceanic modelings. Another way is by forcing a land surface model (LSM) with atmospheric conditions. It provides daily values but with uncertainties associated with the models. We use data assimilation techniques by merging the modeled river discharges by the ORCHIDEE (without human processes currently) LSM and the observations from the Global Runoff Data Centre (GRDC) to obtain optimized discharges over the entire basin. The model systematic errors and human impacts (dam operation, irrigation, etc.) are taken into account by an optimization parameter (with annual variation), which is applied to correct model intermediate variable runoff and drainage over each sub-watershed. The method is illustrated over the Iberian Peninsula with 27 GRDC stations over the period 1979-1989. ORCHIDEE represents a realistic discharge over the north of the Iberian Peninsula with small model systematic errors, while the model overestimates discharges by 30-150% over the south and northeast regions where the blue water footprint is large. The normalized bias has been significantly reduced to less than 30% after assimilation, and the assimilation result is not sensitive to assimilation strategies. This method also corrects the discharge bias for the basins without observations assimilated by extrapolating the correction from adjacent basins. The correction increases the interannual variability in river discharge because of the fluctuation of water usage. The E (P - E) of GLEAM (Global Land Evaporation Amsterdam Model, v3.1a) is lower (higher) than the bias-corrected value, which could be due to the different P forcing and probably the missing processes in the GLEAM model.