The development of symbolic methods for the factorization and integration of linear PDEs, many of the methods being generalizations of the Laplace transformations method, requires the finding of complete generating sets of invariants for the corresponding linear operators and their systems with respect to the gauge transformations \(L -> g(x,y)-^1 O\) \(L\) \(O\) \(g(x,y)^{-1}\). Within the theory of Laplace-like methods, there is no uniform approach to this problem, though some individual invariants for hyperbolic bivariate operators, and complete generating sets of invariants for second- and third-order hyperbolic bivariate ones have been obtained. Here we demonstrate a systematic and much more efficient approach to the same problem by ap...