From 75a4d699cd904d954194fe0ac92c765625c8af1c Mon Sep 17 00:00:00 2001 From: ahmathlete Date: Sat, 13 Apr 2024 11:11:33 +0200 Subject: [PATCH] Rpackage --- README.md | 4 +- RPackage/NAMESPACE | 1 - RPackage/R/cape_f.R | 25 ++++++----- RPackage/R/cape_r.R | 1 - RPackage/R/conv2D_f.R | 7 +--- RPackage/R/conv2D_r.R | 1 - RPackage/R/xcorr2D_f.R | 7 +--- RPackage/R/xcorr2D_r.R | 1 - RPackage/man/AquaFortR-package.Rd | 5 ++- RPackage/man/saturation_vapour_pressure.Rd | 2 +- RPackage/src/AquaFortRmodulec.c | 48 ++++++++++------------ RPackage/src/AquaFortRmodulef.f90 | 1 + 12 files changed, 45 insertions(+), 58 deletions(-) diff --git a/README.md b/README.md index 83c2271..3766b23 100644 --- a/README.md +++ b/README.md @@ -7,7 +7,7 @@ [![License](https://img.shields.io/badge/License-CCBY-blue)](#license) -This is the official repository of the Project AquaFortR: Streamlining Atmospheric Science, Oceanography, Climate, and Water Research with Fortran-accelerated R. +This is the official repository of the project AquaFortR: Streamlining Atmospheric Science, Oceanography, Climate, and Water Research with Fortran-accelerated R. The repository is structured as follows: @@ -38,7 +38,7 @@ swirl::install_course_zip("path/to/AquaFortR_Swirl.zip") ### Book -Materials for the chapter 2 are available here. Please, +Materials for Chapter 2 are available here. Please, revise the path to the shared libraries files in the R-Fortran functions. ## Acknowledgment diff --git a/RPackage/NAMESPACE b/RPackage/NAMESPACE index b1cc5ff..524b40e 100644 --- a/RPackage/NAMESPACE +++ b/RPackage/NAMESPACE @@ -8,5 +8,4 @@ export(mixing_ratio_from_dewpoint) export(saturation_vapour_pressure) export(xcorr2D_f) export(xcorr2D_r) -useDynLib(AquaFortR) useDynLib(AquaFortR, .registration=TRUE) diff --git a/RPackage/R/cape_f.R b/RPackage/R/cape_f.R index 5a4bc8e..0de0010 100644 --- a/RPackage/R/cape_f.R +++ b/RPackage/R/cape_f.R @@ -26,27 +26,26 @@ #' vtc = TRUE #' ) #' @author Klemens Barfus (Original in Fortran), Ahmed Homoudi (Integration in R) -#' @useDynLib AquaFortR #' @export cape_f <- function(t_parcel, dwpt_parcel, mr_parcel, p_profile, t_profile, mr_profile, vtc = TRUE) { - DIM <- c( - length(p_profile), - 4, - as.integer(vtc) - ) + + vtc <-as.integer(vtc) + result <- .Call( c_cape_f, - as.double(t_parcel), - as.double(dwpt_parcel), - as.double(mr_parcel), - as.double(p_profile), - as.double(t_profile), - as.double(mr_profile), - as.integer(DIM) + t_parcel, + dwpt_parcel, + mr_parcel, + p_profile, + t_profile, + mr_profile, + vtc ) names(result) <- c("CAPE", "CIN", "p_LCL", "p_LFC") return(result) } + + diff --git a/RPackage/R/cape_r.R b/RPackage/R/cape_r.R index 5ae2c05..06c4036 100644 --- a/RPackage/R/cape_r.R +++ b/RPackage/R/cape_r.R @@ -26,7 +26,6 @@ #' vtc = TRUE #' ) #' @author Klemens Barfus (Original in Fortran), Ahmed Homoudi (Translation to R) -#' @useDynLib AquaFortR #' @export cape_r <- function(t_parcel, dwpt_parcel, mr_parcel, p_profile, t_profile, mr_profile, diff --git a/RPackage/R/conv2D_f.R b/RPackage/R/conv2D_f.R index 08889bd..c1816b5 100644 --- a/RPackage/R/conv2D_f.R +++ b/RPackage/R/conv2D_f.R @@ -12,16 +12,13 @@ #' b <- matrix(c(5, 6, 7, 8), ncol = 2) #' conv2D_f(a, b) #' @author Ahmed Homoudi -#' @useDynLib AquaFortR #' @export conv2D_f <- function(a, b) { stopifnot(length(dim(a)) == 2 | length(dim(b)) == 2) result <- .Call( c_conv2d_f, - as.integer(dim(a)), - as.double(a), - as.integer(dim(b)), - as.double(b) + a, + b ) return(result) } diff --git a/RPackage/R/conv2D_r.R b/RPackage/R/conv2D_r.R index d3e8992..c9e5179 100644 --- a/RPackage/R/conv2D_r.R +++ b/RPackage/R/conv2D_r.R @@ -12,7 +12,6 @@ #' b <- matrix(c(5, 6, 7, 8), ncol = 2) #' conv2D_f(a, b) #' @author Ahmed Homoudi -#' @useDynLib AquaFortR #' @export conv2D_r <- function(a, b) { stopifnot(length(dim(a)) == 2 | length(dim(b)) == 2) diff --git a/RPackage/R/xcorr2D_f.R b/RPackage/R/xcorr2D_f.R index 1f0b1b7..2552f67 100644 --- a/RPackage/R/xcorr2D_f.R +++ b/RPackage/R/xcorr2D_f.R @@ -12,16 +12,13 @@ #' b <- matrix(c(5, 6, 7, 8), ncol = 2) #' xcorr2D_f(a, b) #' @author Ahmed Homoudi -#' @useDynLib AquaFortR #' @export xcorr2D_f <- function(a, b) { stopifnot(length(dim(a)) == 2 | length(dim(b)) == 2) result <- .Call( c_xcorr2d_f, - as.integer(dim(a)), - as.double(a), - as.integer(dim(b)), - as.double(b) + a, + b ) return(result) } diff --git a/RPackage/R/xcorr2D_r.R b/RPackage/R/xcorr2D_r.R index 1283cdc..bdfce81 100644 --- a/RPackage/R/xcorr2D_r.R +++ b/RPackage/R/xcorr2D_r.R @@ -12,7 +12,6 @@ #' b <- matrix(c(5, 6, 7, 8), ncol = 2) #' xcorr2D_f(a, b) #' @author Ahmed Homoudi -#' @useDynLib AquaFortR #' @export xcorr2D_r <- function(a, b) { stopifnot(length(dim(a)) == 2 | length(dim(b)) == 2) diff --git a/RPackage/man/AquaFortR-package.Rd b/RPackage/man/AquaFortR-package.Rd index d5ddaf9..3fc5409 100644 --- a/RPackage/man/AquaFortR-package.Rd +++ b/RPackage/man/AquaFortR-package.Rd @@ -6,12 +6,13 @@ \alias{AquaFortR-package} \title{AquaFortR: Streamlining Atmospheric Science, Oceanography, Climate, and Water Research with Fortran-accelerated R} \description{ -More about what it does (maybe more than one line) Use four spaces when indenting paragraphs within the Description. +AquaFortR is an educational package for students and researchers in Atmospheric Science, Oceanography, Climate, and Water Research. AquaFortR aims to show that using simple Fortran scripts can fulfil the need to speed up R, especially considering that most data sets in these fields are large discretised arrays representing earth system processes. Fortran is one of the fastest programming languages in number crunching, if not the fastest, and multidimensional arrays are a core part of Fortran. } \seealso{ Useful links: \itemize{ - \item \url{https://AHomoudi.github.io/AquaFortR/} + \item \url{https://AHomoudi.github.io/AquaFortR} + \item Report bugs at \url{https://github.com/AHomoudi/AquaFortR/issues} } } diff --git a/RPackage/man/saturation_vapour_pressure.Rd b/RPackage/man/saturation_vapour_pressure.Rd index edceabb..5bcc90e 100644 --- a/RPackage/man/saturation_vapour_pressure.Rd +++ b/RPackage/man/saturation_vapour_pressure.Rd @@ -15,5 +15,5 @@ saturation_vapour_pressure(temperature, ice = FALSE) Estimation of Saturation vapour pressure [hPa] from temperature [k]. } \author{ -Klemens Barfus +Klemens Barfus (Original in Fortran), translated to R by Ahmed Homoudi } diff --git a/RPackage/src/AquaFortRmodulec.c b/RPackage/src/AquaFortRmodulec.c index d56ae0e..409209c 100644 --- a/RPackage/src/AquaFortRmodulec.c +++ b/RPackage/src/AquaFortRmodulec.c @@ -1,47 +1,45 @@ +#define R_NO_REMAP #include #include #include // for NULL -#include -#include -#include void F77_NAME(xcorr2d_f)(int m, int n, int p, int q, int k, int l, double *a, double *b, double *cc); void F77_NAME(conv2d_f)(int m, int n, int p, int q, int k, int l, double *a, double *b, double *conv); void F77_NAME(cape_f)(double *t_parcel, double *dwpt_parcel, double *mr_parcel, int nlevels, double *p_profile, double *t_profile, double *mr_profile, int vtc, int nresult, double *result); -extern SEXP c_xcorr2d_f(SEXP dim_a, SEXP a, SEXP dim_b, SEXP b) +extern SEXP c_xcorr2d_f(SEXP a, SEXP b) { - int m = INTEGER(dim_a)[0]; - int n = INTEGER(dim_a)[1]; + int m = Rf_nrows(a); + int n = Rf_ncols(a); // - int p = INTEGER(dim_b)[0]; - int q = INTEGER(dim_b)[1]; + int p = Rf_nrows(b); + int q = Rf_ncols(b); // int k = m + p - 1; int l = n + q - 1; SEXP ret; - PROTECT(ret = allocMatrix(REALSXP, k, l)); + PROTECT(ret = Rf_allocMatrix(REALSXP, k, l)); F77_CALL(xcorr2d_f) (m, n, p, q, k, l, REAL(a), REAL(b), REAL(ret)); UNPROTECT(1); return (ret); } -extern SEXP c_conv2d_f(SEXP dim_a, SEXP a, SEXP dim_b, SEXP b) +extern SEXP c_conv2d_f(SEXP a, SEXP b) { - int m = INTEGER(dim_a)[0]; - int n = INTEGER(dim_a)[1]; + int m = Rf_nrows(a); + int n = Rf_ncols(a); // - int p = INTEGER(dim_b)[0]; - int q = INTEGER(dim_b)[1]; + int p = Rf_nrows(b); + int q = Rf_ncols(b); // int k = m + p - 1; int l = n + q - 1; SEXP ret; - PROTECT(ret = allocMatrix(REALSXP, k, l)); + PROTECT(ret = Rf_allocMatrix(REALSXP, k, l)); F77_CALL(conv2d_f) (m, n, p, q, k, l, REAL(a), REAL(b), REAL(ret)); UNPROTECT(1); @@ -50,34 +48,32 @@ extern SEXP c_conv2d_f(SEXP dim_a, SEXP a, SEXP dim_b, SEXP b) extern SEXP c_cape_f(SEXP t_parcel, SEXP dwpt_parcel, SEXP mr_parcel, SEXP p_profile, SEXP t_profile, SEXP mr_profile, - SEXP DIM) + SEXP vtc) { - int nlevels = INTEGER(DIM)[0]; - int nret = INTEGER(DIM)[1]; - int vtc = INTEGER(DIM)[2]; + int nlevels = Rf_length(t_profile); + int vtc2 = INTEGER(vtc)[0]; + int nret = 4; SEXP ret; /* Order: CAPE, CIN, p_LCL, p_LFC */ - PROTECT(ret = allocVector(REALSXP, 4)); + PROTECT(ret = Rf_allocVector(REALSXP, 4)); F77_CALL(cape_f) (REAL(t_parcel), REAL(dwpt_parcel), REAL(mr_parcel), nlevels, - REAL(p_profile), REAL(t_profile), REAL(mr_profile), vtc, + REAL(p_profile), REAL(t_profile), REAL(mr_profile), vtc2, nret, REAL(ret)); UNPROTECT(1); return (ret); } static const R_CallMethodDef CallEntries[] = { - {"c_xcorr2d_f", (DL_FUNC)&c_xcorr2d_f, 4}, - {"c_conv2d_f", (DL_FUNC)&c_conv2d_f, 4}, + {"c_xcorr2d_f", (DL_FUNC)&c_xcorr2d_f, 2}, + {"c_conv2d_f", (DL_FUNC)&c_conv2d_f, 2}, {"c_cape_f", (DL_FUNC)&c_cape_f, 7}, {NULL, NULL, 0}}; -static const R_FortranMethodDef FEntries[] = { - {NULL, NULL, 0}}; void R_init_AquaFortR(DllInfo *dll) { - R_registerRoutines(dll, NULL, CallEntries, FEntries, NULL); + R_registerRoutines(dll, NULL, CallEntries, NULL, NULL); R_useDynamicSymbols(dll, FALSE); } diff --git a/RPackage/src/AquaFortRmodulef.f90 b/RPackage/src/AquaFortRmodulef.f90 index be594c6..1c4a4c6 100644 --- a/RPackage/src/AquaFortRmodulef.f90 +++ b/RPackage/src/AquaFortRmodulef.f90 @@ -4,6 +4,7 @@ module AquaFortRmodule public :: xcorr2d_f public :: conv2d_f + public :: cape_f contains ! Cross-Correlation