diff --git a/py_wake/site/wasp_grid_site.py b/py_wake/site/wasp_grid_site.py
index f180c0863112abe408f5811bf58a7d9e84047656..38779e2d73ab4f49f263b6dea95a055c25ca4f8f 100644
--- a/py_wake/site/wasp_grid_site.py
+++ b/py_wake/site/wasp_grid_site.py
@@ -37,7 +37,12 @@ class WaspGridSite(XRSite):
 
     def _local_wind(self, localWind, ws_bins=None):
         lw = super()._local_wind(localWind.copy(), ws_bins)
-        lw['TI'] = self.interp(self.ds.tke, lw.coords) * (.75 + 3.8 / lw.ws)
+
+        # ti is assumed to be the turbulence intensity given by CFD
+        # (expected value of TI at 15m/s). The Normal Turbulence model
+        # is used to calculate TI at different wind speed,
+        # see footnote 4 at page 24 of IEC 61400-1 (2005)
+        lw['TI'] = self.interp(self.ds.ti15ms, lw.coords) * (.75 + 3.8 / lw.ws)
         return lw
 
     @classmethod
@@ -58,48 +63,6 @@ def load_wasp_grd(path, globstr='*.grd', speedup_using_pickle=True):
 
     globstr: str
         string that is used to glob files if path is a directory.
