Ensemble data assimilation (EDA) has shown great promise in perfect model studies. The flow-dependent background error covariances which EDA provides have been shown to be quite benefical, especially when observations are sparse and/or infrequent. When applied to real forecast systems, with imperfect models, results to date have been mixed. The need to parameterize model error not accounted for in the ensemble itself is thought to be the primary factor limiting the performance of first-generation EDA systems. We have implemented a prototype EDA system for the NCEP Global Forecast System (GFS). Analyses obtained with a 100 member T62L28 ensemble, using real observations for January and February 2004, will be shown and compared with analyses obtained with the NCEP operational analysis system. Sensitivity to ensemble size and the scheme used for parameterizing model error will be discussed, along with prospects for operational implementation.