Previous changeset 48:99e1c1300fcd (2018-06-05) Next changeset 50:927321ed0322 (2018-08-07) |
Commit message:
Uploaded |
modified:
insect_phenology_model.R |
b |
diff -r 99e1c1300fcd -r 1b6864c5b50a insect_phenology_model.R --- a/insect_phenology_model.R Tue Jun 05 07:52:51 2018 -0400 +++ b/insect_phenology_model.R Tue Jun 05 07:52:59 2018 -0400 |
[ |
b'@@ -6,6 +6,7 @@\n make_option(c("--adult_mortality"), action="store", dest="adult_mortality", type="integer", help="Adjustment rate for adult mortality"),\n make_option(c("--adult_accumulation"), action="store", dest="adult_accumulation", type="integer", help="Adjustment of degree-days accumulation (old nymph->adult)"),\n make_option(c("--egg_mortality"), action="store", dest="egg_mortality", type="integer", help="Adjustment rate for egg mortality"),\n+ make_option(c("--end_date"), action="store", dest="end_date", default=NULL, help="End date for custom date interval"),\n make_option(c("--input_norm"), action="store", dest="input_norm", help="30 year normals temperature data for selected station"),\n make_option(c("--input_ytd"), action="store", dest="input_ytd", default=NULL, help="Year-to-date temperature data for selected location"),\n make_option(c("--insect"), action="store", dest="insect", help="Insect name"),\n@@ -24,6 +25,7 @@\n make_option(c("--plot_generations_separately"), action="store", dest="plot_generations_separately", help="Plot Plot P, F1 and F2 as separate lines or pool across them"),\n make_option(c("--plot_std_error"), action="store", dest="plot_std_error", help="Plot Standard error"),\n make_option(c("--replications"), action="store", dest="replications", type="integer", help="Number of replications"),\n+ make_option(c("--start_date"), action="store", dest="start_date", default=NULL, help="Start date for custom date interval"),\n make_option(c("--young_nymph_accumulation"), action="store", dest="young_nymph_accumulation", type="integer", help="Adjustment of degree-days accumulation (egg->young nymph)")\n )\n \n@@ -31,10 +33,11 @@\n args <- parse_args(parser, positional_arguments=TRUE);\n opt <- args$options;\n \n-add_daylight_length = function(temperature_data_frame, num_rows) {\n- # Return a vector of daylight length (photoperido profile) for\n- # the number of days specified in the input_ytd temperature data\n- # (from Forsythe 1995).\n+add_daylight_length = function(temperature_data_frame) {\n+ # Return temperature_data_frame with an added column\n+ # of daylight length (photoperido profile).\n+ num_rows = dim(temperature_data_frame)[1];\n+ # From Forsythe 1995.\n p = 0.8333;\n latitude = temperature_data_frame$LATITUDE[1];\n daylight_length_vector = NULL;\n@@ -63,6 +66,26 @@\n return(data_frame);\n }\n \n+extract_date_interval_rows = function(df, start_date, end_date) {\n+ date_interval_rows = df[df$DATE >= start_date & df$DATE <= end_date];\n+ return(date_interval_rows);\n+}\n+\n+from_30_year_normals = function(norm_data_frame, start_date_doy, end_date_doy, year) {\n+ # The data we want is fully contained within the 30 year normals data.\n+ first_norm_row = which(norm_data_frame$DOY==start_date_doy);\n+ last_norm_row = which(norm_data_frame$DOY==end_date_doy);\n+ # Add 1 to the number of rows to ensure that the end date is included.\n+ tmp_data_frame_rows = last_norm_row - first_norm_row + 1;\n+ tmp_data_frame = get_new_temperature_data_frame(nrow=tmp_data_frame_rows);\n+ j = 0;\n+ for (i in first_norm_row:last_norm_row) {\n+ j = j + 1;\n+ tmp_data_frame[j,] = get_next_normals_row(norm_data_frame, year, i);\n+ }\n+ return (tmp_data_frame);\n+}\n+\n get_file_path = function(life_stage, base_name, life_stage_nymph=NULL, life_stage_adult=NULL) {\n if (!is.null(life_stage_nymph)) {\n lsi = get_life_stage_index(life_stage, life_stage_nymph=life_stage_nymph);\n@@ -123,7 +146,49 @@\n return(list(p_m, p_se, f1_m, f1_se, f2_m, f2_se))\n }\n \n-get_next_normals_row = function(norm_data_frame, year, is_leap_year, index) {\n+get_new_norm_data_frame = function(is_leap_year, input_norm=NULL, nrow=0) {\n+ # The input_norm data has the following 10 columns:\n+ # STATIONID, LATITUDE, LONGITUDE, ELEV_M, NAME, ST, MMDD, DOY, TMIN, TMAX\n+ column_names = c("STATIONID", "LATITUDE","LONGITUDE", "ELEV_M", "NAME", "ST", "MMDD", "DOY"'..b'eded for plots, some of these values are\n+# being reset here since in some case they were set above.\n start_date = data_list[[2]];\n end_date = data_list[[3]];\n-start_doy_ytd = data_list[[4]];\n-end_doy_ytd = data_list[[5]];\n+prepend_end_doy_norm = data_list[[4]];\n+append_start_doy_norm = data_list[[5]];\n is_leap_year = data_list[[6]];\n-total_days = data_list[[7]];\n-total_days_vector = c(1:total_days);\n-location = data_list[[8]];\n+location = data_list[[7]];\n \n+if (is.null(opt$start_date) && is.null(opt$end_date)) {\n+ # We\'re plotting an entire year.\n+ date_interval = FALSE;\n+ # Display the total number of days in the Galaxy history item blurb.\n+ if (processing_year_to_date_data) {\n+ cat("Number of days year-to-date: ", opt$num_days_ytd, "\\n");\n+ } else {\n+ if (is_leap_year) {\n+ num_days = 366;\n+ } else {\n+ num_days = 365;\n+ }\n+ cat("Number of days in year: ", num_days, "\\n");\n+ }\n+} else {\n+ # FIXME: currently custom date fields are free text, but\n+ # Galaxy should soon include support for a date selector\n+ # at which point this tool should be enhanced to use it.\n+ # Validate start_date.\n+ date_interval = TRUE;\n+ # Calaculate the number of days in the date interval rather\n+ # than using the number of rows in the input temperature data.\n+ start_date = validate_date(opt$start_date);\n+ # Validate end_date.\n+ end_date = validate_date(opt$end_date);\n+ if (start_date >= end_date) {\n+ stop_err("The start date must be between 1 and 50 days before the end date when setting date intervals for plots.");\n+ }\n+ # Calculate the number of days in the date interval.\n+ num_days = difftime(end_date, start_date, units=c("days"));\n+ # Add 1 to the number of days to make the dates inclusive. For\n+ # example, if the user enters a date range of 2018-01-01 to\n+ # 2018-01-31, they likely expect the end date to be included.\n+ num_days = num_days + 1;\n+ if (num_days > 50) {\n+ # We need to restrict date intervals since\n+ # plots render tick marks for each day.\n+ stop_err("Date intervals for plotting cannot exceed 50 days.");\n+ }\n+ # Display the total number of days in the Galaxy history item blurb.\n+ cat("Number of days in date interval: ", num_days, "\\n");\n+}\n # Create copies of the temperature data for generations P, F1 and F2 if we\'re plotting generations separately.\n if (plot_generations_separately) {\n temperature_data_frame_P = data.frame(temperature_data_frame);\n@@ -549,7 +892,7 @@\n }\n \n # Get the ticks date labels for plots.\n-ticks_and_labels = get_x_axis_ticks_and_labels(temperature_data_frame, total_days, start_doy_ytd, end_doy_ytd);\n+ticks_and_labels = get_x_axis_ticks_and_labels(temperature_data_frame, prepend_end_doy_norm, append_start_doy_norm, date_interval);\n ticks = c(unlist(ticks_and_labels[1]));\n date_labels = c(unlist(ticks_and_labels[2]));\n # All latitude values are the same, so get the value for plots from the first row.\n@@ -615,6 +958,7 @@\n }\n }\n # Initialize matrices.\n+total_days = dim(temperature_data_frame)[1];\n if (process_eggs) {\n Eggs.replications = matrix(rep(0, total_days*opt$replications), ncol=opt$replications);\n }\n@@ -778,7 +1122,7 @@\n doy = temperature_data_frame$DOY[row];\n # Photoperiod in the day.\n photoperiod = temperature_data_frame$DAYLEN[row];\n- temp.profile = get_temperature_at_hour(latitude, temperature_data_frame, row, total_days);\n+ temp.profile = get_temperature_at_hour(latitude, temperature_data_frame, row);\n mean.temp = temp.profile[1];\n averages.temp = temp.profile[2];\n averages.day[row] = averages.temp;\n@@ -1459,6 +1803,7 @@\n write.csv(temperature_data_frame_F2, file=file_path, row.names=F);\n }\n \n+total_days_vector = c(1:dim(temperature_data_frame)[1]);\n if (plot_generations_separately) {\n for (life_stage in life_stages) {\n if (life_stage == "Egg") {\n' |