Eectiveness of regular sampling grids to collect multivariate contamination data in urban areas is often strongly reduced by buildings and boundary eects. In addition, earlier observations and knowledge on the history of the area may provide valuable information. This paper extends a simulated annealing-based procedure to optimize the sampling scheme, taking sampling constraints and preliminary information into account. A new optimization criterion is formulated that is able to handle multivariate problems. The sampling scheme is optimized using a spatial weight function that allows to distinguish between areas with dierent priorities. A case study in Rotterdam harbour with ยฎve contaminants at two depths showed two sequential sampling stages, in which two weight functions were applied. The ยฎrst stage combined earlier observations and historical knowledge, with emphasis on areas with high priority. The resulting scheme showed a contamination at 17.4% of the samples, with 1.5% heavily contaminated. The second stage used probability maps of exceeding intermediate threshold values to guide additional sampling to possible hotspots. This yielded 26.7% contaminated samples, with 16.7% being heavily contaminated. This included new locations that were not detected during the ยฎrst stage. The proposed method allows to incorporate important preliminary information, and can be used as a valuable tool in environmental decision-making.