Simulation of t copula in Python

I am trying to simulate a t-copula using Python, but my code yields strange results (is not well-behaving):

I followed the approach suggested by Demarta & McNeil (2004) in "The t Copula and Related Copulas", which states:

By intuition, I know that the higher the degrees of freedom parameter, the more the t copula should resemble the Gaussian one (and hence the lower the tail dependency). However, given that I sample from scipy.stats.invgamma.rvs or alternatively from scipy.stats.chi2.rvs, yields higher values for my parameter s the higher my parameter df. This does not made any sense, as I found multiple papers stating that for df--> inf, t-copula --> Gaussian copula.

Here is my code, what am I doing wrong? (I'm a beginner in Python fyi).

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from scipy.stats import invgamma, chi2, t

#Define number of sampling points
n_samples = 1000
df = 10

calib_correl_matrix = np.array([[1,0.8,],[0.8,1]]) #I just took a bivariate correlation matrix here
mu = np.zeros(len(calib_correl_matrix))
s = chi2.rvs(df)
#s = invgamma.pdf(df/2,df/2) 
Z = np.random.multivariate_normal(mu, calib_correl_matrix,n_samples)
X = np.sqrt(df/s)*Z #chi-square method
#X = np.sqrt(s)*Z #inverse gamma method

U = t.cdf(X,df)

My outcomes behave exactly oppisite to what I am (should be) expecting: Higher df create much higher tail-dependency, here also visually:

 U_pd = pd.DataFrame(U)
 fig = plt.gcf()
 fig.set_size_inches(14.5, 10.5)
 pd.plotting.scatter_matrix(U_pd, figsize=(14,10), diagonal = 'kde')
 plt.show()

df=4:

df=100:

It gets especially worse when using the invgamma.rvs directly, even though they should yield the same. For dfs>=30 I often receive a ValueError ("ValueError: array must not contain infs or NaNs")

Thank you very much for your help, much appreciated!

from Simulation of t copula in Python