In this paper, we propose a two-stage algorithm utilizing the Cauchy integral and the matrix-valued adaptive Fourier decomposition (abbreviated as matrix AFD) to identify transfer functions of linear time-invariant (LTI) multi-input multi-output (MIMO) systems in the continuous time case. In recent work of Alpay et al. (2017), a theory of adaptive rational approximation to matrix-valued Hardy space functions on the unit disk was established. The matrix-valued function theory has great potential in applications, in views of the practice of its scalar-valued counterparts. The algorithm and application aspects of the mentioned theory of Alpay et al. (2017) have not been developed. The theory was only written for the unit disk case corresponding to the discrete time systems. The contributions of the present paper are 3-fold. First, we construct an analogous adaptive approximation theory for complex matrix-valued Hardy space functions defined on a half of the complex plane, corresponding to the Laplace transforms of signals of finite energy whose Fourier transforms are supported on a half of the frequency domain. The half plane model corresponds to signals defined in the whole axis range which is an alternative case to signals defined in a compact interval. The second fold contribution lays on maximal selection of the pair (a,P) where a is a point of the right-half plane and P is an orthogonal projection. We show that the optimal selection of P is dependent on a when a is first fixed, that is P=P(a), where P:a→P(a) has an explicit corresponding relation. Due to this relation we reduce the maximal selection of the pair (a,P) to only that of the parameter a. This result can be extended to the compact intervals case as well. The third fold is, with the precise rule from a to P(a), we develop a practical algorithm for the adaptive approximation to the transfer function. Through an example we show that the proposed algorithm is effective in both the noise-free and noisy cases.