The principal aim of the paper is to investigate new integral expression $${int}_0^{infty}x^{s+1}e^{-{sigma}x^2}L_m^{(gamma,delta)};({zeta};{sigma}x^2);L_n^{(alpha,eta)};({xi};{sigma}x^2);J_s;(xy);dx$$, where $y$ is a positive real number; $sigma$, $zeta$ and $xi$ are complex numbers with positive real parts; $s$, $alpha$, $eta$, $gamma$ and $delta$ are complex numbers whose real parts are greater than -1; $J_n(x)$ is Bessel function and $L_n^{(alpha,eta)}$ (${gamma};x$) is generalized Laguerre polynomials. Some integral formulas have been obtained. The Maple implementation has also been examined.