The apparently random pattern of the non-trivial zeroes of the Riemann zeta function (all on the critical line, according to the Riemann hypothesis) has led to the suggestion that they may be related to the spectrum of an operator. It has also been known for some time that the statistical properties of the eigenvalue distribution of an ensemble of random matrices resemble those of the zeroes of the zeta function. With the objective to identify a suitable operator, we start by assuming the Riemann hypothesis and construct a unitary matrix model (UMM) for the zeta function. Our approach, however, could be termed piecemeal, in the sense that, we consider each factor (in the Euler product representation) of the zeta function to get a UMM for each prime, and then assemble these to get a matrix model for the full zeta function. This way we can write the partition function as a trace of an operator. Similar construction works for a family of related zeta functions.