We study matrix elements of Fourier-transformed straight infinite Wilson lines as a way to calculate gauge invariant tree-level amplitudes with off-shell gluons. The off-shell gluons are assigned "polarization vectors" which (in the Feynman gauge) are transverse to their off-shell momenta and define the direction of the corresponding Wilson line operators. The infinite Wilson lines are first regularized to prove the correctness of the method. We have implemented the method in a computer FORM program that can calculate gluonic matrix elements of Wilson line operators automatically. In addition we formulate the Feynman rules that are convenient in certain applications, e.g. proving the Ward identities. Using both the program and the Feynman rules we calculate a few examples, in particular the matrix elements corresponding to gauge invariant and processes. An immediate application of the approach is in the high energy scattering, as in a special kinematic setup our results reduce to the form directly related to Lipatov's vertices. Thus the results we present can be directly transformed into Lipatov's vertices, in particular into and vertices with arbitrary "orientation" of reggeized gluons. Since the formulation itself is not restricted to high-energy scattering, we also apply the method to a decomposition of an ordinary on-shell amplitude into a set of gauge invariant objects.