Gaussian states are the backbone of quantum information protocols with continuous variable systems, whose power relies fundamentally on the entanglement between the different modes. In the case of global pure states, knowledge of the reduced states in a given bipartition of a multipartite quantum system bears information on the entanglement in such bipartition. For Gaussian states, the reduced states are also Gaussian, so there determination requires essentially the experimental determination of their covariance matrix. Here, we develop strategies to determine the covariance matrix of an arbitrary n-mode bosonic Gaussian state through measurement of the total phase acquired when appropriate metaplectic evolutions, associated with quadratic Hamiltonians, are applied. Simply one-mode metaplectic evolutions, such rotations, squeezing and shear transformations, in addition to a single two-mode rotation, allows to determine all the covariance matrix elements of a n-mode bosonic system. All the single-mode metaplectic evolutions are applied conditionally to a state in which an ancilla qubit is entangled with the n-mode system. The ancillary system provides, after measurement, the value of the total phase of each evolution. The proposed method is experimentally friendly to be implemented in the most currently used continuous variable systems.