From db6a09695c0b61d6ab2a9eebe3d72923f08b044c Mon Sep 17 00:00:00 2001 From: Ralf Stubner Date: Fri, 21 Dec 2018 00:17:05 +0100 Subject: [PATCH] Re-add vectorized functions This reverts commit c8161ccd0fc8fe718478cac23c105f134447e048. --- NAMESPACE | 27 + R/vectorisedSEfunctions.R | 1094 ++++++++++++++++++++++++++++++++ man/Vectorised.Rd | 212 +++++++ vignettes/VectorisedSwephR.Rmd | 76 +++ 4 files changed, 1409 insertions(+) create mode 100644 R/vectorisedSEfunctions.R create mode 100644 man/Vectorised.Rd create mode 100644 vignettes/VectorisedSwephR.Rmd diff --git a/NAMESPACE b/NAMESPACE index abe955c9..740e4db2 100644 --- a/NAMESPACE +++ b/NAMESPACE @@ -35,5 +35,32 @@ export(swe_sol_eclipse_when_loc) export(swe_topo_arcus_visionis) export(swe_version) export(swe_vis_limit_mag) +export(vec_azalt) +export(vec_azalt_rev) +export(vec_calc) +export(vec_calc_ut) +export(vec_date_conversion) +export(vec_day_of_week) +export(vec_deltat) +export(vec_deltat_ex) +export(vec_fixstar2) +export(vec_fixstar2_mag) +export(vec_fixstar2_ut) +export(vec_get_planet_name) +export(vec_heliacal_angle) +export(vec_heliacal_pheno_ut) +export(vec_heliacal_ut) +export(vec_julday) +export(vec_lun_eclipse_how) +export(vec_lun_eclipse_when) +export(vec_lun_eclipse_when_loc) +export(vec_pheno) +export(vec_pheno_ut) +export(vec_refrac_extended) +export(vec_revjul) +export(vec_rise_trans_true_hor) +export(vec_sol_eclipse_when_loc) +export(vec_topo_arcus_visionis) +export(vec_vis_limit_mag) importFrom(Rcpp,evalCpp) useDynLib(swephR, .registration = TRUE) diff --git a/R/vectorisedSEfunctions.R b/R/vectorisedSEfunctions.R new file mode 100644 index 00000000..c45d4e36 --- /dev/null +++ b/R/vectorisedSEfunctions.R @@ -0,0 +1,1094 @@ +# Copyright 2018 Ralf Stubner and Victor Reijs +# +# This file is part of swephR. +# +# swephR is free software: you can redistribute it and/or modify it +# under the terms of the GNU Affero General Public License as published by +# the Free Software Foundation, either version 3 of the License, or +# (at your option) any later version. +# +# swephR is distributed in the hope that it will be useful, but +# WITHOUT ANY WARRANTY; without even the implied warranty of +# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +# GNU Affero General Public License for more details. +# +# You should have received a copy of the GNU Affero General Public License +# along with swephR. If not, see . + + +#section2 +##' @title Vectorised swephR functions +##' @description This are the vectorised functions of the +##' ones that are namend swe_xxx. Each input parameter can be a vector. +##' Make sure that the vector length of each input parameter matches +##' the rules around data-frames: so the largest vector length must be an integer +##' multiple of the other vector lengths. +##' @param ipl Body/planet as integer vector (SE$SUN=0, SE$Moon=1, ... SE$PLUTO=9) +##' @param iflag Computation flag as integer vector, many options possible (section 2.3) +##' @param jd_ut UT Julian day number vector (day) +##' @param jd_et ET Julian day number as double vector (day) +##' @param starname Star name as string vector ("" for no star) +##' @param calc_flag Calculation flag as integer vector (refraction direction (SE$TRUE_TO_APP=0 or SE$APP_TO_TRUE=1)) +##' @param coord_flag Coordinate flag as integer vector (reference system (SE$ECL2HOR=0 or SE$EQU2HOR=1)) +##' @param atpress Atmospheric pressure as double vector (hPa) +##' @param attemp Atmospheric temperature as double vector (Celsius) +##' @param athum Atmospheric humidity as double vector (\%) +##' @param atvis Atmospheric visibiliy as double vector (km or -) +##' @param obsage Age of observer as double vector (year) +##' @param obssnellen Aquity of observer as double vector (-) +##' @param obsbin Mono-/bi-nocular observation as double vector (-) +##' @param obsmag Magnification as double vector (-) +##' @param obsaper Aperture of optics as double vector (mm) +##' @param obstrans Transmission as double vector (-) +##' @param ephe_flag Ephemeris flag as integer vector (SE$FLG_JPLEPH=1, SE$FLG_SWIEPH=2 or SE$FLG_MOSEPH=4) +##' @param horhgt Horizon apparent altitude as double vector (deg) +##' @param xin1 Position of body as numeric vector (either ecliptical or equatorial coordinates, depending on coord_flag) +##' @param xin2 Position of body as numeric vector (either ecliptical or equatorial coordinates, depending on coord_flag) +##' @param rsmi Event flag as integer vector (e.g.: SE$CALC_RISE=1, SE$CALC_SET=2,SE$CALC_MTRANSIT=4,SE$CALC_ITRANSIT=8) +##' @param backward backwards search as boolean vector (TRUE) +##' @param ifltype eclipse type as integer vector (e.g.: SE$ECL_CENTRAL=1,SE$ECL_NONCENTRAL=2,SE$ECL_TOTAL=4,SE$ECL_ANNULAR=8,SE$ECL_PARTIAL=16,SE$ECL_ANNULAR_TOTAL=32) +##' @param InAlt object's apparent/topocentric altitude as double vector (depending on calc_flag) (deg) +##' @param lapse_rate lapse rate as double vector (K/m) +##' @param jd_utstart UT Julian day number as double vector (day) +##' @param jd_start Julian day number as double vector (day) +##' @param objectname Name of fixed star or planet as string vector +##' @param event_type Event type as integer vector +##' @param helflag Calculation flag (incl. ephe_flag values) as integer vector +##' @param mag Object's visible magnitude (Vmag) as double vector (-) +##' @param AziO Object's azimuth as double vector (deg) +##' @param AltO Object's altitude as double vecor (deg) +##' @param AziS Sun's azimuth as double (vector deg) +##' @param AziM Moon's azimut as double vecor (deg) +##' @param AltM Moon's altitude as double vector (deg) +##' @param year Astronomical year as integer vector +##' @param month Month as integer vector +##' @param day Day as integer vector +##' @param hour Hour as double vector +##' @param gregflag Calendar type as integer vector (SE$JUL_CAL=0 or SE$GREG_CAL=1) +##' @param cal Calendar type as char vector ("g"[regorian] or "j"[ulian]) +##' @param jd Julian day number as double vector (day) +##' @param longitude Geographic longitude as double vector (deg) +##' @param lat Geographic latitude as double vector (deg) +##' @param height Height as double vector (m) +##' @rdname Vectorised +##' @export +vec_calc_ut <- + function(jd_ut, + ipl, + iflag = 4) { + # default Moshier epheemris + functionvector <- + data.frame(jd_ut, ipl, iflag) + #print(functionvector) + listsize <- nrow(functionvector) + ResultVector <- vector("list", listsize) + for (i in 1:listsize) + { + ResultVector[[i]] <- swe_calc_ut(functionvector$jd_ut[i], + functionvector$ipl[i], + functionvector$iflag[i]) + } + class(ResultVector)<-"swephR.object.ephe" + return(ResultVector) + } + +##' @rdname Vectorised +##' @export +vec_calc <- + function(jd_et, + ipl, + iflag = 4) { + #default Moshier epheemris + functionvector <- + data.frame(jd_et, ipl, iflag) +# ##print(functionvector) + listsize <- nrow(functionvector) + ResultVector <- vector("list", listsize) + for (i in 1:listsize) + { + ResultVector[[i]] <- swe_calc(functionvector$jd_et[i], + functionvector$ipl[i], + functionvector$iflag[i]) + } + class(ResultVector)<-"swephR.object.ephe" + return(ResultVector) + + } + +#section3 +##' @rdname Vectorised +##' @export +vec_get_planet_name <- function(ipl) { + functionvector <- data.frame(ipl) + #print(functionvector) + listsize <- nrow(functionvector) + ResultVector <- vector("list", listsize) + for (i in 1:listsize) + { + ResultVector[[i]] <- swe_get_planet_name(functionvector$ipl[i]) + } + return(ResultVector) +} + +#section4 +##' @rdname Vectorised +##' @export +vec_sol_eclipse_when_loc <- + function(jd_start, + ephe_flag = 4, + # default Moshier epheemris + longitude, + lat, + height = 0, + backward = FALSE) { + #default is forward search + functionvector <- + data.frame(jd_start, ephe_flag, longitude, lat, height, backward) + #print(functionvector) + listsize <- nrow(functionvector) + ResultVector <- vector("list", listsize) + for (i in 1:listsize) + { + geopos <- + c(functionvector$longitude[i], + functionvector$lat[i], + functionvector$height[i]) + ResultVector[[i]] <- swe_sol_eclipse_when_loc( + functionvector$jd_start[i], + functionvector$ephe_flag[i], + geopos, + functionvector$backward[i] + ) + } + class(ResultVector)<-"swephR.eclipse.loc" + return(ResultVector) + } + +#section4 +##' @rdname Vectorised +##' @export +vec_fixstar2_ut <- + function(starname, + jd_ut, + iflag = 4) { + # Default Moshier ephemeris + functionvector <- + data.frame(starname, jd_ut, iflag, stringsAsFactors = FALSE) + #print(functionvector) + listsize <- nrow(functionvector) + ResultVector <- vector("list", listsize) + for (i in 1:listsize) + { + ResultVector[[i]] <- swe_fixstar2_ut(functionvector$starname[i], + functionvector$jd_ut[i], + functionvector$iflag[i]) + } + class(ResultVector)<-"swephR.star.ephe" + return(ResultVector) + } + +##' @rdname Vectorised +##' @export +vec_fixstar2 <- + function(starname, + jd_et, + iflag = 4) { + # Default Moshier ephemeris + functionvector <- + data.frame(starname, jd_et, iflag, stringsAsFactors = FALSE) + #print(functionvector) + listsize <- nrow(functionvector) + ResultVector <- vector("list", listsize) + for (i in 1:listsize) + { + ResultVector[[i]] <- swe_fixstar2(functionvector$starname[i], + functionvector$jd_et[i], + functionvector$iflag[i]) + } + class(ResultVector)<-"swephR.star.ephe" + return(ResultVector) + } + +##' @rdname Vectorised +##' @export +vec_fixstar2_mag <- + function(starname) { + functionvector <- + data.frame(starname, stringsAsFactors = FALSE) + #print(functionvector) + listsize <- nrow(functionvector) + ResultVector <- vector("list", listsize) + for (i in 1:listsize) + { + ResultVector[[i]] <- swe_fixstar2_mag(functionvector$starname[i]) + } + class(ResultVector)<-"swephR.star.mag" + return(ResultVector) + } + +#section6 +##' @rdname Vectorised +##' @export +vec_lun_eclipse_when_loc <- + function(jd_start, + ephe_flag = 4, + # default Moshier epheemris + longitude, + lat, + height = 0, + backward = FALSE) { + #default is forward search + functionvector <- + data.frame(jd_start, ephe_flag, longitude, lat, height, backward) + #print(functionvector) + listsize <- nrow(functionvector) + ResultVector <- vector("list", listsize) + for (i in 1:listsize) + { + geopos <- + c(functionvector$longitude[i], + functionvector$lat[i], + functionvector$height[i]) + ResultVector[[i]] <- swe_lun_eclipse_when_loc( + functionvector$jd_start[i], + functionvector$ephe_flag[i], + geopos, + functionvector$backward[i] + ) + } + class(ResultVector)<-"swephR.eclipse.loc" + return(ResultVector) + } + +##' @rdname Vectorised +##' @export +vec_lun_eclipse_how <- + function(jd_ut, + ephe_flag = 4, + # default Moshier epheemris + longitude, + lat, + height = 0) { + functionvector <- + data.frame(jd_ut, ephe_flag, longitude, lat, height) + #print(functionvector) + listsize <- nrow(functionvector) + ResultVector <- vector("list", listsize) + for (i in 1:listsize) + { + geopos <- + c(functionvector$longitude[i], + functionvector$lat[i], + functionvector$height[i]) + ResultVector[[i]] <- + swe_lun_eclipse_how(functionvector$jd_ut[i], + functionvector$ephe_flag[i], + geopos) + } + class(ResultVector)<-"swephR.eclipse.how" + return(ResultVector) + } + +##' @rdname Vectorised +##' @export +vec_lun_eclipse_when <- + function(jd_start, + ephe_flag = 4, + # default Moshier epheemris + ifltype, + backward = FALSE) { + #default is forward search + functionvector <- + data.frame(jd_start, ephe_flag, ifltype, backward) + #print(functionvector) + listsize <- nrow(functionvector) + ResultVector <- vector("list", listsize) + for (i in 1:listsize) + { + ResultVector[[i]] <- swe_lun_eclipse_when( + functionvector$jd_start[i], + functionvector$ephe_flag[i], + functionvector$ifltype[i], + functionvector$backward[i] + ) + } + class(ResultVector)<-"swephR.eclipse.when" + return(ResultVector) + } + +##' @rdname Vectorised +##' @export +vec_rise_trans_true_hor <- + function(jd_ut, + ipl, + starname = "", + ephe_flag = 4, + # default Moshier epheemris + rsmi, + longitude, + lat, + height = 0, + atpress = 1013.25, + attemp = 15, + horhgt = 0) { + functionvector <- + data.frame( + jd_ut, + ipl, + starname, + ephe_flag, + rsmi, + longitude, + lat, + height, + atpress, + attemp, + horhgt, + stringsAsFactors = FALSE + ) + #print(functionvector) + listsize <- nrow(functionvector) + ResultVector <- vector("list", listsize) + for (i in 1:listsize) + { + geopos <- + c(functionvector$longitude[i], + functionvector$lat[i], + functionvector$height[i]) + ResultVector[[i]] <- swe_rise_trans_true_hor( + functionvector$jd_ut[i], + functionvector$ipl[i], + functionvector$starname[i], + functionvector$ephe_flag[i], + functionvector$rsmi[i], + geopos, + functionvector$atpress[i], + functionvector$attemp[i], + functionvector$horhgt[i] + ) + } + class(ResultVector)<-"swephR.object.event" + return(ResultVector) + } + +##' @rdname Vectorised +##' @export +vec_pheno_ut <- + function(jd_ut, + ipl, + iflag = 4) { + #default Moshier epheemris + functionvector <- + data.frame(jd_ut, ipl, iflag) + #print(functionvector) + listsize <- nrow(functionvector) + ResultVector <- vector("list", listsize) + for (i in 1:listsize) + { + ResultVector[[i]] <- swe_pheno_ut(functionvector$jd_ut[i], + functionvector$ipl[i], + functionvector$iflag[i]) + } + class(ResultVector)<-"swephR.object.pheno" + return(ResultVector) + } + +##' @rdname Vectorised +##' @export +vec_pheno <- + function(jd_et, + ipl, + iflag = 4) { + #default Moshier epheemris + functionvector <- + data.frame(jd_et, ipl, iflag) + #print(functionvector) + listsize <- nrow(functionvector) + ResultVector <- vector("list", listsize) + for (i in 1:listsize) + { + ResultVector[[i]] <- swe_pheno(functionvector$jd_et[i], + functionvector$ipl[i], + functionvector$iflag[i]) + } + class(ResultVector)<-"swephR.object.pheno" + return(ResultVector) + } + +##' @rdname Vectorised +##' @export +vec_azalt <- + function(jd_ut, + coord_flag, + longitude, + lat, + height = 0, + atpress = 1013.25, + attemp = 15, + xin1, + xin2) { + functionvector <- + data.frame(jd_ut, + coord_flag, + longitude, + lat, + height, + atpress, + attemp, + xin1, + xin2) + #print(functionvector) + listsize <- nrow(functionvector) + ResultVector <- vector("list", listsize) + for (i in 1:listsize) + { + geopos <- + c(functionvector$longitude[i], + functionvector$lat[i], + functionvector$height[i]) + xin <- + c(functionvector$xin1[i], + functionvector$xin2[i]) + ResultVector[[i]] <- swe_azalt( + functionvector$jd_ut[i], + functionvector$coord_flag[i], + geopos, + functionvector$atpress[i], + functionvector$attemp[i], + xin + ) + } + class(ResultVector)<-"swephR.coord" + return(ResultVector) + } + +##' @rdname Vectorised +##' @export +vec_azalt_rev <- + function(jd_ut, + coord_flag, + longitude, + lat, + height = 0, + xin1, + xin2) { + functionvector <- + data.frame(jd_ut, + coord_flag, + longitude, + lat, + height, + xin1, + xin2) + #print(functionvector) + listsize <- nrow(functionvector) + ResultVector <- vector("list", listsize) + for (i in 1:listsize) + { + geopos <- + c(functionvector$longitude[i], + functionvector$lat[i], + functionvector$height[i]) + xin <- + c(functionvector$xin1[i], + functionvector$xin2[i]) + ResultVector[[i]] <- swe_azalt_rev(functionvector$jd_ut[i], + functionvector$coord_flag[i], + geopos, + xin) + } + class(ResultVector)<-"swephR.coord" + return(ResultVector) + } + +##' @rdname Vectorised +##' @export +vec_refrac_extended <- + function(InAlt, + height = 0, + atpress = 1013.25, + attemp = 15, + lapse_rate = 0.0065, + calc_flag) { + functionvector <- + data.frame(InAlt, + height, + atpress, + attemp, + lapse_rate, + calc_flag) + #print(functionvector) + listsize <- nrow(functionvector) + ResultVector <- vector("list", listsize) + for (i in 1:listsize) + { + ResultVector[[i]] <- swe_refrac_extended( + functionvector$InAlt[i], + functionvector$height[i], + functionvector$atpress[i], + functionvector$attemp[i], + functionvector$lapse_rate[i], + functionvector$calc_flag[i] + ) + } + class(ResultVector)<-"swephR.refract" + return(ResultVector) + } + +##' @rdname Vectorised +##' @export +vec_heliacal_ut <- + function(jd_utstart, + longitude, + lat, + height = 0, + atpress = 1013.25, + attemp = 15, + athum = 0, + atvis = 0.2, + obsage = 36, + obssnellen = 1.4, + # Default: average acuity + obsbin = 1, + obsmag = 1, + obsaper = 7, + obstrans = 0.8, + objectname, + event_type, + helflag = 260) { + # Default high precision and Moshier ephemeris + functionvector <- + data.frame( + jd_utstart, + longitude, + lat, + height, + atpress , + attemp , + athum, + atvis, + obsage, + obssnellen, + obsbin, + obsmag, + obsaper, + obstrans, + objectname, + event_type, + helflag, + stringsAsFactors = FALSE + ) + #print(functionvector) + listsize <- nrow(functionvector) + ResultVector <- vector("list", listsize) + for (i in 1:listsize) + { + geopos <- + c(functionvector$longitude[i], + functionvector$lat[i], + functionvector$height[i]) + datm <- + c( + functionvector$atpress[i], + functionvector$attemp[i], + functionvector$athum[i], + functionvector$atvis[i] + ) + dobs <- + c( + functionvector$obsage[i], + functionvector$obssnellen[i], + functionvector$obsbin[i], + functionvector$obsmag[i], + functionvector$obsaper[i], + functionvector$obstrans[i] + ) + ResultVector[[i]] <- swe_heliacal_ut( + functionvector$jd_utstart[i], + geopos, + datm, + dobs, + functionvector$objectname[i], + functionvector$event_type[i], + functionvector$helflag[i] + ) + } + class(ResultVector)<-"swephR.heliacal.event" + return(ResultVector) + } + +##' @rdname Vectorised +##' @export +vec_vis_limit_mag <- + function(jd_ut, + longitude, + lat, + height = 0, + atpress = 1013.25, + attemp = 15, + athum = 0, + atvis = 0.2, + obsage = 36, + obssnellen = 1.4, + # Default: average acuity + obsbin = 1, + obsmag = 1, + obsaper = 7, + obstrans = 0.8, + objectname, + helflag = 260) { + # Default high precision and Moshier ephemeris + functionvector <- + data.frame( + jd_ut, + longitude, + lat, + height, + atpress , + attemp , + athum, + atvis, + obsage, + obssnellen, + obsbin, + obsmag, + obsaper, + obstrans, + objectname, + helflag, + stringsAsFactors = FALSE + ) + #print(functionvector) + listsize <- nrow(functionvector) + ResultVector <- vector("list", listsize) + for (i in 1:listsize) + { + geopos <- + c(functionvector$longitude[i], + functionvector$lat[i], + functionvector$height[i]) + datm <- + c( + functionvector$atpress[i], + functionvector$attemp[i], + functionvector$athum[i], + functionvector$atvis[i] + ) + dobs <- + c( + functionvector$obsage[i], + functionvector$obssnellen[i], + functionvector$obsbin[i], + functionvector$obsmag[i], + functionvector$obsaper[i], + functionvector$obstrans[i] + ) + ResultVector[[i]] <- swe_vis_limit_mag( + functionvector$jd_ut[i], + geopos, + datm, + dobs, + functionvector$objectname[i], + functionvector$helflag[i] + ) + } + class(ResultVector)<-"swephR.heliacal.mag" + return(ResultVector) + } + +##' @rdname Vectorised +##' @export +vec_heliacal_pheno_ut <- + function(jd_ut, + longitude, + lat, + height = 0, + atpress = 1013.25, + attemp = 15, + athum = 0, + atvis = 0.2, + obsage = 36, + obssnellen = 1.4, + # Default: average acuity + obsbin = 1, + obsmag = 1, + obsaper = 7, + obstrans = 0.8, + objectname, + event_type, + helflag = 260) { + # Default high precision and Moshier ephemeris + functionvector <- + data.frame( + jd_ut, + longitude, + lat, + height, + atpress , + attemp , + athum, + atvis, + obsage, + obssnellen, + obsbin, + obsmag, + obsaper, + obstrans, + objectname, + event_type, + helflag, + stringsAsFactors = FALSE + ) + #print(functionvector) + listsize <- nrow(functionvector) + ResultVector <- vector("list", listsize) + for (i in 1:listsize) + { + geopos <- + c(functionvector$longitude[i], + functionvector$lat[i], + functionvector$height[i]) + datm <- + c( + functionvector$atpress[i], + functionvector$attemp[i], + functionvector$athum[i], + functionvector$atvis[i] + ) + dobs <- + c( + functionvector$obsage[i], + functionvector$obssnellen[i], + functionvector$obsbin[i], + functionvector$obsmag[i], + functionvector$obsaper[i], + functionvector$obstrans[i] + ) + ResultVector[[i]] <- swe_heliacal_pheno_ut( + functionvector$jd_ut[i], + geopos, + datm, + dobs, + functionvector$objectname[i], + functionvector$event_type[i], + functionvector$helflag[i] + ) + } + class(ResultVector)<-"swephR.heliacal.pheno" + return(ResultVector) + } + +##' @rdname Vectorised +##' @export +vec_topo_arcus_visionis <- + function(jd_ut, + longitude, + lat, + height = 0, + atpress = 1013.25, + attemp = 15, + athum = 0, + atvis = 0.2, + obsage = 36, + obssnellen = 1.4, + # Default: average acuity + obsbin = 1, + obsmag = 1, + obsaper = 7, + obstrans = 0.8, + helflag = 260, + # Default high precision and Moshier ephemeris + mag, + AziO, + AltO, + AziS, + AziM, + AltM) { + functionvector <- + data.frame( + jd_ut, + longitude, + lat, + height, + atpress , + attemp , + athum, + atvis, + obsage, + obssnellen, + obsbin, + obsmag, + obsaper, + obstrans, + helflag, + mag, + AziO, + AltO, + AziS, + AziM, + AltM + ) + #print(functionvector) + listsize <- nrow(functionvector) + ResultVector <- vector("list", listsize) + for (i in 1:listsize) + { + geopos <- + c(functionvector$longitude[i], + functionvector$lat[i], + functionvector$height[i]) + datm <- + c( + functionvector$atpress[i], + functionvector$attemp[i], + functionvector$athum[i], + functionvector$atvis[i] + ) + dobs <- + c( + functionvector$obsage[i], + functionvector$obssnellen[i], + functionvector$obsbin[i], + functionvector$obsmag[i], + functionvector$obsaper[i], + functionvector$obstrans[i] + ) + ResultVector[[i]] <- swe_topo_arcus_visionis( + functionvector$jd_ut[i], + geopos, + datm, + dobs, + functionvector$helflag[i], + functionvector$mag[i], + functionvector$AziO[i], + functionvector$AltO[i], + functionvector$AziS[i], + functionvector$AziM[i], + functionvector$AltM[i] + ) + } + class(ResultVector)<-"swephR.heliacal.av" + return(ResultVector) + } + +##' @rdname Vectorised +##' @export +vec_heliacal_angle <- + function(jd_ut, + longitude, + lat, + height = 0, + atpress = 1013.25, + attemp = 15, + athum = 0, + atvis = 0.2, + obsage = 36, + obssnellen = 1.4, + # Default: average acuity + obsbin = 1, + obsmag = 1, + obsaper = 7, + obstrans = 0.8, + helflag = 260, + # Default high precision and Moshier ephemeris + mag, + AziO, + AziS, + AziM, + AltM) { + functionvector <- + data.frame( + jd_ut, + longitude, + lat, + height, + atpress , + attemp , + athum, + atvis, + obsage, + obssnellen, + obsbin, + obsmag, + obsaper, + obstrans, + helflag, + mag, + AziO, + AziS, + AziM, + AltM + ) + #print(functionvector) + listsize <- nrow(functionvector) + ResultVector <- vector("list", listsize) + for (i in 1:listsize) + { + geopos <- + c(functionvector$longitude[i], + functionvector$lat[i], + functionvector$height[i]) + datm <- + c( + functionvector$atpress[i], + functionvector$attemp[i], + functionvector$athum[i], + functionvector$atvis[i] + ) + dobs <- + c( + functionvector$obsage[i], + functionvector$obssnellen[i], + functionvector$obsbin[i], + functionvector$obsmag[i], + functionvector$obsaper[i], + functionvector$obstrans[i] + ) + ResultVector[[i]] <- swe_heliacal_angle( + functionvector$jd_ut[i], + geopos, + datm, + dobs, + functionvector$helflag[i], + functionvector$mag[i], + functionvector$AziO[i], + functionvector$AziS[i], + functionvector$AziM[i], + functionvector$AltM[i] + ) + } + class(ResultVector)<-"swephR.heliacal.angle" + return(ResultVector) + } + +#section 7 +##' @rdname Vectorised +##' @export +vec_julday <- + function(year, + month, + day, + hour = 12, + #default is midday + gregflag = 1) { + #default (proleptic) Gregorian calendar + functionvector <- + data.frame(year, month, day, hour, gregflag) + #print(functionvector) + listsize <- nrow(functionvector) + ResultVector <- vector("list", listsize) + for (i in 1:listsize) + { + ResultVector[[i]] <- swe_julday( + functionvector$year[i], + functionvector$month[i], + functionvector$day[i], + functionvector$hour[i], + functionvector$gregflag[i] + ) + } + return(ResultVector) + } + +##' @rdname Vectorised +##' @export +vec_date_conversion <- + function(year, + month, + day, + hour = 12, + #default is midday + cal = "g") { + #default (proleptic) Gregorian calendar + functionvector <- + data.frame(year, month, day, hour, cal, + stringsAsFactors = FALSE) + #print(functionvector) + listsize <- nrow(functionvector) + ResultVector <- vector("list", listsize) + for (i in 1:listsize) + { + ResultVector[[i]] <- swe_date_conversion( + functionvector$year[i], + functionvector$month[i], + functionvector$day[i], + functionvector$hour[i], + functionvector$cal[i] + ) + } + class(ResultVector)<-"swephR.jdn" + return(ResultVector) + } + +##' @rdname Vectorised +##' @export +vec_revjul <- + function(jd, + gregflag = 1) { + #default (proleptic) Gregorian calendar + functionvector <- + data.frame(jd, gregflag) + #print(functionvector) + listsize <- nrow(functionvector) + ResultVector <- vector("list", listsize) + for (i in 1:listsize) + { + ResultVector[[i]] <- swe_revjul(functionvector$jd[i], + functionvector$gregflag[i]) + } + class(ResultVector)<-"swephR.calendar" + return(ResultVector) + } + +#section8 +##' @rdname Vectorised +##' @export +vec_deltat_ex <- + function(jd_ut, + ephe_flag = 4) { + # Default Moshier ephemeris + functionvector <- + data.frame(jd_ut, ephe_flag) + #print(functionvector) + listsize <- nrow(functionvector) + ResultVector <- vector("list", listsize) + for (i in 1:listsize) + { + ResultVector[[i]] <- swe_deltat_ex(functionvector$jd_ut[i], + functionvector$ephe_flag[i]) + } + return(ResultVector) + } + +##' @rdname Vectorised +##' @export +vec_deltat <- + function(jd_ut) { + functionvector <- + data.frame(jd_ut) + #print(functionvector) + listsize <- nrow(functionvector) + ResultVector <- vector("list", listsize) + for (i in 1:listsize) + { + ResultVector[[i]] <- swe_deltat(functionvector$jd_ut[i]) + } + return(ResultVector) + } + +#section16 +##' @rdname Vectorised +##' @export +vec_day_of_week <- + function(jd) { + functionvector <- + data.frame(jd) + #print(functionvector) + listsize <- nrow(functionvector) + ResultVector <- vector("list", listsize) + for (i in 1:listsize) + { + ResultVector[[i]] <- swe_day_of_week(functionvector$jd[i]) + } + return(ResultVector) + } diff --git a/man/Vectorised.Rd b/man/Vectorised.Rd new file mode 100644 index 00000000..a10f0e81 --- /dev/null +++ b/man/Vectorised.Rd @@ -0,0 +1,212 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/vectorisedSEfunctions.R +\name{vec_calc_ut} +\alias{vec_calc_ut} +\alias{vec_calc} +\alias{vec_get_planet_name} +\alias{vec_sol_eclipse_when_loc} +\alias{vec_fixstar2_ut} +\alias{vec_fixstar2} +\alias{vec_fixstar2_mag} +\alias{vec_lun_eclipse_when_loc} +\alias{vec_lun_eclipse_how} +\alias{vec_lun_eclipse_when} +\alias{vec_rise_trans_true_hor} +\alias{vec_pheno_ut} +\alias{vec_pheno} +\alias{vec_azalt} +\alias{vec_azalt_rev} +\alias{vec_refrac_extended} +\alias{vec_heliacal_ut} +\alias{vec_vis_limit_mag} +\alias{vec_heliacal_pheno_ut} +\alias{vec_topo_arcus_visionis} +\alias{vec_heliacal_angle} +\alias{vec_julday} +\alias{vec_date_conversion} +\alias{vec_revjul} +\alias{vec_deltat_ex} +\alias{vec_deltat} +\alias{vec_day_of_week} +\title{Vectorised swephR functions} +\usage{ +vec_calc_ut(jd_ut, ipl, iflag = 4) + +vec_calc(jd_et, ipl, iflag = 4) + +vec_get_planet_name(ipl) + +vec_sol_eclipse_when_loc(jd_start, ephe_flag = 4, longitude, lat, + height = 0, backward = FALSE) + +vec_fixstar2_ut(starname, jd_ut, iflag = 4) + +vec_fixstar2(starname, jd_et, iflag = 4) + +vec_fixstar2_mag(starname) + +vec_lun_eclipse_when_loc(jd_start, ephe_flag = 4, longitude, lat, + height = 0, backward = FALSE) + +vec_lun_eclipse_how(jd_ut, ephe_flag = 4, longitude, lat, height = 0) + +vec_lun_eclipse_when(jd_start, ephe_flag = 4, ifltype, + backward = FALSE) + +vec_rise_trans_true_hor(jd_ut, ipl, starname = "", ephe_flag = 4, rsmi, + longitude, lat, height = 0, atpress = 1013.25, attemp = 15, + horhgt = 0) + +vec_pheno_ut(jd_ut, ipl, iflag = 4) + +vec_pheno(jd_et, ipl, iflag = 4) + +vec_azalt(jd_ut, coord_flag, longitude, lat, height = 0, + atpress = 1013.25, attemp = 15, xin1, xin2) + +vec_azalt_rev(jd_ut, coord_flag, longitude, lat, height = 0, xin1, xin2) + +vec_refrac_extended(InAlt, height = 0, atpress = 1013.25, + attemp = 15, lapse_rate = 0.0065, calc_flag) + +vec_heliacal_ut(jd_utstart, longitude, lat, height = 0, + atpress = 1013.25, attemp = 15, athum = 0, atvis = 0.2, + obsage = 36, obssnellen = 1.4, obsbin = 1, obsmag = 1, + obsaper = 7, obstrans = 0.8, objectname, event_type, helflag = 260) + +vec_vis_limit_mag(jd_ut, longitude, lat, height = 0, atpress = 1013.25, + attemp = 15, athum = 0, atvis = 0.2, obsage = 36, + obssnellen = 1.4, obsbin = 1, obsmag = 1, obsaper = 7, + obstrans = 0.8, objectname, helflag = 260) + +vec_heliacal_pheno_ut(jd_ut, longitude, lat, height = 0, + atpress = 1013.25, attemp = 15, athum = 0, atvis = 0.2, + obsage = 36, obssnellen = 1.4, obsbin = 1, obsmag = 1, + obsaper = 7, obstrans = 0.8, objectname, event_type, helflag = 260) + +vec_topo_arcus_visionis(jd_ut, longitude, lat, height = 0, + atpress = 1013.25, attemp = 15, athum = 0, atvis = 0.2, + obsage = 36, obssnellen = 1.4, obsbin = 1, obsmag = 1, + obsaper = 7, obstrans = 0.8, helflag = 260, mag, AziO, AltO, AziS, + AziM, AltM) + +vec_heliacal_angle(jd_ut, longitude, lat, height = 0, + atpress = 1013.25, attemp = 15, athum = 0, atvis = 0.2, + obsage = 36, obssnellen = 1.4, obsbin = 1, obsmag = 1, + obsaper = 7, obstrans = 0.8, helflag = 260, mag, AziO, AziS, AziM, + AltM) + +vec_julday(year, month, day, hour = 12, gregflag = 1) + +vec_date_conversion(year, month, day, hour = 12, cal = "g") + +vec_revjul(jd, gregflag = 1) + +vec_deltat_ex(jd_ut, ephe_flag = 4) + +vec_deltat(jd_ut) + +vec_day_of_week(jd) +} +\arguments{ +\item{jd_ut}{UT Julian day number vector (day)} + +\item{ipl}{Body/planet as integer vector (SE$SUN=0, SE$Moon=1, ... SE$PLUTO=9)} + +\item{iflag}{Computation flag as integer vector, many options possible (section 2.3)} + +\item{jd_et}{ET Julian day number as double vector (day)} + +\item{jd_start}{Julian day number as double vector (day)} + +\item{ephe_flag}{Ephemeris flag as integer vector (SE$FLG_JPLEPH=1, SE$FLG_SWIEPH=2 or SE$FLG_MOSEPH=4)} + +\item{longitude}{Geographic longitude as double vector (deg)} + +\item{lat}{Geographic latitude as double vector (deg)} + +\item{height}{Height as double vector (m)} + +\item{backward}{backwards search as boolean vector (TRUE)} + +\item{starname}{Star name as string vector ("" for no star)} + +\item{ifltype}{eclipse type as integer vector (e.g.: SE$ECL_CENTRAL=1,SE$ECL_NONCENTRAL=2,SE$ECL_TOTAL=4,SE$ECL_ANNULAR=8,SE$ECL_PARTIAL=16,SE$ECL_ANNULAR_TOTAL=32)} + +\item{rsmi}{Event flag as integer vector (e.g.: SE$CALC_RISE=1, SE$CALC_SET=2,SE$CALC_MTRANSIT=4,SE$CALC_ITRANSIT=8)} + +\item{atpress}{Atmospheric pressure as double vector (hPa)} + +\item{attemp}{Atmospheric temperature as double vector (Celsius)} + +\item{horhgt}{Horizon apparent altitude as double vector (deg)} + +\item{coord_flag}{Coordinate flag as integer vector (reference system (SE$ECL2HOR=0 or SE$EQU2HOR=1))} + +\item{xin1}{Position of body as numeric vector (either ecliptical or equatorial coordinates, depending on coord_flag)} + +\item{xin2}{Position of body as numeric vector (either ecliptical or equatorial coordinates, depending on coord_flag)} + +\item{InAlt}{object's apparent/topocentric altitude as double vector (depending on calc_flag) (deg)} + +\item{lapse_rate}{lapse rate as double vector (K/m)} + +\item{calc_flag}{Calculation flag as integer vector (refraction direction (SE$TRUE_TO_APP=0 or SE$APP_TO_TRUE=1))} + +\item{jd_utstart}{UT Julian day number as double vector (day)} + +\item{athum}{Atmospheric humidity as double vector (\%)} + +\item{atvis}{Atmospheric visibiliy as double vector (km or -)} + +\item{obsage}{Age of observer as double vector (year)} + +\item{obssnellen}{Aquity of observer as double vector (-)} + +\item{obsbin}{Mono-/bi-nocular observation as double vector (-)} + +\item{obsmag}{Magnification as double vector (-)} + +\item{obsaper}{Aperture of optics as double vector (mm)} + +\item{obstrans}{Transmission as double vector (-)} + +\item{objectname}{Name of fixed star or planet as string vector} + +\item{event_type}{Event type as integer vector} + +\item{helflag}{Calculation flag (incl. ephe_flag values) as integer vector} + +\item{mag}{Object's visible magnitude (Vmag) as double vector (-)} + +\item{AziO}{Object's azimuth as double vector (deg)} + +\item{AltO}{Object's altitude as double vecor (deg)} + +\item{AziS}{Sun's azimuth as double (vector deg)} + +\item{AziM}{Moon's azimut as double vecor (deg)} + +\item{AltM}{Moon's altitude as double vector (deg)} + +\item{year}{Astronomical year as integer vector} + +\item{month}{Month as integer vector} + +\item{day}{Day as integer vector} + +\item{hour}{Hour as double vector} + +\item{gregflag}{Calendar type as integer vector (SE$JUL_CAL=0 or SE$GREG_CAL=1)} + +\item{cal}{Calendar type as char vector ("g"[regorian] or "j"[ulian])} + +\item{jd}{Julian day number as double vector (day)} +} +\description{ +This are the vectorised functions of the +ones that are namend swe_xxx. Each input parameter can be a vector. +Make sure that the vector length of each input parameter matches +the rules around data-frames: so the largest vector length must be an integer +multiple of the other vector lengths. +} diff --git a/vignettes/VectorisedSwephR.Rmd b/vignettes/VectorisedSwephR.Rmd new file mode 100644 index 00000000..a061a3a2 --- /dev/null +++ b/vignettes/VectorisedSwephR.Rmd @@ -0,0 +1,76 @@ +--- +title: "swephR" +author: "Victor Reijs" +date: "`r Sys.Date()`" +output: rmarkdown::html_vignette +vignette: > + %\VignetteIndexEntry{swephR} + %\VignetteEngine{knitr::rmarkdown} + %\VignetteEncoding{UTF-8} +--- + +```{r setup, include = FALSE} +knitr::opts_chunk$set( + collapse = TRUE, + comment = "#>" +) +``` + +# Introduction +This vignette provides input how the vectorised swephR functions can be used (these function names start with 'vec-').The vectorised functions provide a list of lists as results. +These functions are not finalised and are in testing. A decision STILL has to be made if the reuslts should be a list of matrices or a list of lists. + + +# Some simple steps to do calculations +To compute the position of celestial body or star with SE (Swiss Ephemeris) for two dates (aka using vectorisation), you do the following steps: + + +1. First load `swephR` and constants used in swephR/SE: +```{r} +library(swephR) +data(SE) +``` + +2. Optionaly set the directory path of the ephemeris files (the format of the path depends on your OS), e.g.: +```{r, eval = FALSE} +swe_set_ephe_path("C:\\sweph\\ephe") +``` + +3. For a specific date, compute the Julian day number (in below example: 1 January 2000CE at 12:00 UT (or J2000.0) and 1 January 1BCE at 12:00): +```{r} +year1 <- 2000 +year2 <- 0 +month <- 1 +day <- 1 +hour <- 12 +jdut <- vec_julday(c(year1,year2), month, day, hour, SE$GREG_CAL) +jdut +``` +4. Compute (using Moshier ephemeris) the positions (longitude, latitude, distance, longitude speed and latitude speed) of a planet or other celestial bodies for two dates (in below example: the Sun and Moon): +```{r} +ipl1 <- SE$SUN +ipl2 <- SE$MOON +iflag <- SE$FLG_MOSEPH + SE$FLG_SPEED +result <- vec_calc_ut(c(jdut[[1]],jdut[[2]]), c(ipl1,ipl1,ipl2,ipl2), iflag) +result +``` +or a fixed star on two dates (in below example: Sirius): +```{r} +starname = "sirius" +result <- vec_fixstar2_ut(starname, c(jdut[[1]],jdut[[2]]), iflag) +result +``` + +5. Determine the Julian day number of the Heliacal Rise of Sirius for two epochs: +```{r} +options(digits=15) +result <- vec_heliacal_ut(c(jdut[[1]],jdut[[2]]),0,50,10,1013.25,15,50,0.25,25,1,1,1,5,0.8,starname, + SE$HELIACAL_RISING,SE$HELFLAG_HIGH_PRECISION+SE$FLG_MOSEPH) +result +``` + + +6. At the end of your computations close all files and free up memory: +```{r} +swe_close() +```