From 833ef21c01f7c72a54097e867445c3d21466d048 Mon Sep 17 00:00:00 2001 From: Sven Gastauer Date: Wed, 11 Nov 2020 08:23:48 +0100 Subject: [PATCH] ts.cal Added ts.cal - TS for calibraiton sphere --- NAMESPACE | 1 + R/TScal.R | 112 +++++++++++++++++++++++++++++++++++++++++++++++ man/c_Leroy08.Rd | 8 ++-- man/ts.cal.Rd | 50 +++++++++++++++++++++ 4 files changed, 167 insertions(+), 4 deletions(-) create mode 100644 R/TScal.R create mode 100644 man/ts.cal.Rd diff --git a/NAMESPACE b/NAMESPACE index a29a530..7a50df1 100755 --- a/NAMESPACE +++ b/NAMESPACE @@ -34,6 +34,7 @@ export(rho_c) export(rho_p0) export(rho_smow) export(smoother) +export(ts.cal) import(ggplot2) import(pracma) import(reshape2) diff --git a/R/TScal.R b/R/TScal.R new file mode 100644 index 0000000..a6a0888 --- /dev/null +++ b/R/TScal.R @@ -0,0 +1,112 @@ +#' Theoretical TS of calibration sphere +#' R script of function for band averaged TS of tungsten +#' carbide spheres translated from Matlab code (Version date 31/5/05) +#' written by D. MacLennan (Marine Lab, Aberdeen) +#' updated from R version written by Sascha F?ssler 03/2010 +#' updated Sven Gastauer 04/01/2017 +#' @import ggplot2 +#' @export +#' @param freq numeric frequencies [kHz]; default seq(38,200) +#' @param c numeric ambient soundspeeds [m/s]; default 1500 +#' @param d numeric sphere diameter [cm]; default = 38.1 +#' @param mat string material properties of the aphere Cu = Copper, TC = Tungsten-Carbide; default 'TC' +#' @param water string either sw= seawater or fw = fresh water; will be ignored if rhow (density of seawater is defined); default 'sw +#' @param rhow numeric density of ambient seawater, if NULL parameter water will be used; default NULL +#' @param plot boolean TRUE/FALSE if plot should be part of the output; default TRUE +#' @return dataframe with columns "F","TS","ModF2","c", where F = Frequency [kHz]; TS = TS [dB]; ModF2 = ModF^2 and c soundspeed [m/s], an optional ggplot can be part of the output +#' @examples ts.cal(freq=seq(90,170,by=0.1),c=1500.5,rhow=1.02509,plot=TRUE) +ts.cal <- function(freq=seq(38,200),c=1500,d=38.1,mat="TC",water="sw",rhow=NULL,plot=TRUE){ + require(ggplot2) + a <- d/2 + results <-NA + #select materia l + cc <- switch(mat, + TC = c(6853,4171), + Cu = c(4760,2288)) + if(is.null(rhow)){ + ww <- switch(water, + sw = 1.027, + fw = 1)} + else{ + ww = rhow + } + rho <- switch(mat, + TC = c(14.9/ww), + Cu = c(8.945/ww)) + + for(i in 1:length(c)){ + q<- 2*pi*freq*a/c[i] # ka range + ka=q + + nr <- length(ka) + F <- matrix(0,nr,4) + Lmax <- floor(max(ka))+20 + alpha <- 2*rho*(cc[2]/c[i])^2 + beta <- rho*(cc[1]/c[i])^2-alpha + nn <- 1:(Lmax+1) + n0 <- 2:(Lmax+1) + lh <- nn-0.5 + LL <- 0:Lmax + S <- (floor(LL/2)==LL/2)*1-(floor((LL+1)/2)==(LL+1)/2)*1 + L <- 1:Lmax + L <- matrix(c(L,L,L),length(L),3) + jh <- matrix(0,Lmax+1,3) + djh <- jh + ddjh <- jh + yh <- matrix(0,Lmax+1,1) + + for (jj in 1:nr){ + q <- ka[jj] + q1 <- q*c[i]/cc[1] + q2 <- q*c[i]/cc[2] + qq <- matrix(c(q1,q2,q),length(q),3) + bfac <- sqrt((qq*2/pi)); + Qm <- t(matrix(qq,length(qq),Lmax)) + Bf <- t(matrix(bfac,length(bfac),Lmax+1)) + + # jh/yh are spherical Bessel functions of the first/second kinds + # starting at order L=0 (i.e. order is row number - 1). + jh <- (matrix(c(besselJ(q1,lh),besselJ(q2,lh),besselJ(q,lh)),length(nn),3))/Bf + yh <- besselY(q,lh)/Bf[,3] + djh[1,] <- -jh[1,]+cos(qq) # Zero order derivatives + dyh <- -yh[1]+sin(q) # djh= x*jL'(x); dyh= x*yL'(x) + ddjh[1,] <- -2*djh[1,]-qq*sin(qq) # ddjh = x^2*jL''(x) + djh[n0,] <- -(L+1)*jh[n0,]+Qm*jh[n0-1,] # Derivatives for orders 2 upwards + dyh[n0] <- -(L[,1]+1)*yh[n0]+Qm[,3]*yh[n0-1] + ddjh[n0,] <- ((L+1)*(L+2)-Qm*Qm)*jh[n0,]-2*Qm*jh[n0-1,] + + a2 <- (LL*LL+LL-2)*jh[,2]+ddjh[,2] + a1 <- 2*LL*(LL+1)*(djh[,1]-jh[,1]) + b2 <- a2*(beta*q1^2*jh[,1]-alpha*ddjh[,1])-alpha*a1*(jh[,2]-djh[,2]) + b1 <- q*(a2*djh[,1]-a1*jh[,2]) + x <- -b2*djh[,3]/q + b1*jh[,3] + y <- b2*dyh/q - b1*yh + z <- S*(2*LL+1)*x/(x*x+y*y) + y <- z*y + x <- z*x + + F[jj,3:4] <- -(2/q)*c(sum(y),sum(x)) # The form function + F[jj,1] <- F[jj,3]^2+F[jj,4]^2 # and its modulus squared + F[jj,2]<-10*log10((a/2000)^2*F[jj,1]) + } + tmp <- as.data.frame(cbind(rbind(cbind(freq, F[,2], F[,1])),c[i])) + names(tmp)<-c("F","TS","ModF2","c") + if(i>1 | jj>1){ + results <- rbind(results,tmp) + }else{ + results<-tmp + } + } + if(plot == TRUE){ + #plot(results$TS,col=results$c) + pp<-ggplot(data=results, aes(x=F,y=TS,group=c))+ + geom_line(aes(lty=as.factor(c)))+ + labs(lty = "c [m/s]")+ + theme_classic()+ + theme(legend.position = "top")+ + xlab("Frequency [kHz]")+ylab("TS [dB]") + print(pp) + } + results <- na.omit(results) + return(results) +} diff --git a/man/c_Leroy08.Rd b/man/c_Leroy08.Rd index 825403f..bf7fef6 100755 --- a/man/c_Leroy08.Rd +++ b/man/c_Leroy08.Rd @@ -25,14 +25,14 @@ Returns the sound speed according to Leroy et al (2008). This "newer" equation s } \examples{ # TABLE III in Leroy et al. (2008) -# Common oceans, lat = 30°, P= 80 MPa Z= 7808.13 m, S= 34.7\% +# Common oceans, lat = 30B0, P= 80 MPa Z= 7808.13 m, S= 34.7\% lat=30; Z=7808.13; S=34.7; T=c(1,1.5,2,2.5,3) c_Leroy08(Z,T,S,lat) -# Common oceans, lat = 30°, P= 80 MPa Z= 7808.13 m, T=2 °C +# Common oceans, lat = 30B0, P= 80 MPa Z= 7808.13 m, T=2 B0C c_Leroy08(Z,T=2,S=seq(33.5,35.5,.5),lat) -# Common oceans, = 30°, P= 5 MPa Z= 497.12 m, S= 35\% +# Common oceans, = 30B0, P= 5 MPa Z= 497.12 m, S= 35\% c_Leroy08(Z=497.12,T=seq(-2,20,2),S=35,lat) -# Common oceans, = 30°, P= 5 MPa Z= 497.12 m, T=8 °C +# Common oceans, = 30B0, P= 5 MPa Z= 497.12 m, T=8 B0C c_Leroy08(Z=497.12,T=8,S=seq(33,37,1),lat) } \references{ diff --git a/man/ts.cal.Rd b/man/ts.cal.Rd new file mode 100644 index 0000000..268247a --- /dev/null +++ b/man/ts.cal.Rd @@ -0,0 +1,50 @@ +% Generated by roxygen2: do not edit by hand +% Please edit documentation in R/TScal.R +\name{ts.cal} +\alias{ts.cal} +\title{Theoretical TS of calibration sphere +R script of function for band averaged TS of tungsten +carbide spheres translated from Matlab code (Version date 31/5/05) +written by D. MacLennan (Marine Lab, Aberdeen) +updated from R version written by Sascha F?ssler 03/2010 +updated Sven Gastauer 04/01/2017} +\usage{ +ts.cal( + freq = seq(38, 200), + c = 1500, + d = 38.1, + mat = "TC", + water = "sw", + rhow = NULL, + plot = TRUE +) +} +\arguments{ +\item{freq}{numeric frequencies [kHz]; default seq(38,200)} + +\item{c}{numeric ambient soundspeeds [m/s]; default 1500} + +\item{d}{numeric sphere diameter [cm]; default = 38.1} + +\item{mat}{string material properties of the aphere Cu = Copper, TC = Tungsten-Carbide; default 'TC'} + +\item{water}{string either sw= seawater or fw = fresh water; will be ignored if rhow (density of seawater is defined); default 'sw} + +\item{rhow}{numeric density of ambient seawater, if NULL parameter water will be used; default NULL} + +\item{plot}{boolean TRUE/FALSE if plot should be part of the output; default TRUE} +} +\value{ +dataframe with columns "F","TS","ModF2","c", where F = Frequency [kHz]; TS = TS [dB]; ModF2 = ModF^2 and c soundspeed [m/s], an optional ggplot can be part of the output +} +\description{ +Theoretical TS of calibration sphere +R script of function for band averaged TS of tungsten +carbide spheres translated from Matlab code (Version date 31/5/05) +written by D. MacLennan (Marine Lab, Aberdeen) +updated from R version written by Sascha F?ssler 03/2010 +updated Sven Gastauer 04/01/2017 +} +\examples{ +ts.cal(freq=seq(90,170,by=0.1),c=1500.5,rhow=1.02509,plot=TRUE) +}