We report an improved technique for diffuse foreground minimization from Cosmic Microwave Background (CMB) maps using a new multi-phase iterative internal-linear-combination (ILC) approach in harmonic space. The new procedure consists of two phases. In phase 1, a diffuse foreground cleaned map is obtained by performing a usual ILC operation in the harmonic space in a single iteration over the desired portion of the sky. In phase 2, we obtain the final foreground cleaned map using an iterative ILC approach also in the harmonic space, however, now, during each iteration of foreground minimization, some of the regions of the sky that are not being cleaned in the current iteration, are replaced by the corresponding cleaned portions of the phase 1 cleaned map. The new ILC method nullifies a foreground leakage signal that is otherwise inevitably present in the old and usual harmonic space iterative ILC method. The new method is flexible to handle input frequency maps, irrespective of whether or not they initially have the same instrumental and pixel resolution, by bringing them to a common and maximum possible beam and pixel resolution at the beginning of the analysis. This dramatically reduces data redundancy and hence memory usage and computational cost. During the ILC weight calculation it avoids any need to deconvolve partial sky spherical harmonic coefficients by the beam and pixel window functions, which in strict mathematical sense, is not well-defined for azimuthally symmetric window functions. Using WMAP 9-year and Planck-2015 published frequency maps we obtain a pair of foreground cleaned CMB maps and CMB angular power spectrum. Our power spectrum match well with Planck-2015 results, with some difference. Finally, we show that the weights for ILC foreground minimization have an intrinsic characteristic that it tends to produce a statistically isotropic CMB map as well.