-
-    Returns
-    -------
-    WindResourceGrid: :any:`WindResourceGrid`:
-
-    Examples
-    --------
-    >>> from mowflot.wind_resource import WindResourceGrid
-    >>> path = '../mowflot/tests/data/WAsP_grd/'
-    >>> wrg = WindResourceGrid.from_wasp_grd(path)
-    >>> print(wrg)
-        <xarray.Dataset>
-        Dimensions:            (sec: 12, x: 20, y: 20, z: 3)
-        Coordinates:
-          * sec                (sec) int64 1 2 3 4 5 6 7 8 9 10 11 12
-          * x                  (x) float64 5.347e+05 5.348e+05 5.349e+05 5.35e+05 ...
-          * y                  (y) float64 6.149e+06 6.149e+06 6.149e+06 6.149e+06 ...
-          * z                  (z) float64 10.0 40.0 80.0
-        Data variables:
-            flow_inc           (x, y, z, sec) float64 1.701e+38 1.701e+38 1.701e+38 ...
-            ws_mean            (x, y, z, sec) float64 3.824 3.489 5.137 5.287 5.271 ...
-            meso_rgh           (x, y, z, sec) float64 0.06429 0.03008 0.003926 ...
-            obst_spd           (x, y, z, sec) float64 1.701e+38 1.701e+38 1.701e+38 ...
-            orog_spd           (x, y, z, sec) float64 1.035 1.039 1.049 1.069 1.078 ...
-            orog_trn           (x, y, z, sec) float64 -0.1285 0.6421 0.7579 0.5855 ...
-            power_density      (x, y, z, sec) float64 77.98 76.96 193.5 201.5 183.9 ...
-            rix                (x, y, z, sec) float64 0.0 0.0 0.0 0.0 0.0 0.0 0.0 ...
-            rgh_change         (x, y, z, sec) float64 6.0 10.0 10.0 10.0 6.0 4.0 0.0 ...
-            rgh_spd            (x, y, z, sec) float64 1.008 0.9452 0.8578 0.9037 ...
-            f                  (x, y, z, sec) float64 0.04021 0.04215 0.06284 ...
-            tke                (x, y, z, sec) float64 1.701e+38 1.701e+38 1.701e+38 ...
-            A                  (x, y, z, sec) float64 4.287 3.837 5.752 5.934 5.937 ...
-            k                  (x, y, z, sec) float64 1.709 1.42 1.678 1.74 1.869 ...
-            flow_inc_tot       (x, y, z) float64 1.701e+38 1.701e+38 1.701e+38 ...
-            ws_mean_tot        (x, y, z) float64 5.16 6.876 7.788 5.069 6.85 7.785 ...
-            power_density_tot  (x, y, z) float64 189.5 408.1 547.8 178.7 402.2 546.6 ...
-            rix_tot            (x, y, z) float64 0.0 0.0 0.0 9.904e-05 9.904e-05 ...
-            tke_tot            (x, y, z) float64 1.701e+38 1.701e+38 1.701e+38 ...
-            A_tot              (x, y, z) float64 5.788 7.745 8.789 5.688 7.716 8.786 ...
-            k_tot              (x, y, z) float64 1.725 1.869 2.018 1.732 1.877 2.018 ...
-            elev               (x, y) float64 37.81 37.42 37.99 37.75 37.46 37.06 ...
-
     '''
 
     var_name_dict = {
@@ -114,7 +77,7 @@ def load_wasp_grd(path, globstr='*.grd', speedup_using_pickle=True):
         'Roughness changes': 'rgh_change',
         'Roughness speed': 'rgh_spd',
         'Sector frequency': 'f',
-        'Turbulence intensity': 'tke',
+        'Turbulence intensity': 'ti15ms',
         'Weibull-A': 'A',
         'Weibull-k': 'k',
         'Elevation': 'elev',
@@ -269,8 +232,8 @@ def load_wasp_grd(path, globstr='*.grd', speedup_using_pickle=True):
 
     ######################################################################
 
-    if 'tke' in ds and np.mean(ds['tke']) > 1:
-        ds['tke'] *= 0.01
+    if 'ti15ms' in ds and np.mean(ds['ti15ms']) > 1:
+        ds['ti15ms'] *= 0.01
 
     return ds
 
diff --git a/py_wake/tests/test_sites/test_wasp_grid_site.py b/py_wake/tests/test_sites/test_wasp_grid_site.py
index 2c85c0ecf1b4633907dd3f071256d39249f4843b..d484cdd9ac94fed2a0d57e143752c186663f70ea 100644
--- a/py_wake/tests/test_sites/test_wasp_grid_site.py
+++ b/py_wake/tests/test_sites/test_wasp_grid_site.py
@@ -135,9 +135,9 @@ def test_wasp_resources_grid_point(site):
                                           ws=wt_u, ct=wt_ct, power=wt_p, power_unit='kw')
 
     lw = site.local_wind(x, y, 30, wd=range(0, 360, 30))
-    A_lst, k_lst, f_lst, spd_lst, orog_trn_lst, flow_inc_lst, tke_lst = [
+    A_lst, k_lst, f_lst, spd_lst, orog_trn_lst, flow_inc_lst, ti15ms_lst = [
         site.interp(site.ds[n], lw.coords).sel(i=0).values
-        for n in ['Weibull_A', 'Weibull_k', 'Sector_frequency', 'Speedup', 'Turning', 'flow_inc', 'tke']]
+        for n in ['Weibull_A', 'Weibull_k', 'Sector_frequency', 'Speedup', 'Turning', 'flow_inc', 'ti15ms']]
     pdf_lst = [lambda x, A=A, k=k: k / A * (x / A)**(k - 1) * np.exp(-(x / A)**k) * (x[1] - x[0])
                for A, k in zip(A_lst, k_lst)]
 #     cdf_lst = [lambda x, A=A, k=k: 1 - np.exp(-(x / A) ** k) for A, k in zip(A_lst, k_lst)]
@@ -151,7 +151,7 @@ def test_wasp_resources_grid_point(site):
     npt.assert_array_almost_equal(spd_lst, wasp_spd)
     npt.assert_array_almost_equal(orog_trn_lst, wasp_trn)
     npt.assert_array_almost_equal(flow_inc_lst, wasp_inc)
-    npt.assert_array_almost_equal(tke_lst, np.array(wasp_ti) / 100)
+    npt.assert_array_almost_equal(ti15ms_lst, np.array(wasp_ti) / 100)
 
     # compare pdf, u_mean and aep to wasp
     lw = site.local_wind(x, np.array(y) + 1e-6, 30, wd=np.arange(0, 360, 30), ws=ws)