diff --git a/Documentation/.gitkeep b/Documentation/.gitkeep
new file mode 100644
index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391
diff --git a/Documentation/Data/.gitkeep b/Documentation/Data/.gitkeep
new file mode 100644
index 0000000000000000000000000000000000000000..e69de29bb2d1d6434b8b29ae775ad8c2e48c5391
diff --git a/Documentation/Data/BARREN.png b/Documentation/Data/BARREN.png
new file mode 100644
index 0000000000000000000000000000000000000000..b1199ec348ec606fd75ab6e21528213e10421d09
Binary files /dev/null and b/Documentation/Data/BARREN.png differ
diff --git a/Documentation/Data/Bassin_Map.png b/Documentation/Data/Bassin_Map.png
new file mode 100644
index 0000000000000000000000000000000000000000..3f6ff0775272516553f1b054b5e6cc91546a0686
Binary files /dev/null and b/Documentation/Data/Bassin_Map.png differ
diff --git a/Documentation/Data/CROPLAND.png b/Documentation/Data/CROPLAND.png
new file mode 100644
index 0000000000000000000000000000000000000000..90b247ebaa56bd964021053e435cfb70ba8b132b
Binary files /dev/null and b/Documentation/Data/CROPLAND.png differ
diff --git a/Documentation/Data/Cluster_Map.png b/Documentation/Data/Cluster_Map.png
new file mode 100644
index 0000000000000000000000000000000000000000..c68fcfaf35fc4cad0ac19800753f2d45e911af18
Binary files /dev/null and b/Documentation/Data/Cluster_Map.png differ
diff --git a/Documentation/Data/Commodities.csv b/Documentation/Data/Commodities.csv
new file mode 100644
index 0000000000000000000000000000000000000000..d67cf7f945304fb2601fcd9975a1fb68ca7d600c
--- /dev/null
+++ b/Documentation/Data/Commodities.csv
@@ -0,0 +1,11 @@
+Commodity,Description *not exhaustive list
+Cereal C3,Include wheat and barley
+Cereal C4,"Include maize, millet, sorghum"
+Oil crops,"Include soybeans, coconuts, olives, sunflower, rapeseed, oil palm fruit"
+Pulses,"Include beans, chickpeas, lentils, vetches, lupins, peas"
+Starchy roots,"Include potatoes, sweet potatoes, cassava, taro, yams, yautia"
+Sugar,"Include sugar cane, sugar beet, honey"
+Fruits and vegetables,"Include apples, pears, peaches, raspberries, courgette, butternut"
+Energy crops,Include miscanthus
+Ruminant livestock products,Include dairy and meat product
+Monogastric livestock products,Include meat products
diff --git a/Documentation/Data/Constraints.csv b/Documentation/Data/Constraints.csv
new file mode 100644
index 0000000000000000000000000000000000000000..db4ba4f9c3e0357647ec1f34ce1e672b02c6db39
--- /dev/null
+++ b/Documentation/Data/Constraints.csv
@@ -0,0 +1,4 @@
+Annual change limt,Variable name,Value
+Import,MAX_IMPORT_CHANGE,0.05 * time step
+Demand fraction,ANNUAL_MAX_DEMAND_FRACT_CHANGE,0.5
+Technology level,TECHNOLOGY_CHANGE_ANNUAL_RATE,0.002
diff --git a/Documentation/Data/Conutry_Map.png b/Documentation/Data/Conutry_Map.png
new file mode 100644
index 0000000000000000000000000000000000000000..75bf717e0cb96a5b00508d147e0efd8ce1c129ab
Binary files /dev/null and b/Documentation/Data/Conutry_Map.png differ
diff --git a/Documentation/Data/Country_Map.png b/Documentation/Data/Country_Map.png
new file mode 100644
index 0000000000000000000000000000000000000000..75bf717e0cb96a5b00508d147e0efd8ce1c129ab
Binary files /dev/null and b/Documentation/Data/Country_Map.png differ
diff --git a/Documentation/Data/Croptype.csv b/Documentation/Data/Croptype.csv
new file mode 100644
index 0000000000000000000000000000000000000000..020fa4fd059423b1446a5741e344d5e0cb851ab0
--- /dev/null
+++ b/Documentation/Data/Croptype.csv
@@ -0,0 +1,14 @@
+Crop Types
+Wheat
+Maize
+Rice
+Oil crops
+Pulses
+Starchy roots
+Energy crops
+Set aside
+Monogastrics
+Ruminants
+Fruits and vegetables
+Sugar
+Pasture
diff --git a/Documentation/Data/Eclispe_environment.png b/Documentation/Data/Eclispe_environment.png
new file mode 100644
index 0000000000000000000000000000000000000000..97c28573e1d68246e111472469fb502cf6cde49f
Binary files /dev/null and b/Documentation/Data/Eclispe_environment.png differ
diff --git a/Documentation/Data/Eddie_table1.csv b/Documentation/Data/Eddie_table1.csv
new file mode 100644
index 0000000000000000000000000000000000000000..b15764858c4c7f917148d13eaa0cd32c06740336
--- /dev/null
+++ b/Documentation/Data/Eddie_table1.csv
@@ -0,0 +1,6 @@
+ssp,protectionism,automation,diet,climate,cyber,financial
+SSP1,1,0,0,0,0,0
+SSP2,1,0,0,0,0,0
+SSP3,1,0,0,0,0,0
+SSP4,1,0,0,0,0,0
+SSP5,1,0,0,0,0,0
diff --git a/Documentation/Data/Indices.csv b/Documentation/Data/Indices.csv
new file mode 100644
index 0000000000000000000000000000000000000000..384202617664eb0601cfd67564453a0014921cb4
--- /dev/null
+++ b/Documentation/Data/Indices.csv
@@ -0,0 +1,10 @@
+Indices,Description
+i,spatial locations within a country
+j,"crop and pasture land use types (wheat, maize, rice, oilcrops, pulses, starchy roots, energy crops and pasture) "
+j-crops,"crop land use types, i.e. excluding pasture."
+k,"all agricultural commodities, including animal products"
+k-noncereals,"agricultural crop commodities excluding cereals, i.e. oilcrops, pulses, starchy roots and energy crops"
+k-cereals,"agricultural cereal commodities, i.e. wheat, maize and rice"
+k-feeds,"crops for use as livestock feed, i.e. wheat, maize and oilcrops"
+k-feeds_and_pasturefeed,feed crops plus pasture
+l,"land cover types (cropland, pasture and natural)"
diff --git a/Documentation/Data/Irrigation_constraint.png b/Documentation/Data/Irrigation_constraint.png
new file mode 100644
index 0000000000000000000000000000000000000000..4ad6829758c00faa4341f9d82af068d62e55703e
Binary files /dev/null and b/Documentation/Data/Irrigation_constraint.png differ
diff --git a/Documentation/Data/LandSyMM-modules.png b/Documentation/Data/LandSyMM-modules.png
new file mode 100644
index 0000000000000000000000000000000000000000..41c5a45f9b4443c46643df48ff14dc02172b10a8
Binary files /dev/null and b/Documentation/Data/LandSyMM-modules.png differ
diff --git a/Documentation/Data/LandSyMM_v01-PLUMv2_script_names.pdf b/Documentation/Data/LandSyMM_v01-PLUMv2_script_names.pdf
new file mode 100644
index 0000000000000000000000000000000000000000..56e821f3a152480de0334ae87c72f092f138d02c
Binary files /dev/null and b/Documentation/Data/LandSyMM_v01-PLUMv2_script_names.pdf differ
diff --git a/Documentation/Data/LandSyMM_v01-Page-1.png b/Documentation/Data/LandSyMM_v01-Page-1.png
new file mode 100644
index 0000000000000000000000000000000000000000..6ada3dae49b8f429a603b9ec317a1770efeb6878
Binary files /dev/null and b/Documentation/Data/LandSyMM_v01-Page-1.png differ
diff --git a/Documentation/Data/LandSyMM_v01-Very_simplefied_LandSYMM.png b/Documentation/Data/LandSyMM_v01-Very_simplefied_LandSYMM.png
new file mode 100644
index 0000000000000000000000000000000000000000..0c4fa1d411241b57ffffe28687f687828f813f53
Binary files /dev/null and b/Documentation/Data/LandSyMM_v01-Very_simplefied_LandSYMM.png differ
diff --git a/Documentation/Data/Landuse.csv b/Documentation/Data/Landuse.csv
new file mode 100644
index 0000000000000000000000000000000000000000..e6a98af2ffc1b2f693ccfa347c061d63a8a2a648
--- /dev/null
+++ b/Documentation/Data/Landuse.csv
@@ -0,0 +1,8 @@
+Land use and land cover
+Cropland
+Pasture
+Managed forest
+Unmanaged forest
+Urban
+Barren
+Other natural
diff --git a/Documentation/Data/MANAGED_FOREST.png b/Documentation/Data/MANAGED_FOREST.png
new file mode 100644
index 0000000000000000000000000000000000000000..9aa32360d99762e5ff1b68b594c99dbceba1a1c6
Binary files /dev/null and b/Documentation/Data/MANAGED_FOREST.png differ
diff --git a/Documentation/Data/NATURAL.png b/Documentation/Data/NATURAL.png
new file mode 100644
index 0000000000000000000000000000000000000000..f346d12370148728af2a5184708a85294acc600f
Binary files /dev/null and b/Documentation/Data/NATURAL.png differ
diff --git a/Documentation/Data/PASTURE.png b/Documentation/Data/PASTURE.png
new file mode 100644
index 0000000000000000000000000000000000000000..c5783f2a088263735601a764745cc5c00b88735c
Binary files /dev/null and b/Documentation/Data/PASTURE.png differ
diff --git a/Documentation/Data/Parameters.csv b/Documentation/Data/Parameters.csv
new file mode 100644
index 0000000000000000000000000000000000000000..a6d162cb02494aa4ba7264bf548f79f28fca4d51
--- /dev/null
+++ b/Documentation/Data/Parameters.csv
@@ -0,0 +1,17 @@
+Parameters,Description
+demandk,demand for commodity k (t)
+export_pck,export price for commodity k (\$/t)
+import_pck, import price for commodity k (\$/t) 
+base_lu_costj,base land use costs for crop j (\$/ha)
+f_cost, cost per tonne of fertiliser (\$/t)
+w_cost_indexi,"irrigation cost index at location i (unitless, 0-1 range) "
+w_cost,irrigation cost (\$/m2)
+m_costj,management intensity cost for crop j (\$)
+f_max,"maximum fertiliser rate, i.e. at f=1(t/ha)"
+"w_maxi,j","maximum irrigation rate, i.e. at w=1, for crop j at location i (m2/ha)"
+irrigiation_eff, Irrigation efficiency (%)
+crop_dmk,crop dry matter for commodity k (%)
+suitable_areai,area suitable for agriculture at location i (ha)
+water_availi,water availability at location i (m3/ha)
+min/max_net_importk, minimum and maximum net imports for commodity k (t)
+lc_change_unit_costl,land cover change unit cost for cover type l (\$/ha)
diff --git a/Documentation/Data/Parameters_example.csv b/Documentation/Data/Parameters_example.csv
new file mode 100644
index 0000000000000000000000000000000000000000..dc5bdc1ed7300673ec6edfddf27b132bfc97191b
--- /dev/null
+++ b/Documentation/Data/Parameters_example.csv
@@ -0,0 +1,18 @@
+Parameters,Central value
+"Irrigation cost, w_cost (\$/m2)",0.5
+"Fertiliser cost, f_cost (\$/t)",1800
+"Other intensity cost, m_cost (\$)",600
+"Land cover change cost, lc_change_cost: Natural to agricultural (\$/ha)",60
+"Land cover change cost, lc_change_cost: Managed forest to agricultural (\$/ha)",160
+"Land cover change cost, lc_change_cost: Agricultural land to natural land (\$/ha)",200
+"Land cover change cost, lc_change_cost: Pasture to cropland (\$/ha)",220
+"Land cover change cost, lc_change_cost: Cropland to pasture (\$/ha)",370
+Minimum natural or managed forest cover,10%
+Pasture harvest fraction,50%
+Seed and waste rate,10%
+"Technology yield change rate, ?",0.20%
+Intial price shift factor,1%
+"International market price sensitivity, ?",30%
+"International import tariff, i_tariff",20%
+"Transport cost, t_cost (\$/t)",50%
+"Transportation losses, t_cost",5%
diff --git a/Documentation/Data/Protected_area.png b/Documentation/Data/Protected_area.png
new file mode 100644
index 0000000000000000000000000000000000000000..c3ac5072714c616094300f47492e484902de750b
Binary files /dev/null and b/Documentation/Data/Protected_area.png differ
diff --git a/Documentation/Data/SSP.PNG b/Documentation/Data/SSP.PNG
new file mode 100644
index 0000000000000000000000000000000000000000..1bb3844457b08f1dc80065615c9266d203e428fc
Binary files /dev/null and b/Documentation/Data/SSP.PNG differ
diff --git a/Documentation/Data/TradeBarriers_map.png b/Documentation/Data/TradeBarriers_map.png
new file mode 100644
index 0000000000000000000000000000000000000000..95d68f4874d77428a73b60d29750fc4e775ef4f6
Binary files /dev/null and b/Documentation/Data/TradeBarriers_map.png differ
diff --git a/Documentation/Data/UNMANAGED_FOREST.png b/Documentation/Data/UNMANAGED_FOREST.png
new file mode 100644
index 0000000000000000000000000000000000000000..f63673cebdfa2a1edc5dfec59cff7eb52fa5ad99
Binary files /dev/null and b/Documentation/Data/UNMANAGED_FOREST.png differ
diff --git a/Documentation/Data/URBAN.png b/Documentation/Data/URBAN.png
new file mode 100644
index 0000000000000000000000000000000000000000..685149c924ad22ef8f70c373f561a5fb592afba4
Binary files /dev/null and b/Documentation/Data/URBAN.png differ
diff --git a/Documentation/Data/Variables.csv b/Documentation/Data/Variables.csv
new file mode 100644
index 0000000000000000000000000000000000000000..eff451a74d0490879d833b2774997698d97cf4f9
--- /dev/null
+++ b/Documentation/Data/Variables.csv
@@ -0,0 +1,9 @@
+Variables,Description
+"areai,j",area for crop j at location i (ha)
+"fi,j",fertiliser intensity factor for crop j at location i (unitless)
+"wi,j",irrigation intensity factor for crop j at location i (unitless)
+"mi,j",management intensity factor for crop j at location i (unitless)
+ruminant_feedk,ruminant feed use for commodity k (t)
+monogastric_feedk,monogastric feed use for commodity k (t)
+importk,imports for commodity k (t)
+exportk,exports for commodity k (t)
diff --git a/Documentation/Data/country_groups.csv b/Documentation/Data/country_groups.csv
new file mode 100644
index 0000000000000000000000000000000000000000..0856449476614629ffd3edd7587cf7cbd7ad8bf6
--- /dev/null
+++ b/Documentation/Data/country_groups.csv
@@ -0,0 +1,194 @@
+Country,population,Region,IncomeGroup,PlumGroup
+New Zealand,4.3507,East Asia & Pacific,High income: OECD,Australia & NZ
+Australia,22.03175,East Asia & Pacific,High income: OECD,Australia & NZ
+New Caledonia,,East Asia & Pacific,,Australia & NZ
+China,1337.705,East Asia & Pacific,Upper middle income,China
+Tuvalu,0.009827,East Asia & Pacific,Upper middle income,East Asia & Pacific_other
+Palau,0.02047,East Asia & Pacific,Upper middle income,East Asia & Pacific_other
+Marshall Islands,0.052428,East Asia & Pacific,Upper middle income,East Asia & Pacific_other
+Kiribati,0.102648,East Asia & Pacific,Lower middle income,East Asia & Pacific_other
+Micronesia (Federated States of),0.103619,East Asia & Pacific,Lower middle income,East Asia & Pacific_other
+Tonga,0.103947,East Asia & Pacific,Upper middle income,East Asia & Pacific_other
+Samoa,0.186029,East Asia & Pacific,Lower middle income,East Asia & Pacific_other
+Vanuatu,0.236299,East Asia & Pacific,Lower middle income,East Asia & Pacific_other
+Brunei Darussalam,0.393302,East Asia & Pacific,High income: nonOECD,East Asia & Pacific_other
+Solomon Islands,0.526177,East Asia & Pacific,Lower middle income,East Asia & Pacific_other
+Macau,0.534626,East Asia & Pacific,High income: nonOECD,East Asia & Pacific_other
+Fiji,0.859952,East Asia & Pacific,Upper middle income,East Asia & Pacific_other
+Lao People's Democratic Republic,6.260544,East Asia & Pacific,Lower middle income,East Asia & Pacific_other
+Papua New Guinea,6.847517,East Asia & Pacific,Lower middle income,East Asia & Pacific_other
+French Polynesia,,East Asia & Pacific,,East Asia & Pacific_other
+Democratic People's Republic of Korea,,East Asia & Pacific,,East Asia & Pacific_other
+Indonesia,241.613126,East Asia & Pacific,Lower middle income,Indonesia
+Japan,128.07,East Asia & Pacific,High income: OECD,Japan
+Mongolia,2.712657,East Asia & Pacific,Lower middle income,Mongolia
+Myanmar,51.733013,East Asia & Pacific,Low income,Myanmar
+Singapore,5.076732,East Asia & Pacific,High income: nonOECD,Philippines & Malaysia
+Malaysia,28.1195,East Asia & Pacific,Upper middle income,Philippines & Malaysia
+Philippines,93.038902,East Asia & Pacific,Lower middle income,Philippines & Malaysia
+Republic of Korea,49.410366,East Asia & Pacific,High income: OECD,Republic of Korea
+Thailand,66.692024,East Asia & Pacific,Upper middle income,Thailand
+Cambodia,14.363586,East Asia & Pacific,Low income,Viet Nam & Cambodia
+Viet Nam,86.9325,East Asia & Pacific,Lower middle income,Viet Nam & Cambodia
+Albania,2.913021,Europe & Central Asia,Upper middle income,Eastern Europe
+Armenia,2.963496,Europe & Central Asia,Lower middle income,Eastern Europe
+Slovakia,5.391428,Europe & Central Asia,High income: OECD,Eastern Europe
+Serbia,7.291436,Europe & Central Asia,Upper middle income,Eastern Europe
+Bulgaria,7.395599,Europe & Central Asia,Upper middle income,Eastern Europe
+Hungary,10.000023,Europe & Central Asia,Upper middle income,Eastern Europe
+Czechia,10.47441,Europe & Central Asia,High income: OECD,Eastern Europe
+Romania,20.246871,Europe & Central Asia,Upper middle income,Eastern Europe
+Montenegro,0.619428,Europe & Central Asia,Upper middle income,ex-Yugoslavia
+Kosovo,1.77568,Europe & Central Asia,Lower middle income,ex-Yugoslavia
+Slovenia,2.048583,Europe & Central Asia,High income: OECD,ex-Yugoslavia
+The former Yugoslav Republic of Macedonia,2.062443,Europe & Central Asia,Upper middle income,ex-Yugoslavia
+Republic of Moldova,3.562045,Europe & Central Asia,Lower middle income,ex-Yugoslavia
+Bosnia and Herzegovina,3.835258,Europe & Central Asia,Upper middle income,ex-Yugoslavia
+Croatia,4.417781,Europe & Central Asia,High income: nonOECD,ex-Yugoslavia
+Luxembourg,0.506953,Europe & Central Asia,High income: OECD,France Netherlands & Benlex
+Ireland,4.560155,Europe & Central Asia,High income: OECD,Ireland
+Belgium,10.895586,Europe & Central Asia,High income: OECD,France Netherlands & Benlex
+Netherlands,16.615394,Europe & Central Asia,High income: OECD,France Netherlands & Benlex
+France,65.027512,Europe & Central Asia,High income: OECD,France Netherlands & Benlex
+Switzerland,7.824909,Europe & Central Asia,High income: OECD,Germany Austria & Switzerland
+Austria,8.363404,Europe & Central Asia,High income: OECD,Germany Austria & Switzerland
+Germany,81.77693,Europe & Central Asia,High income: OECD,Germany Austria & Switzerland
+Italy,59.277417,Europe & Central Asia,High income: OECD,Italy
+Georgia,3.926,Europe & Central Asia,Lower middle income,Other former USSR
+Turkmenistan,5.041995,Europe & Central Asia,Upper middle income,Other former USSR
+Kyrgyzstan,5.4479,Europe & Central Asia,Low income,Other former USSR
+Tajikistan,7.581696,Europe & Central Asia,Low income,Other former USSR
+Azerbaijan,9.054332,Europe & Central Asia,Upper middle income,Other former USSR
+Belarus,9.49,Europe & Central Asia,Upper middle income,Other former USSR
+Kazakhstan,16.321581,Europe & Central Asia,Upper middle income,Other former USSR
+Uzbekistan,28.5624,Europe & Central Asia,Lower middle income,Other former USSR
+Poland,38.042794,Europe & Central Asia,High income: OECD,Poland
+Estonia,1.331475,Europe & Central Asia,High income: OECD,Scandinavia
+Latvia,2.097555,Europe & Central Asia,High income: nonOECD,Scandinavia
+Lithuania,3.097282,Europe & Central Asia,High income: nonOECD,Scandinavia
+Russian Federation,142.849449,Europe & Central Asia,High income: nonOECD,Russian Federation
+Greenland,0.056905,Europe & Central Asia,High income: nonOECD,Scandinavia
+Iceland,0.318041,Europe & Central Asia,High income: OECD,Scandinavia
+Norway,4.889252,Europe & Central Asia,High income: OECD,Scandinavia
+Finland,5.363352,Europe & Central Asia,High income: OECD,Scandinavia
+Denmark,5.547683,Europe & Central Asia,High income: OECD,Scandinavia
+Sweden,9.378126,Europe & Central Asia,High income: OECD,Scandinavia
+Portugal,10.5731,Europe & Central Asia,High income: OECD,Spain & Portugal
+Spain,46.576897,Europe & Central Asia,High income: OECD,Spain & Portugal
+Cyprus,1.103685,Europe & Central Asia,High income: nonOECD,Turkey Cyprus & Greece
+Greece,11.121341,Europe & Central Asia,High income: OECD,Turkey Cyprus & Greece
+Turkey,72.310416,Europe & Central Asia,Upper middle income,Turkey Cyprus & Greece
+Ukraine,45.8707,Europe & Central Asia,Lower middle income,Ukraine
+United Kingdom,62.766365,Europe & Central Asia,High income: OECD,United Kingdom
+Argentina,41.222875,Latin America & Caribbean,Upper middle income,Argentina
+Bolivia (Plurinational State of),9.918245,Latin America & Caribbean,Lower middle income,Bolivia and Chile
+Chile,17.015048,Latin America & Caribbean,High income: OECD,Bolivia and Chile
+Brazil,198.614208,Latin America & Caribbean,Upper middle income,Brazil
+Saint Kitts and Nevis,0.052352,Latin America & Caribbean,High income: nonOECD,Caribbean
+Dominica,0.071167,Latin America & Caribbean,Upper middle income,Caribbean
+Antigua and Barbuda,0.087233,Latin America & Caribbean,High income: nonOECD,Caribbean
+Aruba,0.101597,Latin America & Caribbean,High income: nonOECD,Caribbean
+Grenada,0.104677,Latin America & Caribbean,Upper middle income,Caribbean
+Saint Vincent and the Grenadines,0.109316,Latin America & Caribbean,Upper middle income,Caribbean
+Saint Lucia,0.177397,Latin America & Caribbean,Upper middle income,Caribbean
+Barbados,0.279566,Latin America & Caribbean,High income: nonOECD,Caribbean
+Bahamas,0.36083,Latin America & Caribbean,High income: nonOECD,Caribbean
+Trinidad and Tobago,1.328095,Latin America & Caribbean,High income: nonOECD,Caribbean
+Jamaica,2.690824,Latin America & Caribbean,Upper middle income,Caribbean
+Puerto Rico,3.721526,Latin America & Caribbean,High income: nonOECD,Caribbean
+Dominican Republic,9.897983,Latin America & Caribbean,Upper middle income,Caribbean
+Haiti,9.999617,Latin America & Caribbean,Low income,Caribbean
+Cuba,11.308133,Latin America & Caribbean,Upper middle income,Caribbean
+Bermuda,,Latin America & Caribbean,,Caribbean
+Belize,0.321609,Latin America & Caribbean,Upper middle income,Central America
+Panama,3.620506,Latin America & Caribbean,Upper middle income,Central America
+Costa Rica,4.545273,Latin America & Caribbean,Upper middle income,Central America
+Nicaragua,5.737722,Latin America & Caribbean,Lower middle income,Central America
+El Salvador,6.038306,Latin America & Caribbean,Lower middle income,Central America
+Honduras,7.503875,Latin America & Caribbean,Lower middle income,Central America
+Guatemala,14.732261,Latin America & Caribbean,Lower middle income,Central America
+Colombia,45.918101,Latin America & Caribbean,Upper middle income,Colombia
+Mexico,118.617542,Latin America & Caribbean,Upper middle income,Mexico
+Uruguay,3.374414,Latin America & Caribbean,High income: nonOECD,Paraguay & Uruguay
+Paraguay,6.209877,Latin America & Caribbean,Lower middle income,Paraguay & Uruguay
+Ecuador,14.934692,Latin America & Caribbean,Upper middle income,Peru & Ecuador
+Peru,29.373644,Latin America & Caribbean,Upper middle income,Peru & Ecuador
+Suriname,0.518141,Latin America & Caribbean,Upper middle income,Venezuela Guyana & Suriname
+Guyana,0.753362,Latin America & Caribbean,Lower middle income,Venezuela Guyana & Suriname
+Venezuela (Bolivarian Republic of),28.995745,Latin America & Caribbean,Upper middle income,Venezuela Guyana & Suriname
+Algeria,36.036159,Middle East & North Africa,Upper middle income,Algeria
+Egypt,82.040994,Middle East & North Africa,Lower middle income,Egypt
+Iran (Islamic Republic of),74.253373,Middle East & North Africa,Upper middle income,Iran (Islamic Republic of)
+Iraq,30.868156,Middle East & North Africa,Upper middle income,Iraq
+Bahrain,1.261319,Middle East & North Africa,High income: nonOECD,Middle East other
+Qatar,1.765513,Middle East & North Africa,High income: nonOECD,Middle East other
+Oman,2.943747,Middle East & North Africa,High income: nonOECD,Middle East other
+Kuwait,3.059473,Middle East & North Africa,High income: nonOECD,Middle East other
+Lebanon,4.337156,Middle East & North Africa,Upper middle income,Middle East other
+Jordan,6.517912,Middle East & North Africa,Upper middle income,Middle East other
+Israel,7.6236,Middle East & North Africa,High income: OECD,Middle East other
+United Arab Emirates,8.329453,Middle East & North Africa,High income: nonOECD,Middle East other
+Yemen,23.591972,Middle East & North Africa,Lower middle income,Middle East other
+Saudi Arabia,28.090647,Middle East & North Africa,High income: nonOECD,Middle East other
+Morocco,32.107739,Middle East & North Africa,Lower middle income,Morocco
+Malta,0.414508,Europe & Central Asia,High income: nonOECD,Turkey Cyprus & Greece
+Djibouti,0.830802,Sub-Saharan Africa,Lower middle income,East Africa
+Libya,6.265697,Middle East & North Africa,Upper middle income,North Africa other
+Tunisia,10.5471,Middle East & North Africa,Upper middle income,North Africa other
+Canada,34.005274,North America,High income: OECD,Canada
+United States of America,309.346863,North America,High income: OECD,United States of America
+Afghanistan,27.962207,South Asia,Low income,Pakistan & Afghanistan
+Bangladesh,151.616777,South Asia,Low income,Bangladesh
+Maldives,0.367,South Asia,Upper middle income,India  & Sri Lanka
+Sri Lanka,20.118,South Asia,Lower middle income,India  & Sri Lanka
+India,1230.984504,South Asia,Lower middle income,India  & Sri Lanka
+Bhutan,0.720246,South Asia,Lower middle income,Nepal & Butan
+Nepal,26.87591,South Asia,Low income,Nepal & Butan
+Pakistan,170.043918,South Asia,Lower middle income,Pakistan & Afghanistan
+Sao Tome and Principe,0.17088,Sub-Saharan Africa,Lower middle income,Central Africa
+Comoros,0.698695,Sub-Saharan Africa,Low income,Central Africa
+Equatorial Guinea,0.72871,Sub-Saharan Africa,High income: nonOECD,Central Africa
+Congo,4.066078,Sub-Saharan Africa,Lower middle income,Central Africa
+Central African Republic,4.444973,Sub-Saharan Africa,Low income,Central Africa
+Burundi,9.461117,Sub-Saharan Africa,Low income,Central Africa
+Rwanda,10.293669,Sub-Saharan Africa,Low income,Central Africa
+Chad,11.89638,Sub-Saharan Africa,Low income,Central Africa
+Cameroon,20.590666,Sub-Saharan Africa,Lower middle income,Central Africa
+Democratic Republic of the Congo,65.938712,Sub-Saharan Africa,Low income,Democratic Republic of the Congo
+Eritrea,4.689664,Sub-Saharan Africa,Low income,East Africa
+South Sudan,10.056475,Sub-Saharan Africa,Lower middle income,East Africa
+Ethiopia,87.561814,Sub-Saharan Africa,Low income,Ethiopia
+Kenya,40.328313,Sub-Saharan Africa,Low income,Kenya
+Nigeria,159.424742,Sub-Saharan Africa,Lower middle income,Nigeria
+South Africa,50.771826,Sub-Saharan Africa,Upper middle income,South Africa
+Swaziland,1.193148,Sub-Saharan Africa,Lower middle income,Southern Africa other
+Mauritius,1.2504,Sub-Saharan Africa,Upper middle income,Southern Africa other
+Lesotho,2.010586,Sub-Saharan Africa,Lower middle income,Southern Africa other
+Botswana,2.047831,Sub-Saharan Africa,Upper middle income,Southern Africa other
+Namibia,2.193643,Sub-Saharan Africa,Upper middle income,Southern Africa other
+Zambia,13.917439,Sub-Saharan Africa,Lower middle income,Southern Africa other
+Zimbabwe,13.973897,Sub-Saharan Africa,Low income,Southern Africa other
+Malawi,14.769824,Sub-Saharan Africa,Low income,Southern Africa other
+Madagascar,21.079532,Sub-Saharan Africa,Low income,Southern Africa other
+Angola,21.219954,Sub-Saharan Africa,Upper middle income,Central Africa
+Mozambique,24.321457,Sub-Saharan Africa,Low income,Southern Africa other
+Sudan,36.114885,Sub-Saharan Africa,Lower middle income,Sudan
+Uganda,33.149417,Sub-Saharan Africa,Low income,Uganda
+United Republic of Tanzania,45.648525,Sub-Saharan Africa,Low income,United Republic of Tanzania
+Cabo Verde,0.490379,Sub-Saharan Africa,Lower middle income,West Africa
+Gabon,1.541936,Sub-Saharan Africa,Upper middle income,Central Africa
+Guinea-Bissau,1.634196,Sub-Saharan Africa,Low income,West Africa
+Gambia,1.693002,Sub-Saharan Africa,Low income,West Africa
+Mauritania,3.5914,Sub-Saharan Africa,Lower middle income,West Africa
+Liberia,3.95799,Sub-Saharan Africa,Low income,West Africa
+Sierra Leone,5.775902,Sub-Saharan Africa,Low income,West Africa
+Togo,6.390851,Sub-Saharan Africa,Low income,West Africa
+Benin,9.509798,Sub-Saharan Africa,Low income,West Africa
+Guinea,11.012406,Sub-Saharan Africa,Low income,West Africa
+Senegal,12.956791,Sub-Saharan Africa,Lower middle income,West Africa
+Mali,15.167286,Sub-Saharan Africa,Low income,West Africa
+Burkina Faso,15.632066,Sub-Saharan Africa,Low income,West Africa
+Niger,16.29199,Sub-Saharan Africa,Low income,West Africa
+Cote d'Ivoire,20.131707,Sub-Saharan Africa,Lower middle income,West Africa
+Ghana,24.317734,Sub-Saharan Africa,Lower middle income,West Africa
+Somalia,,Sub-Saharan Africa,,East Africa
diff --git a/Documentation/Data/shockProbabilities.csv b/Documentation/Data/shockProbabilities.csv
new file mode 100644
index 0000000000000000000000000000000000000000..3f1bb3d086121a4a8039fd8a750d21070d551098
--- /dev/null
+++ b/Documentation/Data/shockProbabilities.csv
@@ -0,0 +1,6 @@
+ssp,protectionism,automation,diet,climate,cyber,financial
+SSP1,0,1,0,0,0,0
+SSP2,0,1,0,0,0,0
+SSP3,0,1,0,0,0,0
+SSP4,0,1,0,0,0,0
+SSP5,0,1,0,0,0,0
diff --git a/Documentation/Documentation.Rmd b/Documentation/Documentation.Rmd
index 3bc44db5f383a86eab50362b6c9aced309f41f9c..e3765481081644fe558046d51384042d4723803b 100644
--- a/Documentation/Documentation.Rmd
+++ b/Documentation/Documentation.Rmd
@@ -28,7 +28,7 @@ PLUMv2 was developed to inform policy making about interaction between policy me
 
 PLUM was developed by Kerstin Engström team and collaborators in 2016. It is open-source (except GAMS license) and mainly implemented in java and the post processing is usually implemented in the R software.
 ```{r, animation.hook='gifski', width = 2000, height = 1600, delay = 0.1, progress = TRUE, loop = FALSE, echo=FALSE, message=FALSE, warning = FALSE}
-output_path<-"C:/Users/jmaire/Development/plumv2/output/Baseline_test/"
+output_path<-"C:/Users/jmaire/Development/plumv2_juliette_local/output/Baseline_test/"
   for (i in 2011:2012){
     test<-fread((paste0(output_path,i,"/LandCoverFract.txt")))
     test<-test %>%
@@ -491,6 +491,7 @@ $$
 
 ### *2019* {#DEV20190}
 
+
 * Addition of DEMAND_PRICE_IMPORT_AND_PROD_COST 
 
 # **INSTALLATION**{.tabset .tabset-fade .tabset-pills}
diff --git a/Documentation/Documentation.html b/Documentation/Documentation.html
index 3bfabbf2217729bc01eb8e5c3ff15cf6d364585b..9beb6be74384e447eea9de78c87e294ae124670c 100644
--- a/Documentation/Documentation.html
+++ b/Documentation/Documentation.html
@@ -1290,6 +1290,192 @@ window.buildTabsets = function(tocID) {
     }
 });
 </script>
+<style type="text/css">
+.lightable-minimal {
+border-collapse: separate;
+border-spacing: 16px 1px;
+width: 100%;
+margin-bottom: 10px;
+}
+.lightable-minimal td {
+margin-left: 5px;
+margin-right: 5px;
+}
+.lightable-minimal th {
+margin-left: 5px;
+margin-right: 5px;
+}
+.lightable-minimal thead tr:last-child th {
+border-bottom: 2px solid #00000050;
+empty-cells: hide;
+}
+.lightable-minimal tbody tr:first-child td {
+padding-top: 0.5em;
+}
+.lightable-minimal.lightable-hover tbody tr:hover {
+background-color: #f5f5f5;
+}
+.lightable-minimal.lightable-striped tbody tr:nth-child(even) {
+background-color: #f5f5f5;
+}
+.lightable-classic {
+border-top: 0.16em solid #111111;
+width: 100%;
+margin-bottom: 10px;
+margin: 10px 5px;
+}
+.lightable-classic caption {
+color: #222222;
+}
+.lightable-classic td {
+padding-left: 5px;
+padding-right: 5px;
+color: #222222;
+}
+.lightable-classic th {
+padding-left: 5px;
+padding-right: 5px;
+font-weight: normal;
+color: #222222;
+}
+.lightable-classic thead tr:last-child th {
+border-bottom: 0.10em solid #111111;
+}
+.lightable-classic tbody tr:last-child td {
+border-bottom: 0.14em solid #111111;
+}
+.lightable-classic.lightable-hover tbody tr:hover {
+background-color: #F9EEC1;
+}
+.lightable-classic.lightable-striped tbody tr:nth-child(even) {
+background-color: #f5f5f5;
+}
+.lightable-classic-2 {
+border-top: 3px double #111111;
+width: 100%;
+margin-bottom: 10px;
+}
+.lightable-classic-2 caption {
+color: #222222;
+}
+.lightable-classic-2 td {
+padding-left: 5px;
+padding-right: 5px;
+color: #222222;
+}
+.lightable-classic-2 th {
+padding-left: 5px;
+padding-right: 5px;
+font-weight: normal;
+color: #222222;
+}
+.lightable-classic-2 tbody tr:last-child td {
+border-bottom: 3px double #111111;
+}
+.lightable-classic-2 thead tr:last-child th {
+border-bottom: 1px solid #111111;
+}
+.lightable-classic-2.lightable-hover tbody tr:hover {
+background-color: #F9EEC1;
+}
+.lightable-classic-2.lightable-striped tbody tr:nth-child(even) {
+background-color: #f5f5f5;
+}
+.lightable-material {
+min-width: 100%;
+white-space: nowrap;
+table-layout: fixed;
+font-family: Roboto, sans-serif;
+border: 1px solid #EEE;
+border-collapse: collapse;
+margin-bottom: 10px;
+}
+.lightable-material th {
+height: 56px;
+padding-left: 16px;
+padding-right: 16px;
+}
+.lightable-material td {
+height: 52px;
+padding-left: 16px;
+padding-right: 16px;
+border-top: 1px solid #eeeeee;
+}
+.lightable-material.lightable-hover tbody tr:hover {
+background-color: #f5f5f5;
+}
+.lightable-material.lightable-striped tbody tr:nth-child(even) {
+background-color: #f5f5f5;
+}
+.lightable-material.lightable-striped tbody td {
+border: 0;
+}
+.lightable-material.lightable-striped thead tr:last-child th {
+border-bottom: 1px solid #ddd;
+}
+.lightable-material-dark {
+min-width: 100%;
+white-space: nowrap;
+table-layout: fixed;
+font-family: Roboto, sans-serif;
+border: 1px solid #FFFFFF12;
+border-collapse: collapse;
+margin-bottom: 10px;
+background-color: #363640;
+}
+.lightable-material-dark th {
+height: 56px;
+padding-left: 16px;
+padding-right: 16px;
+color: #FFFFFF60;
+}
+.lightable-material-dark td {
+height: 52px;
+padding-left: 16px;
+padding-right: 16px;
+color: #FFFFFF;
+border-top: 1px solid #FFFFFF12;
+}
+.lightable-material-dark.lightable-hover tbody tr:hover {
+background-color: #FFFFFF12;
+}
+.lightable-material-dark.lightable-striped tbody tr:nth-child(even) {
+background-color: #FFFFFF12;
+}
+.lightable-material-dark.lightable-striped tbody td {
+border: 0;
+}
+.lightable-material-dark.lightable-striped thead tr:last-child th {
+border-bottom: 1px solid #FFFFFF12;
+}
+.lightable-paper {
+width: 100%;
+margin-bottom: 10px;
+color: #444;
+}
+.lightable-paper thead tr:last-child th {
+color: #666;
+vertical-align: bottom;
+border-bottom: 1px solid #00000020;
+line-height: 1.15em;
+padding: 10px 5px;
+}
+.lightable-paper td {
+vertical-align: middle;
+border-bottom: 1px solid #00000010;
+line-height: 1.15em;
+padding: 7px 5px;
+}
+.lightable-paper.lightable-hover tbody tr:hover {
+background-color: #F9EEC1;
+}
+.lightable-paper.lightable-striped tbody tr:nth-child(even) {
+background-color: #00000008;
+}
+.lightable-paper.lightable-striped tbody td {
+border: 0;
+}
+</style>
 
 
 <style type="text/css">code{white-space: pre;}</style>
diff --git a/GAMS/IntExtOpt.gms b/GAMS/IntExtOpt.gms
index 4e157e7f7929095687e534d7ada269ca6b7d4d18..6b558131cfd3a70f3fcccec78b30cb27176729b8 100644
--- a/GAMS/IntExtOpt.gms
+++ b/GAMS/IntExtOpt.gms
@@ -45,7 +45,6 @@
  SCALAR pastureDecCost                       price for decreasing pasture;
   
  SCALAR meatEfficency                        efficiency of converting feed and pasture into animal products;
- SCALAR animalFeedFromOtherSources           animal feed from source other than modelled crops;
  SCALAR fertiliserUnitCost                   fert cost at max fert rate;
  SCALAR otherIParam                          yield response to other intensity;
  SCALAR otherICost                           cost of other intensity;
@@ -60,15 +59,11 @@ $load location, suitableLandArea, demand, agriExpansionCost, cropIncCost, pastur
 $load previousArea, previousFertIntensity, previousIrrigIntensity, previousOtherIntensity, previousRuminantFeed, previousMonogastricFeed, previousImportAmount, previousExportAmount
 $load yieldNone, yieldFertOnly, yieldIrrigOnly, yieldBoth, yieldShock
 $load fertParam, irrigParam, otherIParam, exportPrices, importPrices, maxNetImport, minNetImport, unhandledCropRate, setAsideRate, maxLandExpansionRate, subsidyRate
-$load meatEfficency, otherICost, irrigCost, irrigMaxRate, irrigConstraint, fertiliserUnitCost, domesticPriceMarkup, minDemandPerCereal, minDemandPerOilcrop, seedAndWasteRate, animalFeedFromOtherSources
+$load meatEfficency, otherICost, irrigCost, irrigMaxRate, irrigConstraint, fertiliserUnitCost, domesticPriceMarkup, minDemandPerCereal, minDemandPerOilcrop, seedAndWasteRate
 $gdxin    
  
  SCALAR delta use to smooth power function see 7.5 www.gams.com dd docs solversconopt.pdf / 0.00000000001 /;
  
- SCALAR ruminantOtherFeed;
- SCALAR monogastricOtherFeed;
- ruminantOtherFeed = animalFeedFromOtherSources * 0.25;
- monogastricOtherFeed = animalFeedFromOtherSources * 0.75;
  demand(cereal_crop) = demand('cereals') * minDemandPerCereal(cereal_crop);
  demand(oilpulse_crop) = demand('oilcropspulses') * minDemandPerOilcrop(oilpulse_crop);
  
@@ -177,9 +172,9 @@ $gdxin
  
  TOTAL_OIL_PULSE_DEMAND_CONSTRAINT .. sum(oilpulse_crop, net_supply(oilpulse_crop)) =G= demand('oilcropspulses');
  
- RUMINANT_DEMAND_CONSTRAINT .. meatEfficency*(sum(feed_crop, ruminantFeed(feed_crop) * cropDM(feed_crop)) + ruminantOtherFeed) * (1 - seedAndWasteRate('ruminants')) =G= (demand('ruminants') - importAmount('ruminants') + exportAmount('ruminants'));
+ RUMINANT_DEMAND_CONSTRAINT .. meatEfficency * sum(feed_crop, ruminantFeed(feed_crop) * cropDM(feed_crop)) * (1 - seedAndWasteRate('ruminants')) =G= (demand('ruminants') - importAmount('ruminants') + exportAmount('ruminants'));
     
- MONOGASTRICS_DEMAND_CONSTRAINT .. meatEfficency*(sum(feed_crop_less_pasture, monogastricFeed(feed_crop_less_pasture) * cropDM(feed_crop_less_pasture)) + monogastricOtherFeed) * (1 - seedAndWasteRate('monogastrics')) =G= (demand('monogastrics') - importAmount('monogastrics') + exportAmount('monogastrics'));
+ MONOGASTRICS_DEMAND_CONSTRAINT .. meatEfficency * sum(feed_crop_less_pasture, monogastricFeed(feed_crop_less_pasture) * cropDM(feed_crop_less_pasture)) * (1 - seedAndWasteRate('monogastrics')) =G= (demand('monogastrics') - importAmount('monogastrics') + exportAmount('monogastrics'));
     
  TOTAL_NON_PASTURE_FEED_DM_CALC .. totalFeedDM =E= sum(feed_crop_less_pasture, (ruminantFeed(feed_crop_less_pasture) + monogastricFeed(feed_crop_less_pasture)) * cropDM(feed_crop_less_pasture));
  FEED_MIX_CONSTRAINT(feed_crop_less_pasture) .. (ruminantFeed(feed_crop_less_pasture) + monogastricFeed(feed_crop_less_pasture)) * cropDM(feed_crop_less_pasture) =L= totalFeedDM * 0.7;
@@ -254,8 +249,8 @@ $gdxin
 * Production quantities based on smaller area (before unhandledCropArea adjustment applied)
  totalProd(crop) = sum(location, area.l(crop, location) * yield.l(crop, location));
  productionShock(crop) = sum(location, area.l(crop, location) * yield.l(crop, location) * yieldShock(crop, location));
- totalProd('ruminants') = meatEfficency*(sum(feed_crop, ruminantFeed.l(feed_crop) * cropDM(feed_crop)) + ruminantOtherFeed);
- totalProd('monogastrics') = meatEfficency*(sum(feed_crop, monogastricFeed.l(feed_crop) * cropDM(feed_crop)) + monogastricOtherFeed);
+ totalProd('ruminants') = meatEfficency*(sum(feed_crop, ruminantFeed.l(feed_crop) * cropDM(feed_crop)));
+ totalProd('monogastrics') = meatEfficency*(sum(feed_crop, monogastricFeed.l(feed_crop) * cropDM(feed_crop)));
  
 * Cost based on adjusted area
  area.l(crop_less_pasture, location) = area.l(crop_less_pasture, location) / (1.0 - unhandledCropRate);
diff --git a/data/shockProbabilities.csv b/data/shockProbabilities.csv
index 3ce0dffa51084a3b284190d5bd98564d01a0bc8d..251b66199bf30d31f755bc49d046df582ba844bb 100644
--- a/data/shockProbabilities.csv
+++ b/data/shockProbabilities.csv
@@ -1,6 +1,6 @@
 ssp,protectionism,automation,diet,climate,cyber,financial
-SSP1,0,1,0,0,0,0
-SSP2,0,1,0,0,0,0
-SSP3,0,1,0,0,0,0
-SSP4,0,1,0,0,0,0
-SSP5,0,1,0,0,0,0
+SSP1,0,0,1,0,0,0
+SSP2,0,0,1,0,0,0
+SSP3,0,0,1,0,0,0
+SSP4,0,0,1,0,0,0
+SSP5,0,0,1,0,0,0
diff --git a/debug_config.properties b/debug_config.properties
index 717dffaff034bd4bac6453e69eb9d95c661faabf..6f97253c832cd0d8951ab8fb1e968b6b2fcf3d2f 100644
--- a/debug_config.properties
+++ b/debug_config.properties
@@ -1,30 +1,29 @@
 BASE_DIR=..
 
-YIELD_DIR=/Users/peteralexander/Documents/LURG/LPJ/LPJGPLUM_remap6p7_20190225/rcp26
-YIELD_FILENAME=yield.out
+YIELD_DIR=/Users/palexand/Documents/LURG/LPJ/CMIP_emulated_test_forPLUM_20200821180759/GFDL-ESM4_ssp126_v20200821_rmol/sim_LPJ-GUESS
 
-DEBUG_LIMIT_COUNTRIES=false
-DEBUG_COUNTRY_NAME=Mongolia
+END_TIMESTEP=0
 
-IS_CALIBRATION_RUN = true
-GENERATE_NEW_YIELD_CLUSTERS = false
-NUM_YIELD_CLUSTERS=8000
+END_FIRST_STAGE_CALIBRATION=10
+TIMESTEP_SIZE=1
 
-END_TIMESTEP=0
-TIMESTEP_SIZE=10
+IS_CALIBRATION_RUN=true
 
-INTERPOLATE_OUTPUT_YEARS = false
+GENERATE_NEW_YIELD_CLUSTERS=false
 
-CHANGE_YIELD_DATA_YEAR=false
+DEMAND_PRICE_IMPORT_AND_PROD_COST=true
 
-ORIG_LEAST_COST_MIN=true
+YIELD_FILENAME=yield_intpinfs_rmol.out
+IRRIG_MAX_WATER_FILENAME=gsirrigation_intpinfs_rmol.out
 
-DEBUG_JUST_DEMAND_OUTPUT=false
-PRICE_ELASTIC_DEMAND=true
+MIN_FERT_AMOUNT=10
+MID_FERT_AMOUNT=60
+MAX_FERT_AMOUNT=200
+FERT_AMOUNT_PADDING=3
 
-GAMS_COUNTRY_TO_SAVE=Brazil
+LPJG_TIMESTEP_SIZE=10
 
-SSP_SCENARIO=SSP2_v9_130325
-DONT_REBASE_DEMAND=false
+EXTRAPOLATE_YIELD_FERT_RESPONSE=false
 
-DEMAND_RECALC_MAX_ITERATIONS=0
\ No newline at end of file
+DEBUG_LIMIT_COUNTRIES=true
+DEBUG_COUNTRY_NAME=Italy
\ No newline at end of file
diff --git a/src/ac/ed/lurg/ModelConfig.java b/src/ac/ed/lurg/ModelConfig.java
index 015b641d738dea8eaed3bd4a4bf41c074d7e0790..9039bcfd75852521cc5a06453d07afc9af42b98c 100755
--- a/src/ac/ed/lurg/ModelConfig.java
+++ b/src/ac/ed/lurg/ModelConfig.java
@@ -304,7 +304,6 @@ public class ModelConfig {
 	public static final double MAX_INCOME_PROPORTION_FOOD_SPEND = getDoubleProperty("MAX_INCOME_PROPORTION_FOOD_SPEND", 0.7);
 	
 	public static final double PASTURE_HARVEST_FRACTION = getDoubleProperty("PASTURE_HARVEST_FRACTION", 0.5);
-	public static final double ANIMAL_FEED_FROM_OTHER_SOURCES_RATE = getDoubleProperty("ANIMAL_FEED_FROM_OTHER_SOURCES_RATE", 0.127);  // animal nutrition coming from sources other crops modelled
 	public static final double MEAT_EFFICIENCY = getDoubleProperty("MEAT_EFFICIENCY", 1.0);  // 'meat' is includes feed conversion ratio already, this is tech. change or similar
 	public static final double IRRIGIATION_EFFICIENCY = getDoubleProperty("IRRIGIATION_EFFICIENCY", 0.5);
 	public static final int ELLIOTT_BASEYEAR = 2010;
@@ -394,7 +393,7 @@ public class ModelConfig {
 	public static final String GAMS_COUNTRY_TO_SAVE = getProperty("GAMS_COUNTRY_TO_SAVE", "China");;
 
 	public static final double PASTURE_MAX_IRRIGATION_RATE = getDoubleProperty("DEFAULT_MAX_IRRIGATION_RATE", 50.0); // shouldn't need this but some areas crops don't have a value, but was causing them to be selected
-	public static final int LPJG_TIMESTEP_SIZE = 5;
+	public static final int LPJG_TIMESTEP_SIZE = getIntProperty("LPJG_TIMESTEP_SIZE", 5);
 	
 	public static final int NUM_YIELD_CLUSTERS = getIntProperty("NUM_YIELD_CLUSTERS", 8000);
 	public static final long RANDOM_SEED = getIntProperty("RANDOM_SEED", 1974329);  // any number will do
@@ -408,4 +407,6 @@ public class ModelConfig {
 
 	public static final boolean USE_CRAFTY_COUNTRIES = getBooleanProperty("USE_CRAFTY_COUNTRIES", false);
 	public static final String CRAFTY_PRODUCTION_DIR = getProperty("CRAFTY_PRODUCTION_DIR", OUTPUT_DIR + File.separator + "crafty");
+	
+	public static final boolean EXTRAPOLATE_YIELD_FERT_RESPONSE = getBooleanProperty("EXTRAPOLATE_YIELD_FERT_RESPONSE", false);
 }
diff --git a/src/ac/ed/lurg/ModelMain.java b/src/ac/ed/lurg/ModelMain.java
index 6c661ca005d00fd7dffaa0f2d44038e7188fa2be..070cfc7dfc5e71f8999adf19d3adbf90e6fdae7c 100644
--- a/src/ac/ed/lurg/ModelMain.java
+++ b/src/ac/ed/lurg/ModelMain.java
@@ -74,7 +74,7 @@ public class ModelMain {
 
 	public static void main(String[] args) {
 		ModelMain theModel = new ModelMain();
-	    	System.out.println("Working Directory = " + System.getProperty("user.dir"));
+	    System.out.println("Working Directory = " + System.getProperty("user.dir"));
 		theModel.setup();
 		theModel.run();
 	}
diff --git a/src/ac/ed/lurg/country/CountryAgent.java b/src/ac/ed/lurg/country/CountryAgent.java
index d133224b6913b7a3408753b2378f767327009e00..6aa7b277f633055ae16745635f9494792ae428a2 100644
--- a/src/ac/ed/lurg/country/CountryAgent.java
+++ b/src/ac/ed/lurg/country/CountryAgent.java
@@ -38,7 +38,6 @@ public class CountryAgent extends AbstractCountryAgent {
 
 	private GamsRasterOutput previousGamsRasterOutput;
 	private RasterSet<IntegerRasterItem> yieldClusters;
-	private double animalFeedFromOtherSources;
 	private Map<CropType, Double> subsidyRates;
 	private boolean saveGamsGdxFiles;
 
@@ -127,11 +126,6 @@ public class CountryAgent extends AbstractCountryAgent {
 		if (saveGamsGdxFiles && ModelConfig.PRICE_ELASTIC_DEMAND)
 			saveGDXFile("demand");
 		
-		if (currentTimestep.isInitialTimestep()) {
-			double totalAnimalProductDemand = currentProjectedDemand.get(CommodityType.MONOGASTRICS) + currentProjectedDemand.get(CommodityType.RUMINANTS);
-			animalFeedFromOtherSources = totalAnimalProductDemand * ModelConfig.ANIMAL_FEED_FROM_OTHER_SOURCES_RATE;
-		}
-						
 		if (currentProjectedDemand.size() == 0) {
 			LogWriter.printlnError("No demand for country " + country + " so skipping it");
 		}
@@ -221,7 +215,7 @@ public class CountryAgent extends AbstractCountryAgent {
 		}
 
 		GamsCountryInput countryLevelInputs = new GamsCountryInput(country, currentProjectedDemand, currentCountryPrices, importConstraints, 
-				previousGamsRasterOutput.getCropUsageData(), currentMinDemandFract, animalFeedFromOtherSources, subsidyRates);	
+				previousGamsRasterOutput.getCropUsageData(), currentMinDemandFract, subsidyRates);	
 		GamsRasterInput input = new GamsRasterInput(currentTimestep, countryYieldSurfaces, previousGamsRasterOutput.getLandUses(), irrigData, countryLevelInputs);
 
 		return input;
diff --git a/src/ac/ed/lurg/country/GlobalPrice.java b/src/ac/ed/lurg/country/GlobalPrice.java
index bf0ac8987c979c2ab8fba1243d5c6bcf86c1a1fe..eb54a680298e8cf139aee46f160bbed03f4c4d33 100644
--- a/src/ac/ed/lurg/country/GlobalPrice.java
+++ b/src/ac/ed/lurg/country/GlobalPrice.java
@@ -73,7 +73,12 @@ public class GlobalPrice implements Serializable {
 			double transportLossRate = (!ModelConfig.SHOCKS_POSSIBLE) ? ModelConfig.TRANSPORT_LOSSES : ModelConfig.getParameter(Parameter.TRANSPORT_LOSSES, timestep.getYear());
 			double newTransportLosses = newExportAmountBeforeLoss * transportLossRate;
 			double stockChange = newExportAmountBeforeLoss - newTransportLosses - newImports - oldDiff;
-			double updatedStock = stockLevel + stockChange;
+			
+			double updatedStock;
+			if (ModelConfig.IS_CALIBRATION_RUN  && timestep.getTimestep() <= ModelConfig.END_FIRST_STAGE_CALIBRATION)
+				updatedStock = stockLevel;  // don't update stock in inital stage of calibration
+			else
+				updatedStock = stockLevel + stockChange;
 				
 			LogWriter.println(String.format("     imports %.2f, exports %.2f", newImports, newExportAmountBeforeLoss - newTransportLosses));
 			LogWriter.println(String.format("     updatedStock %.2f, previous %.2f (last timestep %.2f), stockChange %.2f, oldDiff %.2f", updatedStock, stockLevel, stockLevel-oldDiff, stockChange, oldDiff));
diff --git a/src/ac/ed/lurg/country/gams/GamsCountryInput.java b/src/ac/ed/lurg/country/gams/GamsCountryInput.java
index b89b7732e68f7ecfa618312a02584f781088d09e..56e8ef9eed73cd23b36da068a236cf196963d308 100644
--- a/src/ac/ed/lurg/country/gams/GamsCountryInput.java
+++ b/src/ac/ed/lurg/country/gams/GamsCountryInput.java
@@ -19,12 +19,11 @@ public class GamsCountryInput {
 	private Map<CropType, CountryPrice> countryPrices;
 	private Map<CropType, CropUsageData> previousCropUsageData;
 	private Map<CommodityType, Map<CropType, Double>> minDemandFractions;
-	private double animalFeedFromOtherSources;
 	private Map<CropType, Double> subsidyRates;
 
 	public GamsCountryInput(CompositeCountry country, Map<CommodityType, Double> projectedDemand, Map<CropType, CountryPrice> countryPrices, 
 			Map<CropType, TradeConstraint> importConstraints, Map<CropType, CropUsageData> previousCropUsageData, 
-			Map<CommodityType, Map<CropType, Double>> minDemandFracts, double animalFeedFromOtherSources, Map<CropType, Double> subsidyRates) {
+			Map<CommodityType, Map<CropType, Double>> minDemandFracts, Map<CropType, Double> subsidyRates) {
 		super();
 		this.country = country;
 		this.projectedDemand = projectedDemand;
@@ -32,7 +31,6 @@ public class GamsCountryInput {
 		this.countryPrices = countryPrices;
 		this.previousCropUsageData = previousCropUsageData;
 		this.minDemandFractions = minDemandFracts;
-		this.animalFeedFromOtherSources = animalFeedFromOtherSources;
 		this.subsidyRates = subsidyRates;
 	}
 		
@@ -82,10 +80,6 @@ public class GamsCountryInput {
 		return minDemandFractions;
 	}
 	
-	public double getAnimalFeedFromOtherSources() {
-		return animalFeedFromOtherSources;
-	}
-	
 	public Map<CropType, Double> getSubsidyRates() {
 		return subsidyRates;
 	}
diff --git a/src/ac/ed/lurg/country/gams/GamsLocationOptimiser.java b/src/ac/ed/lurg/country/gams/GamsLocationOptimiser.java
index 016fa1658df1dfeeda447e4d01a1a2578fa44092..af4e20c6a75a76683bf519b623c553fb65697916 100644
--- a/src/ac/ed/lurg/country/gams/GamsLocationOptimiser.java
+++ b/src/ac/ed/lurg/country/gams/GamsLocationOptimiser.java
@@ -31,6 +31,7 @@ import ac.ed.lurg.types.CommodityType;
 import ac.ed.lurg.types.CropType;
 import ac.ed.lurg.types.LandCoverType;
 import ac.ed.lurg.types.Parameter;
+import ac.ed.lurg.types.YieldType;
 import ac.ed.lurg.utils.LazyHashMap;
 import ac.ed.lurg.utils.LogWriter;
 import ac.ed.lurg.yield.YieldResponsesItem;
@@ -227,13 +228,17 @@ public class GamsLocationOptimiser {
 				double maxIrrig = irrigationItem.getMaxIrrigAmount(crop);
 				
 				if (DEBUG) LogWriter.println(String.format("%d      %15s,\t %.1f,\t %.1f, \t %.1f,\t %.1f,\t %.2f,\t\t [%.2f],\t [%.2f],\t {%.2f}", 
-						locationId, crop.getGamsName(), yresp.getYieldNone(crop), yresp.getYieldFertOnly(crop), yresp.getYieldIrrigOnly(crop), yresp.getYieldMax(crop), yresp.getShockRate(crop), 
-						yresp.getFertParam(crop), yresp.getIrrigParam(crop), maxIrrig));
-
-				setGamsParamValue(yNoneP.addRecord(v), yresp.getYieldNone(crop), 4);
-				setGamsParamValue(y_fert.addRecord(v), yresp.getYieldFertOnly(crop), 4);
-				setGamsParamValue(y_irrig.addRecord(v), yresp.getYieldIrrigOnly(crop), 4);
-				setGamsParamValue(y_both.addRecord(v), yresp.getYieldMax(crop), 4);
+						locationId, crop.getGamsName(), 
+						yresp.getExtrapolatedYield(YieldType.NO_FERT_NO_IRRIG, crop), 
+						yresp.getExtrapolatedYield(YieldType.FERT_MAX_NO_IRRIG, crop), 
+						yresp.getExtrapolatedYield(YieldType.NO_FERT_IRRIG_MAX, crop), 
+						yresp.getExtrapolatedYield(YieldType.FERT_MAX_IRRIG_MAX, crop), 
+						yresp.getShockRate(crop), yresp.getFertParam(crop), yresp.getIrrigParam(crop), maxIrrig));
+
+				setGamsParamValue(yNoneP.addRecord(v), yresp.getExtrapolatedYield(YieldType.NO_FERT_NO_IRRIG, crop), 4);
+				setGamsParamValue(y_fert.addRecord(v), yresp.getExtrapolatedYield(YieldType.FERT_MAX_NO_IRRIG, crop), 4);
+				setGamsParamValue(y_irrig.addRecord(v), yresp.getExtrapolatedYield(YieldType.NO_FERT_IRRIG_MAX, crop), 4);
+				setGamsParamValue(y_both.addRecord(v), yresp.getExtrapolatedYield(YieldType.FERT_MAX_IRRIG_MAX, crop), 4);
 				setGamsParamValue(y_shock.addRecord(v), yresp.getShockRate(crop), 4);
 				setGamsParamValue(fert_p.addRecord(v), yresp.getFertParam(crop), 4);
 				setGamsParamValue(irrig_p.addRecord(v), yresp.getIrrigParam(crop), 4);
@@ -310,7 +315,6 @@ public class GamsLocationOptimiser {
 		addScalar(inDB, "unhandledCropRate", ModelConfig.UNHANDLED_CROP_RATE, 3);
 		addScalar(inDB, "setAsideRate", ModelConfig.SETASIDE_RATE, 5);
 		addScalar(inDB, "domesticPriceMarkup", ModelConfig.DOMESTIC_PRICE_MARKUP, 3);
-		addScalar(inDB, "animalFeedFromOtherSources", countryInput.getAnimalFeedFromOtherSources(), 2);		
 		double maxExpansion = 1.0;
 		if (!ModelConfig.IS_CALIBRATION_RUN && countryInput.getCountry().getName().equals("China")) {
 			maxExpansion = ModelConfig.MAX_CHINA_LAND_EXPANSION_RATE;
diff --git a/src/ac/ed/lurg/country/gams/GamsRasterOptimiser.java b/src/ac/ed/lurg/country/gams/GamsRasterOptimiser.java
index b6602c587c6356b96f3672708019d982c47bd8a2..75ff7437b552ea2461cb28b7485c0842684bc48e 100644
--- a/src/ac/ed/lurg/country/gams/GamsRasterOptimiser.java
+++ b/src/ac/ed/lurg/country/gams/GamsRasterOptimiser.java
@@ -290,8 +290,8 @@ public class GamsRasterOptimiser {
 			RasterKey key = entry.getKey();
 			CropType crop = CropType.WHEAT;
 
-			if (yresp != null && (yresp.getYieldMax(crop) < yresp.getYieldFertOnly(crop) || yresp.getYieldMax(crop) < yresp.getYieldIrrigOnly(crop))) {
-				logWarningWithCoord("Inconsistency F only:" + yresp.getYieldFertOnly(crop) + ", I only" + yresp.getYieldIrrigOnly(crop) + ", max " + yresp.getYieldMax(crop) + " at ", key, yieldRaster);
+			if (yresp != null && (yresp.getYield(YieldType.FERT_MAX_IRRIG_MAX, crop) < yresp.getYield(YieldType.FERT_MAX_NO_IRRIG, crop) || yresp.getYield(YieldType.FERT_MAX_IRRIG_MAX, crop) < yresp.getYield(YieldType.NO_FERT_IRRIG_MAX, crop))) {
+				logWarningWithCoord("Inconsistency F only:" + yresp.getYield(YieldType.FERT_MAX_NO_IRRIG, crop) + ", I only" + yresp.getYield(YieldType.NO_FERT_IRRIG_MAX, crop) + ", max " + yresp.getYield(YieldType.FERT_MAX_IRRIG_MAX, crop) + " at ", key, yieldRaster);
 			}
 		}
 		
diff --git a/src/ac/ed/lurg/utils/cluster/Cluster.java b/src/ac/ed/lurg/utils/cluster/Cluster.java
index d96f36409648942bd9b48409334db8e2a6cf67f6..1ad395fb9d24e00d7c0d599fa1d67a1001307128 100644
--- a/src/ac/ed/lurg/utils/cluster/Cluster.java
+++ b/src/ac/ed/lurg/utils/cluster/Cluster.java
@@ -55,7 +55,9 @@ public class Cluster<K, P extends ClusteringPoint<K>> {
     protected double distanceFromCentroid(ClusteringPoint<K> p) {
     	double squaredTotal=0;
     	for (K key : centroid.getAllClusteringKeys()) {
-       		squaredTotal += Math.pow(centroid.getClusteringValue(key)-p.getClusteringValue(key), 2);
+    		double a = centroid.getClusteringValue(key);
+    		double b = p.getClusteringValue(key);
+    		squaredTotal += Math.pow(a-b, 2);
     	}
     	
         return Math.sqrt(squaredTotal);
diff --git a/src/ac/ed/lurg/yield/YieldClusterPoint.java b/src/ac/ed/lurg/yield/YieldClusterPoint.java
index f8f5cfc3f07cfa0a30945401b82c0a2b3545362d..7b3882eb2f70f7e23173bd9d18872f6ab639eeff 100644
--- a/src/ac/ed/lurg/yield/YieldClusterPoint.java
+++ b/src/ac/ed/lurg/yield/YieldClusterPoint.java
@@ -5,6 +5,7 @@ import java.util.Collection;
 
 import ac.ed.lurg.landuse.IrrigationItem;
 import ac.ed.lurg.types.CropType;
+import ac.ed.lurg.types.YieldType;
 import ac.ed.lurg.utils.LogWriter;
 import ac.ed.lurg.utils.cluster.ClusteringPoint;
 import ac.sac.raster.RasterKey;
@@ -34,12 +35,12 @@ public class YieldClusterPoint implements ClusteringPoint<String> {
 		this.rasterKey = rasterKey;
 		
 		// not sure if we be better to get a reference to the YieldResponsesItem, rather than caching these values?
-		this.wheatMin = yields.getYieldNone(CropType.WHEAT);
-		this.riceMax = yields.getYieldMax(CropType.RICE);
-		this.maizeMin = yields.getYieldNone(CropType.MAIZE);
-		this.pasture = yields.getYieldNone(CropType.PASTURE);
-		this.wheatMax = yields.getYieldMax(CropType.WHEAT);
-		this.maizeMax = yields.getYieldMax(CropType.MAIZE);
+		this.wheatMin = yields.getYield(YieldType.NO_FERT_NO_IRRIG, CropType.WHEAT);
+		this.riceMax = yields.getYield(YieldType.FERT_MAX_IRRIG_MAX, CropType.RICE);
+		this.maizeMin = yields.getYield(YieldType.NO_FERT_NO_IRRIG, CropType.MAIZE);
+		this.pasture = yields.getYield(YieldType.NO_FERT_NO_IRRIG, CropType.PASTURE);
+		this.wheatMax = yields.getYield(YieldType.FERT_MAX_IRRIG_MAX, CropType.WHEAT);
+		this.maizeMax = yields.getYield(YieldType.FERT_MAX_IRRIG_MAX, CropType.MAIZE);
 		this.irrigWaterAval = irrigItem.getIrrigConstraint();
 		this.protectedArea = protectedArea;
 	}
@@ -50,7 +51,11 @@ public class YieldClusterPoint implements ClusteringPoint<String> {
 
 	@Override
 	public double getClusteringValue(String key) {
-
+		double v = getClusteringValueNaN(key);		
+		return Double.isNaN(v) ? 0.0 : v;  // LPJ with emulator provides NaN yields.  How to cluster these is tricky, but assuming they are 0 
+	}
+	
+	private double getClusteringValueNaN(String key) {
 		switch (key) {
 			case PASTURE:  return pasture;
 			case WHEAT_MIN:  return wheatMin;
@@ -60,9 +65,9 @@ public class YieldClusterPoint implements ClusteringPoint<String> {
 			case MAIZE_MAX:  return maizeMax;
 			case IRRIG_WATER_AVAL:  return irrigWaterAval;
 			case PROTECTED_AREA:  return protectedArea;
+			default:
+				throw new RuntimeException("YieldClusterPoint.getClusteringValue: got unknown value " + key);
 		}
-		LogWriter.printlnError("YieldClusterPoint.getClusteringValue: got unknown value " + key);
-		return Double.NaN;
 	}
 
 	@Override
diff --git a/src/ac/ed/lurg/yield/YieldResponse.java b/src/ac/ed/lurg/yield/YieldResponse.java
index 86e1e9b2a66110c20f554b4cfd9a50422a2d7230..aefcc740829559ce3c95a6dd44321a49041f7e73 100644
--- a/src/ac/ed/lurg/yield/YieldResponse.java
+++ b/src/ac/ed/lurg/yield/YieldResponse.java
@@ -5,9 +5,11 @@ import java.util.Map;
 
 import ac.ed.lurg.ModelConfig;
 import ac.ed.lurg.types.YieldType;
+import ac.ed.lurg.utils.LogWriter;
 
 public class YieldResponse {
 	private Map<YieldType, Double> yields = new HashMap<YieldType, Double>();
+	private Map<YieldType, Double> yieldsExtrapolated;
 	private double yieldShock;
 	private double fertParm;
 	private double irrigParm;
@@ -15,7 +17,7 @@ public class YieldResponse {
 	private static double defaultParam;
 	
 	static {
-		defaultParam =  calcParam(0, 0.8, 1, 0, 0.5, 1.0);  // default assumes 80% at mid point
+		defaultParam =  calcParam(0, 0.8, 1, 0.5);  // default assumes 80% at mid point
 	}
 
 	
@@ -36,14 +38,59 @@ public class YieldResponse {
 		return d.doubleValue();
 	}		
 	
+	public double getExtrapolatedYield(YieldType type) {
+		if (yieldsExtrapolated == null)
+			calcParams();
+		
+		Double d = yieldsExtrapolated.get(type);
+		
+		if (d == null)
+			return Double.NaN;
+		
+		return d.doubleValue();
+	}
+	
 	public double getFertParam() {
-		if (fertParm == 0) {
-			//	fertParm = calcParam (0, 0.5, 1, 5.0/200, 50.0/200, 200.0/200); // we do have MID fert data, but looks wrong 
-			fertParm = calcParam(getYield(YieldType.NO_FERT_NO_IRRIG), getYield(YieldType.FERT_MID_NO_IRRIG), getYield(YieldType.FERT_MAX_NO_IRRIG), 
-					ModelConfig.MIN_FERT_AMOUNT/ModelConfig.MAX_FERT_AMOUNT, ModelConfig.MID_FERT_AMOUNT/ModelConfig.MAX_FERT_AMOUNT, ModelConfig.MAX_FERT_AMOUNT/ModelConfig.MAX_FERT_AMOUNT);
+		if (fertParm == 0)
+			calcParams();
+
+		return fertParm;
+	}
+	
+	private double calcParams() {
+		double yMinNoI = getYield(YieldType.NO_FERT_NO_IRRIG);
+		double yMidNoI = getYield(YieldType.FERT_MID_NO_IRRIG);
+		double yMaxNoI = getYield(YieldType.FERT_MAX_NO_IRRIG);
+		double fMid = (ModelConfig.MID_FERT_AMOUNT-ModelConfig.MIN_FERT_AMOUNT)/(ModelConfig.MAX_FERT_AMOUNT-ModelConfig.MIN_FERT_AMOUNT);
+		
+		if (ModelConfig.EXTRAPOLATE_YIELD_FERT_RESPONSE) {
+			double yMinI = getYield(YieldType.NO_FERT_IRRIG_MAX);
+			double yMidI = getYield(YieldType.FERT_MID_IRRIG_MAX);
+			double yMaxI = getYield(YieldType.FERT_MAX_IRRIG_MAX);
+
+			double asymptoteYieldIncNoI = yieldAsymptoteSolve(yMidNoI-yMinNoI, yMaxNoI-yMinNoI, fMid);
+			fertParm = -1/fMid * Math.log(1-(yMidI-yMinI)/asymptoteYieldIncNoI);
+			double asymptoteYieldIncI = (yMaxI-yMinI) / (1 - Math.exp(-fertParm));
+			//LogWriter.println(String.format("fertParm from solver: asymptoteYieldIncNoI: %.2f, asymptoteYieldIncI: %.2f, fertParm %.4f", asymptoteYieldIncNoI, asymptoteYieldIncI, fertParm));
 			
-			//LogWriter.println(String.format("%s %s %s", getYield(YieldType.NO_FERT_NO_IRRIG), getYield(YieldType.FERT_MID_NO_IRRIG), getYield(YieldType.FERT_MAX_NO_IRRIG)));
+			if (Double.isNaN(asymptoteYieldIncI) || Double.isNaN(fertParm)) {
+				LogWriter.println("Not finding a suitable extrapolating solution.  Defaulting to old style");
+			}
+			else {
+				yieldsExtrapolated = new HashMap<YieldType, Double>();
+				yieldsExtrapolated.put(YieldType.NO_FERT_NO_IRRIG, yMinNoI);
+				yieldsExtrapolated.put(YieldType.FERT_MAX_NO_IRRIG, asymptoteYieldIncNoI+yMinI);
+				yieldsExtrapolated.put(YieldType.NO_FERT_IRRIG_MAX, yMinI);
+				yieldsExtrapolated.put(YieldType.FERT_MAX_IRRIG_MAX, asymptoteYieldIncI+yMinI);
+			}
 		}
+		
+		// old style, if not wanting extrapolate or data doesn't allow
+		if (yieldsExtrapolated == null) {
+			fertParm = calcParam(yMinNoI, yMidNoI, yMaxNoI, fMid);
+			yieldsExtrapolated = yields;
+		}
+		
 		return fertParm;
 	}
 	
@@ -55,11 +102,11 @@ public class YieldResponse {
 		return irrigParm;
 	}
 	
-	private static double calcParam(double yMin, double yMid, double yMax, double xMin, double xMid, double xMax) {
+	private static double calcParam(double yMin, double yMid, double yMax, double xMid) {
 		if (yMid <= yMin || yMax <= yMid)
 			return defaultParam;
-
-		return -Math.log(1 - ((yMid-yMin)/(yMax-yMin))) * (xMax-xMin) / (xMid-xMin); 					// 1-exp(-x) function
+		
+		return -Math.log(1 - ((yMid-yMin)/(yMax-yMin))) / xMid; 					// 1-exp(-x) function
 	//	return Math.log((yMid - yMin)/(yMax - yMin)) / Math.log((xMid - xMin)/(xMax - xMin)) / xMax;  	// Cobb-Douglas
 	}
 	
@@ -70,4 +117,48 @@ public class YieldResponse {
 	public double getShock() {
 		return yieldShock;
 	}
+	
+	private static double f(double x, double yM, double yH, double fM) {
+	   	return 1-yM/x - Math.pow(1-yH/x, fM);
+	}
+
+	private double yieldAsymptoteSolve (double yM, double yH, double fM) {
+
+		double a = yH;
+		double b = yH * 20;
+		double epsilon = 0.0001;
+		
+		if (yM < 0 || yH < 0 || fM * yH > yM) {
+			LogWriter.printlnError(String.format("yieldAsymptoteSolve: less than zero yield or not diminishing yield increases to N: yM=%.3f, yH=%.3f, fM=%.3f", yM, yH, fM));
+			return Double.NaN;
+		}
+		
+		double m = Double.NaN;
+		double fm;
+		double fa = f(a, yM, yH, fM);
+		double fb = f(b, yM, yH, fM);
+
+		if (fa * fb > 0) {
+			LogWriter.printlnError(String.format("yieldAsymptoteSolve: not bracketing root: yM=%.3f, yH=%.3f, fM=%.3f", yM, yH, fM));
+			return Double.NaN;
+		}
+		
+		while ((b-a) > epsilon)
+		{
+			m = (a+b)/2;           // Mid point
+			fm = f(m, yM, yH, fM);
+
+			if ( fm * fa > 0 )
+			{  // f(a) and f(m) have same signs: move a
+				a = m;
+				fa = fm;
+			}
+			else
+			{  // f(a) and f(m) have different signs: move b
+				b = m;
+				fb = fm;
+			}
+		}
+		return m;
+	}
 }
\ No newline at end of file
diff --git a/src/ac/ed/lurg/yield/YieldResponsesItem.java b/src/ac/ed/lurg/yield/YieldResponsesItem.java
index e78fa2ff79316084cd243949d99aa0d0a9265720..e6723643f8a33dee46399e90f8eb44eaeef77d1b 100644
--- a/src/ac/ed/lurg/yield/YieldResponsesItem.java
+++ b/src/ac/ed/lurg/yield/YieldResponsesItem.java
@@ -29,26 +29,11 @@ public class YieldResponsesItem implements RasterItem {
 		return d;
 	}	
 
-	public double getYieldNone(CropType crop) {
-		return getYield(YieldType.NO_FERT_NO_IRRIG, crop);
+	public double getExtrapolatedYield(YieldType yieldType, CropType crop) {
+		YieldResponse yr = getYieldResponseForCrop(crop);
+		return yr.getExtrapolatedYield(yieldType);
 	}	
 	
-	public double getYieldMidFertOnly(CropType crop) {
-		return getYield(YieldType.FERT_MID_NO_IRRIG, crop);
-	}
-
-	public double getYieldFertOnly(CropType crop) {
-		return getYield(YieldType.FERT_MAX_NO_IRRIG, crop);
-	}
-	
-	public double getYieldIrrigOnly(CropType crop) {
-		return getYield(YieldType.NO_FERT_IRRIG_MAX, crop);
-	}
-	
-	public double getYieldMax(CropType crop) {
-		return getYield(YieldType.FERT_MAX_IRRIG_MAX, crop);
-	}
-	
 	public double getFertParam(CropType crop) {
 		return getYieldResponseForCrop(crop).getFertParam();
 	}