At the previous NWP conference, we reported on an incremental 3DVAR data assimilation system developed for small-scale nonhydrostatic models, in particular, the ARPS model. Within the system, a cost function is defined as a sum of background and observation constraints. The background field is typically provided by previous model forecasts. For special applications, e.g., for local-scale Doppler radar data retrievals, the background can also be provided by a single environmental sounding. In non-cycled mode, as is often unavoidable in research settings, the background can be from forecasts of an external model. In the system, the background covariances are modeled using recursive spatial filters. The system currently handles all traditional observations as well as raw Doppler velocity measurements. In this paper, we will report on the further development of this 3DVAR system, in particular, the addition of (weak) equation constraints to the cost function. These equations are based on nonhydrostatic equations of motion and the three dimensional mass continuity equation. These equations are important in providing the coupling among the velocity components and with thermodynamic fields. The paper will present the 3DVAR formulation and the results of numerical experiments.