diff --git a/R/scattnlay.r b/R/scattnlay.r index 53e8bf7..e36c3a9 100644 --- a/R/scattnlay.r +++ b/R/scattnlay.r @@ -1,7 +1,7 @@ Layer <- setClass("Layer", representation( m="complex", - d="numeric" + r="numeric" )) Scatterer <- setClass("Scatterer", @@ -13,7 +13,7 @@ Scatterer <- setClass("Scatterer", setMethod("show", signature(object="Layer"),function(object){ cat("Layer Data \n") - cat("Diameter:", object@d," nm\n") + cat("Diameter:", 2*object@r," nm\n") cat("Epsilon:", object@m,"\n") }) @@ -24,12 +24,12 @@ setMethod("show", signature(object="Scatterer"),function(object){ setGeneric("na<-",function(x,value) standardGeneric("na<-")) setGeneric("lambda<-",function(x,value) standardGeneric("lambda<-")) -setGeneric("d<-",function(x,value) standardGeneric("d<-")) +setGeneric("r<-",function(x,value) standardGeneric("r<-")) setGeneric("m<-",function(x,value) standardGeneric("m<-")) setReplaceMethod("na","Scatterer", function(x,value) {x@na <- value; validObject(x); x}) setReplaceMethod("lambda","Scatterer", function(x,value) {x@lambda <- value; validObject(x); x}) -setReplaceMethod("d","Layer", function(x,value) {x@d <- value; validObject(x); x}) +setReplaceMethod("r","Layer", function(x,value) {x@r <- value; validObject(x); x}) setReplaceMethod("m","Layer", function(x,value) {x@m <- value; validObject(x); x}) setGeneric("scattnlay", function(object) {standardGeneric("scattnlay")}) diff --git a/oil_scat.R b/oil_scat.R new file mode 100644 index 0000000..4d69ee2 --- /dev/null +++ b/oil_scat.R @@ -0,0 +1,27 @@ +library(Rscattnlay) +library(ggplot2) +library(reshape2) + +S <- Scatterer() # create a scatterer object + +dr <- Layer() # and a layer representing the droplet +na(S) <- 1.33 # set the ambient index +lambda(S) <- 600 + +m(dr) <- 1.431+0i; + +scatnlay_ar <- function(r){ + r(dr) <- r + St <- S+dr; # PUT THE STACK HERE + Q <- scattnlay(St); +} + +dr_size <- c(1:100) + +output <- sapply(dr_size,scatnlay_ar) + +data_mat <- data.frame(t(output),dr_size,'oil') +colnames(data_mat) <- c("Qext", "Qsca", "Qabs", "Qbk", "Qpr", "g", "Albedo", "nmax","Size","Layers") + +scat_plot <- ggplot(data=data_mat, aes(x=Size,y=Qsca)) +scat_plot + geom_line() + theme_minimal(base_size=22) + xlim(c(0,100)) + xlab("Radius (nm)") \ No newline at end of file diff --git a/test_2.R b/test_2.R new file mode 100644 index 0000000..507999c --- /dev/null +++ b/test_2.R @@ -0,0 +1,81 @@ +library(Rscattnlay) +library(ggplot2) +library(reshape2) + +S <- Scatterer() # create a scatterer object + +np <- Layer() # and a layer representing the np +lipid <- Layer() # and other layers +inner_layer <- Layer() +water <- Layer() + +na(S) <- 1.33 # set the ambient index +r(np) <- 35 # and the size (radius) of the np in nm + +r(water) <- 40 +m(water) <- 1.33+0i + +r(inner_layer) <- 38 +m(inner_layer) <- 1.33+0i + +r(lipid) <- 40 # and the size (radius) of the other layers +m(lipid) <- 1.45+0i #fixed value + +# load up silver data +palik_ag_vis <- read.table("palik_ag_vis_hb.csv", header=TRUE) +colnames(palik_ag_vis) <- c("lambda","n","k","eps_real","eps_imag") + +lambda_palik <- palik_ag_vis$lambda #convert to nm +n_palik <- palik_ag_vis$n +k_palik <- palik_ag_vis$k + +#interpolation as palik data is sparse +n = 512 #number of points +spl_n <- approx(lambda_palik,n_palik,n=n) +spl_k <- approx(lambda_palik,k_palik,n=n) +lambda = spl_n$x + +scatnlay_ar <-function(k){ + m(np) <- spl_n$y[k]+spl_k$y[k]*(0+1i); + lambda(S) <- lambda[k]; + St <- S+np+water; # PUT THE STACK HERE + Q <- scattnlay(St); +} + +scatnlay_ar_water_lipid <-function(k){ + m(np) <- spl_n$y[k]+spl_k$y[k]*(0+1i); + lambda(S) <- lambda[k]; + St <- S+np+inner_layer+lipid; # PUT THE STACK HERE + Q <- scattnlay(St); +} + +output <- sapply(1:n,scatnlay_ar) +output_lipid <- sapply(1:n,scatnlay_ar_water_lipid) +data_mat_np <- data.frame(t(output),lambda,'water') +data_mat_lipid <- data.frame(t(output_lipid),lambda,'lipid') + +colnames(data_mat_np) <- c("Qext", "Qsca", "Qabs", "Qbk", "Qpr", "g", "Albedo", "nmax","Lambda","Layers") +colnames(data_mat_lipid) <- c("Qext", "Qsca", "Qabs", "Qbk", "Qpr", "g", "Albedo", "nmax","Lambda","Layers") + +data_mat <- rbind(data_mat_np,data_mat_lipid) + +ext_water <- subset(data_mat_np, Lambda < 700, select= c("Qext","Lambda")) +ext_water$Lambda[which.max(ext_water$Qext)] + +ext_lipid <- subset(data_mat_lipid, Lambda < 700, select= c("Qext","Lambda")) +ext_lipid$Lambda[which.max(ext_lipid$Qext)] + +dif_Qext <- data.frame(ext_lipid$Qext-ext_water$Qext) ## might be able to use aggregate here +dif_spectra <- cbind(ext_lipid$Lambda,dif_Qext) +colnames(dif_spectra) <- c("lambda","D_Qext") + +palik_ag_melt <- melt(subset(palik_ag_vis,select = c("lambda","n","k")), id=c("lambda"),variable.name = "type", + value.name = "Index") +ext_plot <- ggplot(data=data_mat, aes(x=Lambda,y=Qext, group=as.factor(Layers), color=as.factor(Layers))) +ext_plot + geom_line() + theme_minimal(base_size=22) + xlim(c(350,550))+ ylim(c(0,10)) + labs(colour = "Experiment") + xlab("Wavelength (nm)") + +scat_plot <- ggplot(data=data_mat, aes(x=Lambda,y=Qsca, group=as.factor(Layers), color=as.factor(Layers))) +scat_plot + geom_line() + theme_minimal(base_size=22) + xlim(c(350,550))+ ylim(c(0,7.5)) + labs(colour = "Experiment") + xlab("Wavelength (nm)") + +dif_plot <- ggplot(data=dif_spectra, aes(x=lambda,y=D_Qext)) +dif_plot + geom_line() + theme_minimal(base_size=22) + xlim(c(350,550)) + xlab("Wavelength (nm)") diff --git a/test_3.R b/test_3.R new file mode 100644 index 0000000..175449c --- /dev/null +++ b/test_3.R @@ -0,0 +1,53 @@ +#scattnlay() +library(Rscattnlay) +library(ggplot2) +library(reshape2) + +S <- Scatterer() # create a scatterer object + +np <- Layer() # and a layer representing the np +na(S) <- 1.33 # set the ambient index +r(np) <- 55 # and the size (radius) of the np in nm + +# load up silver data +palik_ag_vis <- read.table("palik_ag_vis_hb.csv", header=TRUE) +colnames(palik_ag_vis) <- c("lambda","n","k","eps_real","eps_imag") + +lambda_palik <- palik_ag_vis$lambda #convert to nm +n_palik <- palik_ag_vis$n +k_palik <- palik_ag_vis$k + +#interpolation as palik data is sparse +n = 500#number of points +spl_n <- approx(lambda_palik,n_palik,n=n) +spl_k <- approx(lambda_palik,k_palik,n=n) +lambda = spl_n$x + +scatnlay_ar <-function(k){ + m(np) <- spl_n$y[k]+spl_k$y[k]*(0+1i); + lambda(S) <- lambda[k]; + St <- S+np; # PUT THE STACK HERE + Q <- scattnlay(St); +} + +output <- sapply(1:n,scatnlay_ar) +data_mat_np <- data.frame(t(output),lambda,'water') + +Data1 <- data.frame(read.table("/media/mbajb/data/Linux/Comsol Projects/Nottingham Simulations/Data8.txt",skip=5)) +colnames(Data1) <- c("Lambda","Frequency","Qscat") +Data1$scal <- Data1$Qscat/max(Data1$Qscat) +Data1$Lambda <- Data1$Lambda*1e9 + +colnames(data_mat_np) <- c("Qext", "Qsca", "Qabs", "Qbk", "Qpr", "g", "Albedo", "nmax","Lambda","Layers") + +data_mat <- subset(data_mat_np,Lambda<700 & Lambda >350) +#data_mat$Qsca <- data_mat$Qsca/max(data_mat$Qsca) + +palik_ag_melt <- melt(subset(palik_ag_vis,select = c("lambda","n","k")), id=c("lambda"),variable.name = "type", + value.name = "Index") + +scat_plot <- ggplot(data=data_mat, aes(x=Lambda,y=Qext)) +scat_plot + geom_line() + theme_minimal(base_size=22) + xlim(c(300,500)) + xlab("Wavelength (nm)") + +#scat_plot2 <- ggplot(data=Data1, aes(x=Lambda,y=scal)) +#scat_plot2 + geom_line() + theme_minimal(base_size=22) + xlim(c(350,700))+ ylim(c(0,1.2)) + xlab("Wavelength (nm)")