forked from trynthink/scout
-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathplots.R
1467 lines (1366 loc) · 74.8 KB
/
plots.R
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
548
549
550
551
552
553
554
555
556
557
558
559
560
561
562
563
564
565
566
567
568
569
570
571
572
573
574
575
576
577
578
579
580
581
582
583
584
585
586
587
588
589
590
591
592
593
594
595
596
597
598
599
600
601
602
603
604
605
606
607
608
609
610
611
612
613
614
615
616
617
618
619
620
621
622
623
624
625
626
627
628
629
630
631
632
633
634
635
636
637
638
639
640
641
642
643
644
645
646
647
648
649
650
651
652
653
654
655
656
657
658
659
660
661
662
663
664
665
666
667
668
669
670
671
672
673
674
675
676
677
678
679
680
681
682
683
684
685
686
687
688
689
690
691
692
693
694
695
696
697
698
699
700
701
702
703
704
705
706
707
708
709
710
711
712
713
714
715
716
717
718
719
720
721
722
723
724
725
726
727
728
729
730
731
732
733
734
735
736
737
738
739
740
741
742
743
744
745
746
747
748
749
750
751
752
753
754
755
756
757
758
759
760
761
762
763
764
765
766
767
768
769
770
771
772
773
774
775
776
777
778
779
780
781
782
783
784
785
786
787
788
789
790
791
792
793
794
795
796
797
798
799
800
801
802
803
804
805
806
807
808
809
810
811
812
813
814
815
816
817
818
819
820
821
822
823
824
825
826
827
828
829
830
831
832
833
834
835
836
837
838
839
840
841
842
843
844
845
846
847
848
849
850
851
852
853
854
855
856
857
858
859
860
861
862
863
864
865
866
867
868
869
870
871
872
873
874
875
876
877
878
879
880
881
882
883
884
885
886
887
888
889
890
891
892
893
894
895
896
897
898
899
900
901
902
903
904
905
906
907
908
909
910
911
912
913
914
915
916
917
918
919
920
921
922
923
924
925
926
927
928
929
930
931
932
933
934
935
936
937
938
939
940
941
942
943
944
945
946
947
948
949
950
951
952
953
954
955
956
957
958
959
960
961
962
963
964
965
966
967
968
969
970
971
972
973
974
975
976
977
978
979
980
981
982
983
984
985
986
987
988
989
990
991
992
993
994
995
996
997
998
999
1000
# ============================================================================
# Load required packages and files
# ============================================================================
# Define function to load required packages
package_loader <- function(pkg_list) {
options(warn=-1) # Suppress sometimes misleading package load warning messages
# For Windows users, install packages in a directory with administrator priveleges
if (Sys.info()[1]=="Windows"){
# Find Windows user name
user_name = Sys.getenv("USERNAME")
# Create directory for package install (default used by R GUI)
dir_path = file.path('C:', 'Users', user_name, 'Documents','R', 'win-library')
dir.create(dir_path, showWarnings = FALSE, recursive = TRUE)
}else{dir_path = NULL}
# Run through and install/load required packages
for(pkg_name in pkg_list){
# require returns TRUE invisibly if it was able to load package from either the
# default R library path or, in the case of a Windows user, a personal path that the user
# has administrator rights for
if (!require(pkg_name, character.only = TRUE, quietly = TRUE, warn.conflicts=FALSE) &&
!require(pkg_name, lib.loc = dir_path, character.only = TRUE,
quietly = TRUE, warn.conflicts=FALSE)){
# If package was not able to be loaded then download/install using default mirror, notify user
message("Installing R package ", pkg_name)
install.packages(pkg_name, repos='http://cran.us.r-project.org', lib = dir_path, quiet=TRUE)
# Try alternate install mirror if first mirror didn't install package to library directory
if (!require(pkg_name, lib.loc = dir_path, character.only = TRUE,
quietly = TRUE, warn.conflicts=FALSE)){
install.packages(pkg_name, repos='http://cran.cnr.berkeley.edu', lib = dir_path, quiet=TRUE)
}
# Load package after installing
library(pkg_name, character.only = TRUE, quietly=TRUE, warn.conflicts=FALSE, lib.loc = dir_path)
}
}
}
# Load indicated packages
package_loader(c("RColorBrewer", "rjson", "WriteXLS", "stringr", "TeachingDemos", "scales"))
# Get current working directory path
base_dir = getwd()
# Import uncompeted ECM energy, carbon, and cost data
uncompete_results<-fromJSON(file = file.path(base_dir, 'supporting_data','ecm_prep.json'))
# Import competed ECM energy, carbon, and cost data
compete_results_ecms<-fromJSON(file = file.path(base_dir, 'results','ecm_results.json'))
#.Import competed energy, carbon, and cost data summed across all ECMs
compete_results_agg<-fromJSON(file = file.path(base_dir, 'results','agg_results.json'))
# Combine aggregate and individual-level ECM results
compete_results<-c(compete_results_agg, compete_results_ecms)
# Read in global metadata file
glob_dat<-fromJSON(file = file.path(base_dir, 'glob_run_vars.json'))
# ============================================================================
# Set high-level variables needed across multiple plot types
# ============================================================================
# Set ECM adoption scenarios
adopt_scenarios <- glob_dat$'adopt_schemes'
# Set ECM competition scenarios
comp_schemes <- c('uncompeted', 'competed')
# Set full list of ECM names from results file
meas_names <- names(compete_results)
# Set list of ECM names excluding 'All ECMs' (representing results summed across all ECMs),
# the high-level information about site-source energy calculations (stored in 'Energy Output Type' key),
# and information about output breakout granularity (stored in 'Output Resolution')
meas_names_no_all <- meas_names[
(meas_names!= "All ECMs" & meas_names!= "Energy Output Type" &
meas_names!= "Output Resolution" & meas_names!="On-site Generation")]
# Order the list of ECM names excluding 'All ECMs'
meas_names_no_all <- meas_names_no_all[order(meas_names_no_all)]
# Combine the 'All ECMs' name with the ordered list of the individual ECM names
meas_names <- c('All ECMs', meas_names_no_all)
# Set years in modeling time horizon and reorganize in ascending order
years<-row.names(as.matrix(
compete_results[[meas_names[1]]]$'Markets and Savings (Overall)'[[
adopt_scenarios[1]]]$'Baseline Energy Use (MMBtu)'))
# Set an intended starting and ending year for the plotted results
start_yr = 2015
end_yr = max(years)
# Set the years to take a 'snapshot' of certain results in;
# if the measure set was prepared with public cost adders, set
# this year to 2020, as these adders are reflective of current/
# near term utility conditions from EPA data; also flag for the particular
# type of public cost adder (low/high), if any
if (grepl("PHC", meas_names[2], fixed=TRUE) == TRUE){
snap_yr = "2020"
snap_yr_set = c("2020")
if (grepl("low", meas_names[2], fixed=TRUE) == TRUE){
phc_flag = "low"
}else{
phc_flag = "high"
}
}else{
snap_yr = "2050"
snap_yr_set = c("2030", "2050")
phc_flag = FALSE}
# Filter and order the year range
years<-years[(years>=start_yr)&(years<=end_yr)]
years<-years[order(years)]
# Set list of possible climate zones and associated colors/legend entries for aggregate savings
# and cost effectiveness plots
# Set various region names lists
aia_reg_names <- c("AIA CZ1", "AIA CZ2", "AIA CZ3", "AIA CZ4", "AIA CZ5")
emm_reg_names_detail <- c(
'TRE', 'FRCC', 'MISW', 'MISC', 'MISE', 'MISS',
'ISNE', 'NYCW', 'NYUP', 'PJME', 'PJMW', 'PJMC',
'PJMD', 'SRCA', 'SRSE', 'SRCE', 'SPPS', 'SPPC',
'SPPN', 'SRSG', 'CANO', 'CASO', 'NWPP', 'RMRG', 'BASN')
emm_reg_names <- c(
'Northwest', 'Great Basin', 'California', 'Rocky Mountains',
'Upper Midwest', 'Lower Midwest', 'Lakes/Mid-Atl.', 'Texas',
'Southwest', 'Southeast', 'Northeast')
state_names_detail <- c(
'AL', 'AK', 'AZ', 'AR', 'CA', 'CO', 'CT', 'DE', 'DC', 'FL',
'GA', 'HI', 'ID', 'IL', 'IN', 'IA', 'KS', 'KY', 'LA', 'ME',
'MD', 'MA', 'MI', 'MN', 'MS', 'MO', 'MT', 'NE', 'NV', 'NH',
'NJ', 'NM', 'NY', 'NC', 'ND', 'OH', 'OK', 'OR', 'PA', 'RI',
'SC', 'SD', 'TN', 'TX', 'UT', 'VT', 'VA', 'WA', 'WV', 'WI',
'WY')
state_names <- c(
'New England', 'Mid Atlantic', 'East North Central',
'West North Central', 'South Atlantic', 'East South Central',
'West South Central', 'Mountain', 'Pacific')
# Pull czones from first available ECM
czones_test = toString(
compete_results[[meas_names[2]]]$'Filter Variables'$'Applicable Regions')
# Check for AIA regional breakouts; if not there, check for EMM regional breakouts;
# if not there, assume state breakouts
if (grepl("AIA", czones_test, fixed=TRUE) == TRUE){
# Set region legend names
czones_out<-aia_reg_names
# Set region legend colors
czones_out_col<-brewer.pal(length(czones_out), "Dark2")
czones_out_lgnd<-czones_out
# Set alternate region flags
emm_flag = FALSE
emm_det_flag = FALSE
state_flag = FALSE
state_det_flag = FALSE
}else if (
any(grepl(paste(c(emm_reg_names, emm_reg_names_detail),
collapse="|"), czones_test) == TRUE)){
# Set region legend names
if (compete_results_agg$`Output Resolution`%in% c(
"detail", "detail_reg", "detail_reg_bldg", "detail_reg_fuel")){
czones_out<-emm_reg_names_detail
emm_det_flag = TRUE
}else{
czones_out<-emm_reg_names
emm_det_flag = FALSE
}
# Set region legend colors; in the EMM case, extend the color list to
# accomodate the larger number of region categories
czones_out_col<-colorRampPalette(brewer.pal(8, "Dark2"))(length(czones_out))
czones_out_lgnd<-czones_out
# Set alternate region flags
emm_flag = TRUE
state_flag = FALSE
state_det_flag = FALSE
}else{
if (compete_results_agg$`Output Resolution`%in% c(
"detail", "detail_reg", "detail_reg_bldg", "detail_reg_fuel")){
czones_out<-state_names_detail
state_det_flag = TRUE
}else{
czones_out<-state_names
state_det_flag = FALSE
}
# Set region legend colors; in the state case, extend the color list to
# accomodate the larger number of region categories
czones_out_col<-colorRampPalette(brewer.pal(8, "Dark2"))(length(czones_out))
czones_out_lgnd<-czones_out
# Set alternate region flags
emm_flag = FALSE
emm_det_flag = FALSE
state_flag = TRUE
}
# Set list of possible building classes and associated colors for aggregate savings plot
if (compete_results_agg$`Output Resolution`%in%c(
"detail", "detail_bldg", "detail_reg_bldg", "detail_bldg_fuel")){
bclasses_out_agg<-c('Single Family Homes', 'Multi Family Homes',
'Hospitals', 'Large Offices', 'Small/Medium Offices',
'Retail', 'Hospitality', 'Education', 'Assembly/Other',
'Warehouses')
bclasses_out_finmets<-c(
list(c('Single Family Homes', 'Multi Family Homes')),
list(c('Hospitals', 'Large Offices', 'Small/Medium Offices',
'Retail', 'Hospitality', 'Education', 'Assembly/Other',
'Warehouses')))
}else{
bclasses_out_agg<-c('Residential (New)', 'Residential (Existing)',
'Commercial (New)', 'Commercial (Existing)')
bclasses_out_finmets<-c(
list(c('Residential (New)', 'Residential (Existing)')),
list(c('Commercial (New)', 'Commercial (Existing)')))
}
bclasses_out_agg_col<-colorRampPalette(brewer.pal(8, "Dark2"))(length(bclasses_out_agg))
# Set list of possible building classes and associated shapes/legend entries for
# cost effectiveness plot
bclasses_out_finmets_shp<-c(21, 22)
bclasses_out_finmets_lgnd<-c('Residential', 'Commercial')
# Set list of possible end uses and associated colors/legend entries for
# aggregate savings plot
euses_out_agg<-c('Heating (Equip.)', 'Cooling (Equip.)', 'Heating (Env.)', 'Cooling (Env.)',
'Ventilation', 'Lighting', 'Water Heating', 'Refrigeration', 'Cooking', 'Electronics', 'Other')
euses_out_agg_col<-brewer.pal(length(euses_out_agg), "Paired")
# Set list of possible end use names from the raw data and associated colors/legend entries for
# cost effectiveness plot
euses_out_finmets<-c(
list(c('Heating (Equip.)', 'Cooling (Equip.)', 'Ventilation')),
list(c('Heating (Env.)', 'Cooling (Env.)')),
list('Lighting'),
list('Water Heating'), list('Refrigeration'), list('Cooking'),
list('Computers and Electronics'), list('Other'))
euses_out_finmets_col<-brewer.pal(length(euses_out_finmets), "Dark2")
euses_out_finmets_lgnd<-c(
'HVAC', 'Envelope', 'Lighting', 'Water Heating', 'Refrigeration',
'Cooking', 'Electronics', 'Other')
# Determine axis label parameters to use across plots
# Determine site vs. source energy and (if source) captured vs. fossil
# site-source conversion calculation methods, to write to label
if (compete_results_agg$`Energy Output Type`[1] == "site"){
e_site_source = "Site Energy"
ss_calc = ""
}else if (compete_results_agg$`Energy Output Type`[1] == "captured"){
e_site_source = "Primary Energy"
ss_calc = ", CE S-S"
}else{
e_site_source = "Primary Energy"
ss_calc = ", FF S-S"
}
# Case without TSV metrics - set conversion factors, units, and
# axis scaling factors for aggregate savings plots
if (compete_results_agg$`Energy Output Type`[2] == "NA"){
# If results are EMM-resolved and public health cost adders
# are not present, convert quads to TWh for reporting;
# if results are state-resolved, used TBtu for reporting
if ((emm_flag == TRUE)&(phc_flag == FALSE)){
e_conv_emm_st = 293.07
e_axis_units = "TWh"
}else if (state_flag == TRUE){
e_conv_emm_st = 1000
e_axis_units = "TBtu"
}else{
e_conv_emm_st = 1
e_axis_units = "Quads"
}
# No further conversions needed for CO2 or cost
c_conv_emm = 1
c_axis_units = "Mt"
cs_conv_emm = 1
cs_axis_units = "Billion $"
# Set scaling for aggregate savings plots
cex_ann_ax = 0.85
cex_cum_ax = 0.85
# No further TSV metric assumptions text to append
append_txt = ")"
# Case with TSV metrics - set conversion factors, units, scaling factors
# for aggregate savings plots, and text to append to axis label indicating
# TSV metric assumptions
}else{
# When TSV metrics are present, assume that GW (hourly savings, abbreviated in
# the output as "Hr.") or GWh (multi-hour savings, abbreviated in the output as
# "Prd.") is the correct energy unit to use (savings will generally be smaller)
e_conv_emm_st = 293071.07
if (compete_results_agg$`Energy Output Type`[1] == "Prd."){
e_axis_units = "GWh"
}else{
e_axis_units = "GW"
}
# When TSV metrics are present, convert million metric tons to thousand metric
# tons and billion $ to million $ (savings will generally be smaller)
c_conv_emm = 1000
c_axis_units = "Tt"
cs_conv_emm = 1000
cs_axis_units = "Million $"
# Set scaling for aggregate savings plots (a bit lower than when there are no
# TSV metrics, to accomodate additional axis label text for TSV assumptions)
cex_ann_ax = 0.75
cex_cum_ax = 0.75
# Set text to append to axis label indicating key elements of the TSV metrics
# assumptions
append_txt = paste(paste(", ",
compete_results_agg$`Energy Output Type`[2], " ",
compete_results_agg$`Energy Output Type`[3], " ",
compete_results_agg$`Energy Output Type`[4], " ",
compete_results_agg$`Energy Output Type`[5], ")", sep=""))
}
# ============================================================================
# Set high-level variables needed to generate individual ECM plots
# ============================================================================
# Set names of plot files and axes
# Plot of individual ECM energy, carbon, and cost totals
file_names_ecms <- c('Total Energy', 'Total CO2', 'Total Cost')
plot_names_ecms <- c('Total Energy', expression("Total"~ CO[2]), 'Total Cost')
# Set individual ECM plot y axis labels
plot_axis_labels_ecm<-c(
paste(e_site_source, ' (', e_axis_units, ss_calc, append_txt, sep=""),
as.expression(bquote(CO[2] ~ "Emissions ("*.(c_axis_units)*.(append_txt))),
paste('Energy Cost (', cs_axis_units, append_txt, sep=""))
# Set colors for uncompeted baseline, efficient and low/high results
plot_col_uc_base = "gray60"
plot_col_uc_eff = "gray80"
plot_col_uc_lowhigh = "gray90"
# Set variable names to use in accessing all uncompeted energy, carbon, and cost results from JSON data
var_names_uncompete <- c('energy', 'carbon', 'cost')
results_folder_names <- c('energy', 'co2', 'cost')
# Set output units for each variable type
var_units <- c(e_axis_units, c_axis_units, cs_axis_units)
# Set variable names to use in accessing competed baseline energy, carbon, and cost results from
# JSON data. Note that each variable potentially has a '(low)' and '(high)' variant in the JSON.
var_names_compete_base_m <- c(
'Baseline Energy Use (MMBtu)', 'Baseline CO₂ Emissions (MMTons)', 'Baseline Energy Cost (USD)')
var_names_compete_base_l <- c(
'Baseline Energy Use (low) (MMBtu)', 'Baseline CO₂ Emissions (low) (MMTons)',
'Baseline Energy Cost (low) (USD)')
var_names_compete_base_h <- c(
'Baseline Energy Use (high) (MMBtu)', 'Baseline CO₂ Emissions (high) (MMTons)',
'Baseline Energy Cost (high) (USD)')
# Set variable names to use in accessing competed efficient energy, carbon, and cost results from
# JSON data. Note that each variable potentially has a '(low)' and '(high)' variant in the JSON.
var_names_compete_eff_m <- c(
'Efficient Energy Use (MMBtu)', 'Efficient CO₂ Emissions (MMTons)', 'Efficient Energy Cost (USD)')
var_names_compete_eff_l <- c(
'Efficient Energy Use (low) (MMBtu)', 'Efficient CO₂ Emissions (low) (MMTons)',
'Efficient Energy Cost (low) (USD)')
var_names_compete_eff_h <- c(
'Efficient Energy Use (high) (MMBtu)', 'Efficient CO₂ Emissions (high) (MMTons)',
'Efficient Energy Cost (high) (USD)')
# ============================================================================
# Set high-level variables needed to generate aggregated savings plots
# ============================================================================
# File names for aggregated ECM savings plots
plot_names_agg <- c('Total Energy Savings', 'Total Avoided CO2', 'Total Cost Savings')
# Titles for aggregated ECM savings plots
plot_titles_agg <- c('Energy Savings', expression("Avoided"~ CO[2]), 'Cost Savings')
# Set aggregate annual/cumulative savings plot y axis labels
plot_axis_labels_agg_ann<-c(
paste('Annual ', e_site_source, ' Savings (', e_axis_units, ss_calc, append_txt, sep=""),
as.expression(bquote("Annual Avoided"~ CO[2] ~ "Emissions ("*.(c_axis_units)*.(append_txt))),
paste('Annual Energy Cost Savings (', cs_axis_units, append_txt, sep=""))
plot_axis_labels_agg_cum<-c(
paste('Cumulative ', e_site_source, ' Savings (', e_axis_units, ss_calc, append_txt, sep=""),
as.expression(bquote("Cumulative Avoided"~ CO[2] ~ "Emissions ("*.(c_axis_units)*.(append_txt))),
paste('Cumulative Energy Cost Savings (', cs_axis_units, append_txt, sep=""))
# Set names of variables used to retrieve savings data to aggregate
var_names_compete_save <- c(
'Energy Savings (MMBtu)', 'Avoided CO₂ Emissions (MMTons)', 'Energy Cost Savings (USD)')
# Set names of variables to filter aggregated savings
filter_var<-c('Climate Zone', 'Building Class', 'End Use')
# Set transparent background for legend
transparent_back<-alpha("white", 0.5)
# ============================================================================
# Set high-level variables needed to generate ECM cost effectiveness plots
# ============================================================================
# Names for ECM cost effectiveness plots
plot_names_finmets = c('Cost Effective Energy Savings', 'Cost Effective Avoided CO2',
'Cost Effective Operation Cost Savings')
# X axis labels for cost effectiveness plots
plot_axis_labels_finmets_x<-c(
paste(snap_yr, " ", e_site_source, ' Savings (Competed, ', e_axis_units, ss_calc, append_txt, sep=""),
as.expression(bquote(.(snap_yr) ~ "Avoided" ~ CO[2] ~ "Emissions (Competed, "*.(c_axis_units)*.(append_txt))),
paste(snap_yr, ' Energy Cost Savings (Competed, ', cs_axis_units, append_txt, sep=""))
# Y axis labels for cost effectiveness plots
plot_axis_labels_finmets_y<-c(
"IRR (%)", "Payback (years)", "CCE ($/MMBtu saved)", expression("CCC ($/t" ~ CO[2] ~ " avoided)"))
# Financial metric titles
plot_title_labels_finmets<-c(
"Internal Rate of Return (IRR)", "Simple Payback", "Cost of Conserved Energy (CCE)",
expression("Cost of Conserved" ~ CO[2] ~ "(CCC)"))
# Default plot limits for each financial metric
plot_lims_finmets <- c(
list(c(-50, 150)), list(c(0, 25)), list(c(-50, 150)), list(c(-500, 1000)))
# Cost effectiveness threshold lines for each financial metric
plot_ablines_finmets <- c(0, 5, 13, 73)
# Financial metric type and key names for retrieving JSON data on each
fin_metrics <- c('IRR (%)', 'Payback (years)', 'Cost of Conserved Energy ($/MMBtu saved)',
'Cost of Conserved CO₂ ($/MTon CO₂ avoided)')
# ============================================================================
# Set high-level variables needed to generate XLSX data
# ============================================================================
# Define cost-effectiveness metrics column names
xlsx_names_finmets <- c(paste("IRR (%),", snap_yr),
paste("Payback (years),", snap_yr),
paste("CCE ($/MMBtu saved),", snap_yr),
paste("CCC ($/tCO2 avoided),", snap_yr))
# ============================================================================
# Loop through all adoption scenarios, plotting variables, ECMs, and
# competition scenarios, gather data needed for plots, and implement plots
# ============================================================================
# Loop through all adoption scenarios
for (a in 1:length(adopt_scenarios)){
# Set plot colors for competed baseline, efficient, and low/high results
# (varies by adoption scenario); also set Excel summary data file name for adoption scenario
if (adopt_scenarios[a] == "Technical potential"){
# Set plot colors
plot_col_c_base = "midnightblue"
plot_col_c_eff = "lightskyblue"
# # Set Excel summary data file name
xlsx_file_name = file.path(base_dir, 'results', 'plots', 'tech_potential', "Summary_Data-TP.xlsx")
}else{
# Set plot colors
plot_col_c_base = "red3"
plot_col_c_eff = "pink"
# # Set Excel summary data file name
xlsx_file_name = file.path(base_dir, 'results', 'plots', 'max_adopt_potential',
"Summary_Data-MAP.xlsx")
}
# Preallocate list for variable names to be used later to export data
# to xlsx-formatted Excel files
xlsx_var_name_list <- vector("character", length=length(var_names_uncompete))
# Loop through all plotting variables
for (v in 1:length(var_names_uncompete)){
# Finalize column names for the annual total energy, carbon, or cost results that are written
# to the XLSX file
xlsx_names_years = matrix(NA, length(years))
for (yn in (1:length(xlsx_names_years))){
# Append variable units to each year name
xlsx_names_years[yn] = paste(years[yn], " (", var_units[v], ")", sep="")
}
# Define all column names for XLSX file
col_names_xlsx<- c('ECM Name', 'Results Scenario', 'Climate Zones', 'Building Classes', 'End Uses',
xlsx_names_finmets, xlsx_names_years)
# Initialize data frame to write to Excel worksheet (note: number of rows equals to
# number of ECMs * number of results scenarios (baseline/efficient + competed/uncompeted) plus
# two additional rows to accommodate baseline/efficient competed results summed across all ECMs)
xlsx_data<-data.frame(matrix(ncol = length(col_names_xlsx),
nrow = (length(meas_names)*4 + 2)))
# Set column names for the worksheet
colnames(xlsx_data) = col_names_xlsx
# Set a factor to convert the results data to final plotting units for given variable
# (quads for energy, Mt for CO2, and billion $ for cost)
if ((var_names_uncompete[v] == "energy") | (var_names_uncompete[v] == "cost")){
unit_translate = 1/1000000000 # converts energy from MBtu -> quads or cost from $ -> billion $
# Layer on a second conversion factor to handle cases where alternate EMM regions are used
# (energy goes to TWh and/or TSV metrics are used (energy goes to GW/GWh, cost goes to million $))
if (var_names_uncompete[v] == "energy"){
unit_translate = unit_translate * e_conv_emm_st
}else if (var_names_uncompete[v] == "cost"){
unit_translate = unit_translate * cs_conv_emm
}
}else{
unit_translate = 1 # CO2 results data are already imported in Mt units
# Layer on a second conversion factor to handle cases where TSV metrics are used
# (carbon goes to thousand metric tons)
unit_translate = unit_translate * c_conv_emm
}
# Initialize a matrix for storing individual ECM energy, carbon, or cost totals (3 rows
# accommodate mean, low, and high totals values; 4 columns accommodate 2 outputs
# (baseline and efficient) x 2 adoption scenarios)
results <- matrix(list(), 3, 2*length(comp_schemes))
# Initialize a matrix for storing aggregated ECM savings results (3 rows accommodate
# 3 filtering variables (climate zone, building class, end use); 3 columns accomodate
# the possible name categories for each filtering variable, the aggregated savings data
# associated with each of those categories, and the colors associated with those categories
results_agg <- matrix(list(), length(filter_var), 3)
# Initialize climate zone names, annual by-climate aggregated savings data, and line colors
results_agg[1,1:3] <- c(list(czones_out),
list(rep(list(matrix(0, length(years))), length(czones_out))),
list(czones_out_col))
# Initialize building class names, annual by-building aggregated savings data and line colors
results_agg[2,1:3] <- c(list(bclasses_out_agg),
list(rep(list(matrix(0, length(years))), length(bclasses_out_agg))),
list(bclasses_out_agg_col))
# Initialize end use names, annual by-end use aggregated savings data and line colors
results_agg[3,1:3] <- c(list(euses_out_agg),
list(rep(list(matrix(0, length(years))), length(euses_out_agg))),
list(euses_out_agg_col))
# Initialize a matrix for storing individual ECM financial metrics/savings and filter
# variable data
results_finmets <- matrix(list(), length(meas_names_no_all), (length(fin_metrics) + 4))
# Initialize uncertainty flag for the adoption scenario
uncertainty = FALSE
# Set the file name for the plot based on the adoption scenario and plotting variable
if (adopt_scenarios[a] == 'Technical potential'){
# ECM energy, carbon, and cost totals
plot_file_name_ecms = file.path(
base_dir, 'results', 'plots', 'tech_potential', results_folder_names[v],
paste(file_names_ecms[v],"-TP", sep=""))
# Aggregate energy, carbon, and cost savings
plot_file_name_agg = file.path(
base_dir, 'results', 'plots', 'tech_potential', results_folder_names[v],
paste(plot_names_agg[v],"-TP", sep=""))
# ECM cost effectiveness
plot_file_name_finmets = file.path(
base_dir, 'results', 'plots', 'tech_potential', results_folder_names[v],
paste(plot_names_finmets[v],"-TP.pdf", sep=""))
}else{
# ECM energy, carbon, and cost totals
plot_file_name_ecms = file.path(
base_dir, 'results', 'plots', 'max_adopt_potential', results_folder_names[v],
paste(file_names_ecms[v],"-MAP", sep = ""))
# Aggregate energy, carbon, and cost savings
plot_file_name_agg = file.path(
base_dir, 'results', 'plots', 'max_adopt_potential', results_folder_names[v],
paste(plot_names_agg[v],"-MAP", sep=""))
# ECM cost effectiveness
plot_file_name_finmets = file.path(
base_dir, 'results', 'plots', 'max_adopt_potential', results_folder_names[v],
paste(plot_names_finmets[v],"-MAP.pdf", sep=""))
}
# Open PDF device for plotting each ECM's energy, carbon, or cost totals
pdf(paste(plot_file_name_ecms, "-byECM.pdf", sep=""),width=15.5,height=14)
# Set number of rows and columns per page in PDF output
par(mfrow=c(4,4))
# Reconfigure space around each side of the plot for best fit
par(mar=c(5.1,5.1,3.1,2.1))
# Loop through all ECMs
for (m in 1:length(meas_names)){
# Add ECM name to Excel worksheet data frame
# Set appropriate starting row for the ECM's data
if (meas_names[m] == "All ECMs"){
# For the 'All ECMs' data, start at the beginning of the file
row_ind_start = 1
}else{
# Otherwise, leave the first 4 rows for the 'All ECMs' total uncompeted/
# competed energy, carbon, and cost data, plus another set of rows for 'All ECMs'
# energy, carbon, and cost savings totals, first summarized across all climate
# zones, building types, and end uses (1 row) and subsequently broken out by
# climate zone (5, 25, or 48 rows, depending on regionality in the data),
# building type (4 rows), and end use (11 rows)
row_ind_start = (m-1)*4 + 1 + (1 + length(czones_out) + 4 + 11)
}
xlsx_data[row_ind_start:(row_ind_start + 3), 1] = meas_names[m]
# Set applicable climate zones, end uses, and building classes for ECM
# and add to Excel worksheet data frame
czones = toString(compete_results[[meas_names[m]]]$'Filter Variables'$'Applicable Regions')
bldg_types = toString(compete_results[[meas_names[m]]]$'Filter Variables'$'Applicable Building Classes')
end_uses = toString(compete_results[[meas_names[m]]]$'Filter Variables'$'Applicable End Uses')
xlsx_data[row_ind_start:(row_ind_start + 3), 3] = czones
xlsx_data[row_ind_start:(row_ind_start + 3), 4] = bldg_types
xlsx_data[row_ind_start:(row_ind_start + 3), 5] = end_uses
# If there are more than three end use names, set a single end use name of 'Multiple' such that
# the end use name label will fit easily within each plot region
if (str_count(end_uses, ",") > 1){
end_uses = "Multiple"
}
# Find the index for accessing the item in the list of uncompeted results that corresponds
# to energy, carbon, or cost total data for the current ECM. Note: competed energy, carbon, and
# cost totals are accessed directly by ECM name, and do not require an index
for (uc in 1:length(uncompete_results)){
if (uncompete_results[[uc]]$'name' == meas_names[m]){
uc_name_ind = uc
}
}
# Set the appropriate database of ECM financial metrics data (used across both competition schemes)
results_database_finmets = compete_results[[meas_names[m]]]$'Financial Metrics'
# Loop through all competition schemes
for (cp in 1:length(comp_schemes)){
# Add name of results scenario (baseline/efficient + competed/uncompeted)
# to Excel worksheet data frame
xlsx_data[(row_ind_start + (cp-1)*2):(row_ind_start + (cp-1)*2 + 1), 2] = c(
paste("Baseline ", comp_schemes[cp], sep = ""), paste("Efficient ", comp_schemes[cp], sep = ""))
# Set matrix for temporarily storing finalized baseline and efficient results
r_temp <- matrix(NA, 6, length(years))
# Find data for uncompeted energy, carbon, and/or cost; exclude the 'All ECMs' case
# (only competed data may be summed across all ECMs)
if ((comp_schemes[cp] == "uncompeted")&(meas_names[m] != "All ECMs")){
# Set the appropriate database of uncompeted energy, carbon, or cost totals
# (access keys vary based on plotted variable)
if (var_names_uncompete[v] != "cost"){
results_database =
uncompete_results[[uc_name_ind]]$'markets'[[adopt_scenarios[a]]]$
'master_mseg'[[var_names_uncompete[v]]]$'total'
}else{
results_database =
uncompete_results[[uc_name_ind]]$'markets'[[adopt_scenarios[a]]]$
'master_mseg'[[var_names_uncompete[v]]]$'energy'$'total'
}
# Order the uncompeted ECM energy, carbon, or cost totals by year and determine low/high
# bounds on each total value (if applicable)
for (yr in 1:length(years)){
r_temp[1:3, yr] = results_database$'baseline'[years[yr]][[1]]
# Set mean, low, and high values for case with ECM input/output uncertainty
if (length(results_database$'efficient'[years[1]][[1]]) > 1){
# Take mean of list of values from uncompeted results
r_temp[4, yr] = mean(results_database$'efficient'[years[yr]][[1]])
# Take 5th/95th percentiles of list of values from uncompeted results
r_temp[5:6, yr] = quantile(results_database$'efficient'[years[yr]][[1]], c(0.05, 0.95))
uncertainty = TRUE
# Set mean, low, and high values for case without ECM input/output uncertainty
# (all values equal to mean value)
}else{
r_temp[4:6, yr] = results_database$'efficient'[years[yr]][[1]]
}
}
# Find data for competed energy, carbon, and/or cost
}else{
# Set the appropriate database of ECM competed energy, carbon, or cost totals
results_database = compete_results[[meas_names[m]]]$
'Markets and Savings (Overall)'[[adopt_scenarios[a]]]
# Set the appropriate database of ECM competed energy, carbon, or cost totals
# broken down by climate zone, building class, and end use; exclude the 'All ECMs'
# case (breakdowns are not calculated for this case)
if (meas_names[m] != "All ECMs"){
# Set the appropriate database of ECM competed energy, carbon , or cost savings,
# which are broken out by each of the three filtering variables (climate zone,
# building class, end use)
results_database_agg = compete_results[[meas_names[m]]]$
'Markets and Savings (by Category)'[[adopt_scenarios[a]]][[var_names_compete_save[v]]]
# Set the appropriate database of ECM data for categorizing cost effectiveness outcomes
# based on climate zone, building type, and end use
results_database_filters = compete_results[[meas_names[m]]]$'Filter Variables'
}
# Order competed ECM energy, carbon, or cost totals by year and determine low/high bounds
# on each total value (if applicable)`
for (yr in 1:length(years)){
base_uncertain_check = sapply(c("Baseline", "low"), grepl, names(results_database))
if (length(which((base_uncertain_check[,1]&base_uncertain_check[,2])==TRUE)>0)) {
# Take mean value output directly from competed results
r_temp[1, yr] = results_database[[var_names_compete_base_m[v]]][years[yr]][[1]]
# Take 'low' value output directly from competed results (represents 5th percentile)
r_temp[2, yr] = results_database[[var_names_compete_base_l[v]]][years[yr]][[1]]
# Take 'high' value output directly from competed results (represents 5th percentile)
r_temp[3, yr] = results_database[[var_names_compete_base_h[v]]][years[yr]][[1]]
# Flag output uncertainty in the current plot
uncertainty = TRUE
# Set mean, low, and high values for case without ECM input/output uncertainty
# (all values equal to mean value)
}else{
# Take 'low' value output directly from competed results (represents 5th percentile)
r_temp[1:3, yr] = results_database[[var_names_compete_base_m[v]]][years[yr]][[1]]
}
# Set mean, low, and high values for case with ECM input/output uncertainty
eff_uncertain_check = sapply(c("Efficient", "low"), grepl, names(results_database))
if (length(which((eff_uncertain_check[,1]&eff_uncertain_check[,2])==TRUE)>0)) {
# Take mean value output directly from competed results
r_temp[4, yr] = results_database[[var_names_compete_eff_m[v]]][years[yr]][[1]]
# Take 'low' value output directly from competed results (represents 5th percentile)
r_temp[5, yr] = results_database[[var_names_compete_eff_l[v]]][years[yr]][[1]]
# Take 'high' value output directly from competed results (represents 95th percentile)
r_temp[6, yr] = results_database[[var_names_compete_eff_h[v]]][years[yr]][[1]]
# Flag output uncertainty in the current plot
uncertainty = TRUE
# Set mean, low, and high values for case without ECM input/output uncertainty
# (all values equal to mean value)
}else{
r_temp[4:6, yr] = results_database[[var_names_compete_eff_m[v]]][years[yr]][[1]]
}
# Find the amount of the ECM's energy, carbon, or cost savings that can be
# attributed to each of the three aggregate savings filtering variables
# (climate zone, building class, end use) and add those savings to the
# aggregated total for each filter variable across all ECMs; exclude the 'All ECMs'
# case (filter variable breakdowns are not calculated for this case)
if (meas_names[m] != "All ECMs"){
# Loop through the three variables used to filter aggregated savings, each
# represented by a row in the 'results_agg' matrix
for (fv in (1:nrow(results_agg))){
# Set the name categories associated with the given savings filter
# variable (e.g., 'Residential (New)', 'Residential (Existing)', etc. for
# building class)
fv_opts <- results_agg[fv, 1][[1]]
# Initialize a vector used to store the ECM's energy, carbon, or cost
# savings that are attributable to the given filter variable; this
# vector must be as long as the number of category names for the
# filter variable, so a savings total can be stored for each category
add_val = matrix(0, length(fv_opts))
# Retrieve the savings data for the ECM that is attributable to each
# filter variable name category. The data are stored in a three-level
# nested dict with climate zone at the top level, building class at
# the middle level, and end use at the bottom level. All three of these
# levels must be looped through to retrieve the savings data.
# Loop through all climate zones
for (levone in (1:length(results_agg[1, 1][[1]]))){
# Set the climate zone name to use in proceeding down to the
# building class level of the dict
czone = results_agg[1, 1][[1]][[levone]]
# Reduce the dict to the building class level
r_agg_temp = results_database_agg[[czone]]
# Loop through all building classes
for (levtwo in (1:length(results_agg[2, 1][[1]]))){
# Set the building class name to use in proceeding down to the
# end use level of the dict
bldg = results_agg[2, 1][[1]][[levtwo]]
# Reduce the dict to the end use level
r_agg_temp = results_database_agg[[czone]][[bldg]]
# Loop through all end uses
for (levthree in (1:length(results_agg[3, 1][[1]]))){
# Set the end use name to use in retrieving data values
euse = results_agg[3, 1][[1]][[levthree]]
# Reset the predefined 'Electronics' end use name (short for
# later use in plot legends) to the longer 'Computers and Electronics'
# name used in the dict
if (euse == "Electronics"){
euse = "Computers and Electronics"
}
# Reduce the dict to the level of the retrieved data values
r_agg_temp = results_database_agg[[czone]][[bldg]][[euse]]
# If data values exist, add them to the ECM's energy/carbon/cost
# savings-by-filter variable vector initialized above
if (length(r_agg_temp)!=0){
# Determine which index to use in adding the retrieved data to
# the ECM's energy/carbon/cost savings-by-filter variable vector
if (fv == 1){
index = levone
}else if (fv == 2){
index = levtwo
}else{
index = levthree
}
# Add retrieved data to ECM's savings-by-filter variable vector;
# handle case where end use savings are further split out by
# fuel type (may be electric/non-electric or more detailed)
if (years[yr]%in%names(r_agg_temp)){
add_val[index] = add_val[index] + r_agg_temp[years[yr]][[1]]
}else{
fts = names(r_agg_temp)
for (fuel in fts){
if (length(r_agg_temp[[fuel]])!=0){
add_val[index] = add_val[index] +
r_agg_temp[[fuel]][years[yr]][[1]]
}
}
}
}
}
}
}
# Add ECM's savings-by-filter variable vector data to the aggregated total for
# each filter variable across all ECMs
for (fvo in 1:length(fv_opts)){
results_agg[fv, 2][[1]][[fvo]][yr] = results_agg[fv, 2][[1]][[fvo]][yr] + add_val[fvo]
}
}
# If cycling through the year in which snapshots of ECM cost effectiveness are taken,
# retrieve the ECM's competed financial metrics, savings, and filter variable data needed
# to develop those snapshots for the cost effectiveness plots
if (years[yr] == snap_yr){
# Retrieve ECM competed portfolio-level and consumer-level financial metrics data
for (fm in 1:length(fin_metrics)){
# Multiply IRR fractions in JSON data by 100 to convert to final % units
if ((fin_metrics[fm]) == "IRR (%)"){
unit_translate_finmet = 100
}else{
unit_translate_finmet = 1
}
# Consumer-level data are NOT keyed by adoption scenario
results_finmets[(m - 1), fm][[1]] =
results_database_finmets[[fin_metrics[fm]]][[years[yr]]][[1]] * unit_translate_finmet
# Replace all 99900 values with 999 (proxy for NaN)
results_finmets[(m - 1), fm][[1]][results_finmets[(m - 1), fm][[1]] == 99900] = 999
}
# Write ECM cost effectiveness metrics data to XLSX sheet
xlsx_data[row_ind_start:(row_ind_start + 3),
6: (6 + length(plot_title_labels_finmets) - 1)] =
as.matrix(results_finmets[(m - 1),1:4])
# Retrieve, ECM energy, carbon, or cost savings data, convert to final units
results_finmets[(m - 1), (length(fin_metrics)+1)][[1]] =
results_database[[var_names_compete_save[v]]][years[yr]][[1]] * unit_translate
# Determine the outline color, shape, and fill color parameters needed to
# distinguish ECM points on the cost effectiveness plots by their climate
# zone, building type, and end use categories
# # Determine appropriate ECM point outline color for applicable climate
# # zones
# # Set ECM's applicable climate zones
# czones = results_database_filters$'Applicable Climate Zones'
# # Match applicable climate zones to climate zone names used in plotting
# czone_match = which(czones_out %in% czones)
# # If more than one climate zone name was matched, set the outline color
# # to gray, representative of 'Multiple' applicable climate zones; otherwise
# # set to the point outline color appropriate to the matched climate zone
# if (length(czone_match)>1){
# results_finmets[(m - 1), 6] = "gray50"
# }else{
# results_finmets[(m - 1), 6] = czones_out_col[czone_match]
# }
# Determine appropriate ECM point shape for applicable building type
# Set ECM's applicable building type
bldg = results_database_filters$'Applicable Building Classes'
# Match applicable building classes to building type names used in plotting
bldg_match = matrix(NA, length(bclasses_out_finmets))
for (b in 1:length(bclasses_out_finmets)){
if (length(which(bldg%in%bclasses_out_finmets[[b]])>0)){
bldg_match[b] = b
}
}
# If more than one building type name was matched, set the point shape to
# a triangle, representative of 'Multiple' applicable building types;
# otherwise set to the point shape appropriate for the matched building type
if (length(bldg_match[is.finite(bldg_match)])>1){
results_finmets[(m - 1), 7] = 24
}else{
results_finmets[(m - 1), 7] = bclasses_out_finmets_shp[bldg_match[is.finite(bldg_match)]]
}
# Determine appropriate ECM point fill color for applicable end uses
# Set ECM's applicable end uses
euse = results_database_filters$'Applicable End Uses'
# Match applicable end uses to end use names used in plotting
euse_match = matrix(NA, length(euses_out_finmets))
for (e in 1:length(euses_out_finmets)){
if (length(which(euse%in% euses_out_finmets[[e]])>0)){
euse_match[e] = e
}
}
# If more than one end use name was matched, set the point fill color to
# gray, representative of 'Multiple' applicable end uses; otherwise set
# to the point fill color appropriate for the matched end use
if (length(euse_match[is.finite(euse_match)])>1){
results_finmets[(m - 1), 8] ="gray50"
}else{
results_finmets[(m - 1), 8] = euses_out_finmets_col[euse_match[is.finite(euse_match)]]
}
}
}
}
}
# Set the column start and stop indexes to use in updating the matrix of ECM
# energy, carbon or cost totals
col_ind_start = ((cp-1)*(length(comp_schemes))) + 1
col_ind_end = col_ind_start + 1 # note this accommodates baseline and efficient outcomes
# Update results matrix with mean, low, and high baseline and efficient outcomes
results[, col_ind_start:col_ind_end] = rbind(
cbind(list(r_temp[1,]), list(r_temp[4,])),
cbind(list(r_temp[2,]), list(r_temp[5,])),
cbind(list(r_temp[3,]), list(r_temp[6,])))
}
# Set uncompeted and competed baseline energy, carbon, or cost totals for given
# adoption scenario, plotting variable, and ECM (mean and low/high values
# for competed case)
base_uc = unlist(results[1, 1]) * unit_translate
base_c_m = unlist(results[1, 3]) * unit_translate
base_c_l = unlist(results[2, 3]) * unit_translate
base_c_h = unlist(results[3, 3]) * unit_translate
# Set uncompeted and competed efficient energy, carbon, or cost totals for
# adoption scenario, plotting variable, and ECM (mean and low/high values)
eff_uc_m = unlist(results[1, 2]) * unit_translate
eff_uc_l = unlist(results[2, 2]) * unit_translate
eff_uc_h = unlist(results[3, 2]) * unit_translate
eff_c_m = unlist(results[1, 4]) * unit_translate
eff_c_l = unlist(results[2, 4]) * unit_translate
eff_c_h = unlist(results[3, 4]) * unit_translate
# Add annual ECM energy, carbon, or cost totals to XLSX worksheet data frame
xlsx_data[row_ind_start:(row_ind_start + 3),
(6 + (length(plot_title_labels_finmets))):ncol(xlsx_data)] =
rbind(base_uc, eff_uc_m, base_c_m, eff_c_m)
# Find the min. and max. values in the ECM energy, carbon, or cost totals data
# to be plotted
min_val = min(c(base_uc, base_c_m, base_c_l, base_c_h,
eff_uc_m, eff_uc_l, eff_uc_h,
eff_c_m, eff_c_l, eff_c_h))
max_val = max(c(base_uc, base_c_m, base_c_l, base_c_h,
eff_uc_m, eff_uc_l, eff_uc_h,
eff_c_m, eff_c_l, eff_c_h))
# Set limits of y axis for plot based on min. and max. values in data
ylims = pretty(c(min_val-0.05*max_val, max_val+0.05*max_val))
# Initialize the plot with uncompeted baseline ECM energy, carbon, or cost totals
plot(years, base_uc, typ='l', lwd=5, ylim = c(min(ylims), max(ylims)),
xlab=NA, ylab=NA, col=plot_col_uc_base, main = meas_names[m], xaxt="n", yaxt="n")
# Determine legend parameters based on whether uncertainty is present in totals
if (uncertainty == TRUE){
# Set legend names for a plot with uncertainty in the ECM totals
legend_param = c(
"Baseline (Uncompeted)", "Baseline (Competed)",
"Efficient (Uncompeted)", "Efficient (Competed)",
"Baseline (Competed, 5th/95th PCT)", "Efficient (Uncompeted, 5th/95th PCT)",
"Efficient (Competed, 5th/95th PCT)")
col_param = c(plot_col_uc_base, plot_col_c_base, plot_col_uc_eff,
plot_col_c_eff, plot_col_c_base,
plot_col_uc_eff, plot_col_c_eff)
lwd_param = c(5, 3.5, 3, 2, rep(1, 3))
lty_param = c(rep(1, 4), rep(6, 3))
}else{
# Set legend names for a plot with no uncertainty in the ECM totals
legend_param = c("AEO Baseline (Uncompeted)", "AEO Baseline (Competed)",
"Efficient (Uncompeted)", "Efficient (Competed)")
col_param = c(plot_col_uc_base, plot_col_c_base, plot_col_uc_eff, plot_col_c_eff)
lwd_param = c(5, 3.5, 3, 2)
lty_param = rep(1, 4)
}
# Add low/high bounds on uncompeted and competed baseline and efficient
# ECM energy, carbon, or cost totals, if applicable
if (uncertainty == TRUE){
lines(years, eff_uc_l, lwd=1, lty=6, col=plot_col_uc_eff)
lines(years, eff_uc_h, lwd=1, lty=6, col=plot_col_uc_eff)
lines(years, base_c_l, lwd=1, lty=6 , col=plot_col_c_base)
lines(years, base_c_h, lwd=1, lty=6 , col=plot_col_c_base)
lines(years, eff_c_l, lwd=1, lty=6 , col=plot_col_c_eff)
lines(years, eff_c_h, lwd=1, lty=6 , col=plot_col_c_eff)
}
# Add mean uncompeted efficient ECM energy, carbon, or cost totals
lines(years, eff_uc_m, lwd=3, col=plot_col_uc_eff)
# Add mean competed baseline results
lines(years, base_c_m, lwd=3.5, col=plot_col_c_base)
# Add mean competed efficient results
lines(years, eff_c_m, lwd=2, col=plot_col_c_eff)
# Add x and y axis labels
mtext("Year", side=1, line=3.5, cex=0.925)
mtext(plot_axis_labels_ecm[v], side=2, line=3.65, cex=0.925)
# Add tick marks and labels to bottom and left axes
axis(side=1, at=pretty(c(min(years), max(years))),
labels=pretty(c(min(years), max(years))), cex.axis = 1.15)
axis(side=2, at=ylims, labels = ylims, cex.axis = 1.15, las=1)
# Add tick marks to top and right axes
axis(side=3, at=pretty(c(min(years), max(years))), labels = NA)
axis(side=4, at=ylims, labels = NA)
# Annotate total savings in a snapshot years for the 'All ECMs' case;
# otherwise, annotate the applicable ECM end uses
if (meas_names[m] == "All ECMs"){
# Annotate the plot with snapshot year total savings figure
# Find x and y values for annotation
for (s in 1:length(snap_yr_set)){
# Set the x value for annotation point
xval_snap = as.numeric(snap_yr_set[s])
# Set x value for annotation label
if (as.numeric(snap_yr_set[s])>2045){
# For later years, set the label to the left of the x val;
# otherwise, set the label under the x val
xlab_pos = 2
}else{
xlab_pos = 1
}
# Set y values for annotation point
yval_snap_eff = eff_c_m[which(years==snap_yr_set[s])]
yval_snap_base = base_c_m[which(years== snap_yr_set[s])]
# Draw line segment connecting snapshot year baseline and efficient results
points(xval_snap, yval_snap_base, col="forestgreen", pch = 1, cex=1.5, lwd=2.5)
points(xval_snap, yval_snap_eff, col="forestgreen", pch = 1, cex=1.5, lwd=2.5)
segments(xval_snap, yval_snap_eff, xval_snap, yval_snap_base, col="forestgreen", lty=3)
# Add snapshot year savings figure
text(xval_snap, yval_snap_eff - (yval_snap_eff - min(ylims))/7,
paste(toString(sprintf("%.1f", yval_snap_base-yval_snap_eff)),
toString(var_units[v]), sep=" "), pos = xlab_pos, col="forestgreen")
}
}else{
# Add ECM end use labels
text(min(years), max(ylims), labels=paste("End Uses: ", end_uses, sep=""), col="gray50",
pos=4, cex=0.93)
}
# If needed, add disclaimer about the use of public health costs
if (var_names_uncompete[v] == "cost" & phc_flag != FALSE){
# Set increment to use in plotting disclaimer text
ylim_inc = (max(ylims) - min(ylims))*0.05
if (phc_flag == "low"){
text(min(years), min(ylims), labels=
"* Reflects low EPA health benefits - ignore results after 2025",
col="red", pos=4, cex=0.85)
}else{
text(min(years), min(ylims), labels=