Multispectral (MS) images with high spatial resolution (HR) can be obtained by fusing MS images and panchromatic (PAN) image, the HR MS images have an important significance in image interpretation and classification, etc. In the present paper, a new image fusion method based on gradient consistency constraint for MS/PAN images is developed. The method is based on maximum a posteriori (MAP) framework. The relationship of desired HR MS images and PAN image is formulated by gradient consistency constraint. Observation model of MS images and the Huber-Markov priori are combined to solve the fused image by gradient descent algorithm. In the proposed method, gradient consistency constraint is introduced, and defect of band number restriction is overcomed in conventional model-based fusion methods. Iterative step for every band is solved adaptively, and spectral characteristics of each band are fully taken into account, so it not only ensures the spectral information fidelity, but also improves the integration degree of spatial information of fused image. The proposed method has been tested using IKONOS and WorldView-2 images. It is compared with GS, AIHS and AMBF fusion methods from both qualitative and quantitative aspects. Experimental results show that the proposed method can better preserve spectral information while enhance spatial resolution, and it has broader applicability and better fusion result than other methods.