The purpose of this article is to introduce a new Lévy process, termed Variance Gamma++ process, to model the dynamic of assets in illiquid markets. Such a process has the mathematical tractability of the Variance Gamma process and is obtained applying the self-decomposability of the gamma law. Compared to the Variance Gamma model, it has an additional parameter representing the measure of the trading activity. We give a full characterization of the Variance Gamma++ process in terms of its characteristic triplet, characteristic function and transition density. In addition, we provide efficient path simulation algorithms, both forward and backward in time. We also obtain an efficient "integral-free" explicit pricing formula for European options. These results are instrumental to apply Fourier-based option pricing and maximum likelihood techniques for the parameter estimation. Finally, we apply our model to illiquid markets, namely to the calibration of European power future market data. We accordingly evaluate exotic derivatives using the Monte Carlo method and compare these values to those obtained using the Variance Gamma process and give an economic interpretation of the obtained results. Finally, we illustrate an extension to the multivariate framework.