---
title: "Programming-Notes"
author: "Naomi Tague"
date: "January 28, 2016"
output: html_document
---

```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = TRUE)
```

## Code from inclass examples 

Working with data types and functions

```{r powergen}
# create a function that computes power generated from a reservoir
power_gen = function(height, flow, rho=1000, g=9.8, Keff=0.8) {
result = rho * height * flow * g * Keff
return(result)
} 

# examples of function use
power_gen(20,1)
power_gen(height=20, flow=1)
power.estimate = power_gen(height=20,flow=1)
power.estimate

```

## Scoping
```{r scoping, echo=FALSE}

# this will fail
#power.estimate = power_gen(Keff=0.6)
power.estimate = power_gen(height=20, flow=0.6, Keff=0.6)

# this will fail
#Keff 

```
## Nicely documented function 
```{r speed }
#' Power Required by Speed
#'
#' This function determines the power required to keep a vehicle moving
#' a given speed
#' @param cdrag coefficient due to drag default=0.3
#' @param crolling coefficient due to rolling/friction default=0.015
#' @param v vehicle speed (m/2)
#' @param m vehicle mass (kg)
#' @param A area of front of vehicle (m2)
#' @param g acceleration due to gravity (m/s) default=9.8
#' @param pair (kg/m3) default =1.2
#' @return power (W)
power = function(cdrag=0.3, crolling=0.015,pair=1.2,g=9.8,V,m,A) {
P = crolling*m*g*V + 1/2*A*pair*cdrag*V**3
return(P)
}

v=seq(from=0, to=100, by=10)
plot(v, power(V=0.447*v, m=31752, A=25))
lines(v, power(V=0.447*v, m=61752, A=25)) 
```
## Data types in R (that you can exchange with functions)
```{r datafunctions }
mth.names=c("Jan","Feb","Mar","Apr","May","Jun","Jul","Aug","Sep","Oct"
,"Nov","Dec")
reservoir.operation = data.frame(month=mth.names)
reservoir.operation$height = seq(from=32, to=10, by=-2)
reservoir.operation$flowrate = rnorm(n=12, mean=3, sd=0.25)

barplot( power_gen(height=reservoir.operation$height, flow=reservoir.operation$flow), names=reservoir.operation$month)
```

## Using other functions within a function 

```{r morefunctions }

#'
#' This function computes total power generation from a reservoir given its height and flow rate into turbines and number of days (and secs) within those days that the turbines are in operation
#' @param rho Density of water (kg/m3) Default is 1000
#' @param g Acceleration due to gravity (m/sec2) Default is 9.8
#' @param Keff Turbine Efficiency (0-1) Default is 0.8
#' @param height height of water in reservoir (m)
#' @param flow flow rate (m3/sec)
#' @param number of days
#' @param secs in days Default is 86400
#' @author Naomi
#' @examples power_gen(20, 1, 10)
#' @return Power generation (MW)

power_gen_total = function(height, flow, days, secs=86400, rho=1000, g=9.8, Keff=0.8) {

result = rho * height * flow * g * Keff
result = result * days * secs
total =  sum(result)/1e6

return(total)
}

 power_gen_total(reservoir.operation$height, reservoir.operation$flowrate, days=30)
```
## Using lists to return multiple data types from a function 

```{r listreturns }

#' computes profit from  price for forest plot  and Mg/C in that plot
#' @param  price ($)
#' @param  carbon (MgC)
#' @return list with mean, min, and max prices
compute_carbonvalue = function(price, carbon) {

cost.per.carbon = price/carbon
a = mean(cost.per.carbon)
b = max(cost.per.carbon)
c = min(cost.per.carbon)

result = list(avg=a, min=c, max=b)
return(result)
}
# example application, using some data that we 'generate'
obs=data.frame(prices=c(23,44,60,4,2,33,59), forestc=c(59,88,100,10,8,79,300) )
obs
forestvalue = compute_carbonvalue(obs$prices, obs$forestc)
forestvalue
```

## Function that uses factors (and has some wierd comments but ignore that for now
```{r diversity}

#' Describe diversity based on a list of species
#'
#' Compute a species diversity index
#' @param species list of species (names, or code)
#' @return list with the following items
#' \describe{
#' \item{num}{ Number of distinct species}
#' \item{simpson}{Value of simpson diversity index}
#' \item{dominant}{Name of the most frequently occuring species}
#' }
#' @examples
#' compute_diversity(c("butterfly","butterfly","mosquito","butterfly","ladybug","ladybug")))
#' @references
#' http://www.tiem.utk.edu/~gross/bioed/bealsmodules/simpsonDI.html

compute_diversity = function(species) {

species = as.factor(species)
tmp = (summary(species)/sum(summary(species))) ** 2
diversity = sum(tmp)
nspecies = length(summary(species))
tmp = which.max(summary(species))
dominant = names(summary(species)[tmp])
return(list(num=nspecies, simpson=diversity, dominant=dominant))
}
```