Remote retrieval of Normalized Difference Vegetation Index (NDVI) over the Earth’s surface is a critical component of monitoring the surface processes of our planet. NDVI is a widely used and useful indicator of vegetation health and quantity ­ however its retrieval using satellite data is hindered by the frequent presence of clouds in the Earth’s atmosphere. Zeng et al. (2016) developed a novel technique that estimates a surface's Bidirectional Reflectance Distribution Function (BRDF) with a Ross­Li­Maignan (RLM) BRDF model from a set of observations. This method, the Changing­Weight Iterative (CWI) method, uses iterative a posteriori estimation of observation errors to reduce the impact of cloud-contaminated measurements in the sample. Its performance was compared to two conventional methods, ordinary­least squares (OLS) and Li­Gao BRDF-fitting. The three different BRDF­fitting methods were compared in a numerical experiment. 6,000 surface types covering a broad range of surface types were modeled using the canopy radiative transfer model PROSAIL. For each surface, sets of pseudo-observations of the surface’s red and NIR band reflectance were generated using realistic sun­target view geometries from the MODIS and MERSI satellite sensors. The effects of cloud­contamination were simulated by adding different numbers of cloud­contaminated observation to the sample, with varying degrees of contamination. The RLM BRDF model was fitted to these samples using the three different methods to estimate the BRDF model parameters. These were subsequently used to calculate a NDVI composite value. Each method’s estimate was compared to a reference ­value generated by PROSAIL. Results for the 6,000 surfaces confirmed that the CWI method is more noise­resistant than OLS and Li­Gao in situations with many observations (i.e. a large sample), and resulted in estimates that more closely matched the reference value from PROSAIL, compared to the conventional Li­Gao and OLS methods. In scenarios of low­cloud contamination, all three methods failed to detect and significantly suppress the impact of noisy observations, which was expected from existing literature. For a large­sized sample of 13 pseudo­observations studied for the validation site Mongu, Zambia, the CWI method was observed to have a very accurate performance, for up to 5 contaminated observations in the sample. With smaller sized samples of 8 and 10 for two other validation sites, it was found that the RMSE of the CWI method would suddenly increase approximately tenfold when the number of contaminated observations increased beyond 2 and 3, respectively. After these ’tipping points’, the Li­Gao method was more accurate and outperformed CWI. The CWI method therefore performed promisingly when given a large enough sample size, and in these cases it was more accurate than the conventional Li-Gao and OLS methods. However, when it fails to correctly identify noisy observations, its accuracy could decrease suddenly, which should be taken into consideration for operational use. Since the results of the experiment were averaged over 6,000 different sampling points of the PROSAIL model's parameter space, it is suggested that the conclusions apply to a wide range of surface types found all over the Earth.