diff --git a/DESCRIPTION b/DESCRIPTION
index b0343f856f6c7115e9f91c078cb9565320731f95..8c06d06c2f69025fc229e2b899f24fc0ac929bfd 100644
--- a/DESCRIPTION
+++ b/DESCRIPTION
@@ -1,6 +1,6 @@
 Package: gmGeostats
-Version: 0.10-7
-Date: 2020-09-28
+Version: 0.10-7.9003
+Date: 2020-10-05
 Title: Geostatistics for Compositional Analysis
 Authors@R: c(person(given = "Raimon", 
         family = "Tolosana-Delgado", 
diff --git a/NAMESPACE b/NAMESPACE
index f607aa58cfb8a93004bf7a79cf911a11992b132c..98271219ce69c0f9b152e345b78ac09091d3c50e 100644
--- a/NAMESPACE
+++ b/NAMESPACE
@@ -177,6 +177,7 @@ export(gridOrder_gstat)
 export(gridOrder_sp)
 export(gsi.DS)
 export(gsi.EVario2D)
+export(gsi.TurningBands)
 export(gsi.gstatCokriging2compo)
 export(gsi.gstatCokriging2rmult)
 export(gsi.produceV)
diff --git a/NEWS.md b/NEWS.md
index 6b948574052ca32421d81281e544dbc7d50a02f0..c6d5443fbc1caadb70f79b97964cbbb9bfaecf8f 100644
--- a/NEWS.md
+++ b/NEWS.md
@@ -1,3 +1,8 @@
+# gmGeostats 0.10.7.9003
+
+* (2021-04-08) exported turning bands for direct usage
+* (2020-10-01) bug in turning bands spherical variogram corrected
+
 # gmGeostats 0.10.7
 
 * (2020-09-27) updated references in generalised diagonalisations; improved/reimplemented `as.logratioVariogram()` and `as.logratioVariogramAnisotropy()` methods; minor edits at `plot.logratioVariogramAnisotropy()` and `variogramModelPlot.logratioVariogram()`
diff --git a/R/geostats.R b/R/geostats.R
index 45d1e542c2d4105ba4cba70c22a82ffa5f749520..32f495ff92db037ee20874af4fe8d0807a26ff97 100644
--- a/R/geostats.R
+++ b/R/geostats.R
@@ -87,6 +87,7 @@ gsi.calcCgram <- function(X,Y,vgram,ijEqual=FALSE) {
 #' 
 #' @return an array with (npoint, nvar, nsim)-elements, being npoint=nrow(X)
 #' and nvar = nr of variables in vgram
+#' @export
 #' @useDynLib gmGeostats CMVTurningBands
 gsi.TurningBands <- function(X,vgram,nBands,nsim=NULL) {
   # checks and preps
diff --git a/R/gmSimulation.R b/R/gmSimulation.R
index 1e5068b92c20d6497d1e8fcc4bbf089b822c728a..88db6fdd35c9f028b9214bd584f0a5f68177b177 100644
--- a/R/gmSimulation.R
+++ b/R/gmSimulation.R
@@ -249,7 +249,7 @@ gsi.DS4CoDa <- function(n, f, t, n_realiz, nx_TI, ny_TI, nx_SimGrid, ny_SimGrid,
 #' image_cokriged(SimGrid, ivar="V1", breaks=o1$breaks, col=o1$col)
 #' image_cokriged(SimGrid, ivar="V2", breaks=o2$breaks, col=o2$col)
 #' image_cokriged(SimGrid, ivar="mask", breaks=c(-0.0001, 0.5, 1.001))
-#' res = gsi.DS(n=5, f=0.75, t=0.05, n_realiz=2, dim_TI=c(10,7),  dim_SimGrid=c(10,7), 
+#' res = gsi.DS(n=5, f=0.75, t=0.05, n_realiz=100, dim_TI=c(10,7),  dim_SimGrid=c(10,7), 
 #'        TI_input=as.matrix(TI_input), SimGrid_input=as.matrix(SimGrid), 
 #'        ivars_TI = c("V1", "V2"), SimGrid_mask="mask", invertMask=TRUE)
 #' image_cokriged(cbind(xy_TI, getStackElement(res,1)), ivar="V1", breaks=o1$breaks, col=o1$col)
diff --git a/R/zzz.R b/R/zzz.R
index 8c88d2610d2acf67277cda05524c9269c8b9f158..2a0ae9ba58e6a588a064a3d4e575487008531d09 100644
--- a/R/zzz.R
+++ b/R/zzz.R
@@ -200,6 +200,7 @@ setClassUnion(name="GridOrNothing", members = c("NULL", "GridTopology"))
 
 
 .onLoad <- function(libname, pkgname){
+
   ## set package options ---- 
     # grid organisation
     gridOrder = list(refpoint="topleft", cycle=1:2)
@@ -215,7 +216,8 @@ setClassUnion(name="GridOrNothing", members = c("NULL", "GridTopology"))
   invisible()  
 }
 
-.onUnload <- function(libname, pkgname){
+.onUnload <- function(libpath){
+  library.dynam.unload("gmGeostats", libpath)
   ## remove package options ---- 
   options(gmGeostats=NULL)
   
diff --git a/src/gmGeostats.c b/src/gmGeostats.c
index 58d8e05a1a6c87d89fd7a82c21dcf0a52ba63096..4b689abfeac8b844bd89ce9b73bdb1e0fac5697a 100644
--- a/src/gmGeostats.c
+++ b/src/gmGeostats.c
@@ -857,6 +857,14 @@ extern void anaBackwardC(const int *dimX,
 
 
 
+/* ======================================
+ * EXPORT FUNCTIONS TO R
+ * ======================================
+ */
+
+
+
+
 static R_CMethodDef cMethods[] = {
   {"CcalcCgram", (DL_FUNC) &CcalcCgram, 15, CcalcCgram_t},
   {"CMVTurningBands", (DL_FUNC) &CMVTurningBands, 11, CMVTurningBands_t},
@@ -866,11 +874,19 @@ static R_CMethodDef cMethods[] = {
   {NULL, NULL, 0}
 };
 
+//{"CMVTurningBands2", (DL_FUNC) &CMVTurningBands2, 11, CMVTurningBands2_t},
+
 
-void R_init_compositions(DllInfo *info)
+/*
+void R_init_gmGeostats(DllInfo *info)
 {
   R_registerRoutines(info, cMethods, NULL, NULL, NULL);
   R_useDynamicSymbols(info, FALSE);
   R_forceSymbols(info, TRUE);
 }
+*/
+
+
+
+