A method is presented for estimating the error covariance of the errors in the model equations in observation space. Estimating model errors in this systematic way opens up the possibility to use data assimilation for systematic model improvement at the level of the model equations, which would be a huge step forward. This model error covariance is perhaps the hardest covariance matrix to estimate. It represents how the missing physics and errors in parametrizations manifest themselves at the scales the model can resolve. A new element is that we use an efficient particle filter to avoid the need to estimate the error covariance of the state as well, which most other data assimilation methods do require. Starting from a reasonable first estimate, the method generates new estimates iteratively during the data assimilation run, and the method is shown to converge to the correct model error matrix. We also investigate the influence of the accuracy of the observation error covariance on the estimation of the model error covariance and show that, when the observation errors are known, the model error covariance can be estimated well, but, as expected and perhaps unavoidably, the diagonal elements are estimated too low when the observation errors are estimated too high, and vice versa.