Merge branch 'patrick' into Camille_Bis

.gitignore

@@ -102,3 +102,4 @@ target
 org.*
 /src/main/java/miscellaneous/Essai.java
 >>>>>>> master
+/exploration/R/.Rhistory

data/input/fishTryRealBV_CC.xml

@@ -229,18 +229,13 @@
 		<outputPath>data/output/</outputPath>
 		<fileNameFluxes>nutrientFluxes</fileNameFluxes>
 		<lengthAtHatching>2.0</lengthAtHatching>
-		<linfVonBertForFemale>70.0</linfVonBertForFemale>
-		<linfVonBertForMale>70.0</linfVonBertForMale>
+		<linfVonBertForFemale>80.0</linfVonBertForFemale>
+		<linfVonBertForMale>80.0</linfVonBertForMale>
 		<lFirstMaturityForFemale>55.0</lFirstMaturityForFemale>
 		<lFirstMaturityForMale>40.0</lFirstMaturityForMale>
 		<processes>
 			<processesAtBegin>
-				<species.PopulateBasinNetwork>
-					<synchronisationMode>ASYNCHRONOUS</synchronisationMode>
-					<nbSIPerBasin>200</nbSIPerBasin>
-					<initialLength>20.0</initialLength>
-					<nbFishPerSI>2500</nbFishPerSI>
-				</species.PopulateBasinNetwork>
+		
 			</processesAtBegin>
 
 			<processesEachStep>
@@ -252,6 +247,15 @@
 					<temporisation>0</temporisation>
 				</species.PlopProcess>
 
+				<species.PopulateBasinNetworkSeveralTimes>
+					<synchronisationMode>ASYNCHRONOUS</synchronisationMode>
+					<nbSIPerBasin>200</nbSIPerBasin>
+					<nbFishPerSI>2500</nbFishPerSI>
+					<initialLength>20.0</initialLength>
+					<timesOfPopulate>5</timesOfPopulate>
+					<populateSeason>SPRING</populateSeason>
+				</species.PopulateBasinNetworkSeveralTimes>
+				
 				<species.Age>
 					<synchronisationMode>ASYNCHRONOUS</synchronisationMode>
 				</species.Age>
@@ -261,8 +265,8 @@
 					<tempMinGrow>3.0</tempMinGrow>
 					<tempMaxGrow>26.0</tempMaxGrow>
 					<tempOptGrow>17.0</tempOptGrow>
-					<kOptForFemale>0.20323011454056628</kOptForFemale>
-					<kOptForMale>0.3004992503785338</kOptForMale>
+					<kOptForFemale>0.29</kOptForFemale>
+					<kOptForMale>0.21</kOptForMale>
 					<sigmaDeltaLVonBert>0.2</sigmaDeltaLVonBert>
 				</species.Grow>
 
@@ -287,7 +291,7 @@
 
 				<species.Survive>
 					<synchronisationMode>ASYNCHRONOUS</synchronisationMode>
-					<tempMinSurvivalSpawnerInRiv>10.0</tempMinSurvivalSpawnerInRiv>
+					<tempMinSurvivalSpawnerInRiv>9.0</tempMinSurvivalSpawnerInRiv>
 					<tempMaxSurvivalSpawnerInRiv>23.0</tempMaxSurvivalSpawnerInRiv>
 					<!-- <tempOptSurvivalSpawnerInRiv>20.0</tempOptSurvivalSpawnerInRiv> -->
 					<survivalProbOptSpawnerInRiv>1.0</survivalProbOptSpawnerInRiv>
@@ -312,7 +316,7 @@
 				<species.ReproduceAndSurviveAfterReproductionWithDiagnose>
 					<synchronisationMode>ASYNCHRONOUS</synchronisationMode>
 					<reproductionSeason>SPRING</reproductionSeason>
-					<tempMinRep>10.443039492695249</tempMinRep>
+					<tempMinRep>9.8</tempMinRep>
 					<tempMaxRep>26.0</tempMaxRep>
 					<tempOptRep>20.0</tempOptRep>
 					<eta>2.4</eta>

data/input/obsTryReal.xml

@@ -31,6 +31,30 @@
 					<variableName>getFishEffective</variableName>
 				</fr.cemagref.observation.observers.jfreechart.TemporalSerieChart>
 
+				<fr.cemagref.observation.observers.jfreechart.TemporalSerieChart>
+					<graphType>LINE</graphType>
+					<title>Abundance of spawners</title>
+					<xAxisLabel>Time (season)</xAxisLabel>
+					<yAxisLabel>number of fishes</yAxisLabel>
+					<variableName>getSpawnerEffective</variableName>
+				</fr.cemagref.observation.observers.jfreechart.TemporalSerieChart>
+
+				<fr.cemagref.observation.observers.jfreechart.TemporalSerieChart>
+					<graphType>LINE</graphType>
+					<title>Abundance of male spawners</title>
+					<xAxisLabel>Time (season)</xAxisLabel>
+					<yAxisLabel>number of fishes</yAxisLabel>
+					<variableName>getMaleSpawnerEffective</variableName>
+				</fr.cemagref.observation.observers.jfreechart.TemporalSerieChart>
+
+				<fr.cemagref.observation.observers.jfreechart.TemporalSerieChart>
+					<graphType>LINE</graphType>
+					<title>Abundance of female spawners</title>
+					<xAxisLabel>Time (season)</xAxisLabel>
+					<yAxisLabel>number of fishes</yAxisLabel>
+					<variableName>getFemaleSpawnerEffective</variableName>
+				</fr.cemagref.observation.observers.jfreechart.TemporalSerieChart>
+				
 				<miscellaneous.TemporalRangeSerieChart>
 					<title>Range distribution</title>
 					<xAxisLabel>Time (season)</xAxisLabel>
@@ -59,7 +83,7 @@
 					<graphType>LINE</graphType>
 					<title>Statistic for male spawners age</title>
 					<xAxisLabel>Time (season)</xAxisLabel>
-					<yAxisLabel>age (year)</yAxisLabel>
+					<yAxisLabel></yAxisLabel>
 					<variableName>computeMaleSpawnerForFirstTimeSummaryStatistic</variableName>
 				</fr.cemagref.observation.observers.jfreechart.TemporalSerieChart>
 
@@ -67,7 +91,7 @@
 					<graphType>LINE</graphType>
 					<title>Statistic for female spawners age</title>
 					<xAxisLabel>Time (season)</xAxisLabel>
-					<yAxisLabel>age (year)</yAxisLabel>
+					<yAxisLabel></yAxisLabel>
 					<variableName>computeFemaleSpawnerForFirstTimeSummaryStatistic</variableName>
 				</fr.cemagref.observation.observers.jfreechart.TemporalSerieChart>
 

exploration/scriptR/BruchPrg.R

+library(openxlsx)
+# ============================================================
+
+dataBruch = read.xlsx("BDalosesBruch.xlsx")
+dataBruch$`M.gonades.(g)` = as.numeric(dataBruch$`M.gonades.(g)`)
+
+head(dataBruch)
+
+
+tapply(dataBruch$`Lf.(cm)`, dataBruch[,c('Année', 'Sexe')],min, na.rm = TRUE)
+
+tapply(dataBruch$`Lf.(cm)`, dataBruch[,c('Année', 'Sexe')],min, na.rm = TRUE)
+tapply(dataBruch$`Lf.(cm)`, dataBruch[,c('Année', 'Sexe')],quantile, na.rm = TRUE, probs=.05)
+
+tapply(dataBruch$`Lf.(cm)`, dataBruch[,c('Année', 'Sexe')],max, na.rm = TRUE)
+tapply(dataBruch$`Lf.(cm)`, dataBruch[,c('Année', 'Sexe')],quantile, na.rm = TRUE, probs=.95)
+
+tapply(dataBruch$`Lf.(cm)`, dataBruch[,c('Année', 'Sexe')],quantile, na.rm = TRUE, probs=.5)
+
+sel = dataBruch$Année==2013 & dataBruch$Sexe =='F'
+
+hist(dataBruch$`Lt.(cm)`[sel])
+
+abline(v=quantile(dataBruch$`Lt.(cm)`[sel], probs = 0.05))
+
+sel = dataBruch$Sexe=='M'
+lm (dataBruch$`Lt.(cm)`[sel]~dataBruch$`Lf.(cm)`[sel])
+
+summary(lm (dataBruch$`Lt.(cm)`~dataBruch$`Lf.(cm)`))
+summary(lm (dataBruch$`Lt.(cm)`~dataBruch$`Lf.(cm)` * dataBruch$Sexe))
+
+# ====================================================
+# fecundity
+# Taverny 1991
+# ====================================================
+(41*172895+33*202902+74*186424)/(41+33+74)
+(41*98390+33*110386+74*104325)/(41+33+74)
+
+# ============================================================
+# maximal production of recruit in GR3D for the Garonne basin
+# ============================================================
+bj=-log(1.7e-3) /.33
+cj= 4.1e-4 / (84810*.5356)
+(alphaj = bj*exp(-bj*.33)/(cj*(1-exp(-bj*.33))))
+
+# ================================================
+# growth in GR3D
+# ================================================
+
+# ------------------------------------------------
+# temperature effect on growth
+# ------------------------------------------------
+temperatureEffect= function(temp, Tmin, Topt, Tmax){
+#  if (temp<=Tmin | temp >= Tmax)
+#    return(0)
+#  else
+    response=(temp-Tmin)*(temp-Tmax)/((temp-Tmin)*(temp-Tmax)-(temp-Topt)^2)
+    
+    response[temp<=Tmin | temp >= Tmax] = 0
+      return(response)
+
+  }
+temperature=seq(8,30,.1)
+# temperature effect on spawner survival 
+plot(temperature, temperatureEffect(temperature, 10, 20, 23), type='l')
+# temperature effect on recruit survival 
+lines(temperature, temperatureEffect(temperature, 9.75, 20, 26), , col='red')
+
+
+lines(temperature, temperatureEffect(temperature, 9.75, 20, 26) * temperatureEffect(temperature, 10, 20, 23), type='l', col='green')
+lines(temperature, temperatureEffect(temperature, 9.75, 20, 26) * temperatureEffect(temperature, 10, 20, 23) * exp(-.4*5), type='l', col='blue')
+tempData=read.csv("/home/patrick.lambert/Documents/workspace/GR3D/data/input/reality/SeasonTempBVFacAtlant1801_2100_newCRU_RCP85.csv", sep=";")
+sel = tempData$NOM=="Garonne" & tempData$Year>=2008 & tempData$Year<=2018
+plot(tempData$Year[sel], tempData$Winter[sel], type='l')
+Tref=colMeans(tempData[sel, c("Winter", "Spring", "summer", "Autumn")])
+points(Tref,  temperatureEffect(Tref, 9.75, 20, 26), col="red")
+text(Tref, temperatureEffect(Tref, 9.75, 20, 26),  c("Winter", "Spring", "Summer", "Autumn"), pos=1)
+
+mean( temperatureEffect(Tref, 9.75, 20, 26))
+# ----------------------------------------------
+# growth simulation
+# ----------------------------------------------
+vonBertalaffyGrowth = function(age, L0, Linf, K){
+  t0=log(1-L0/Linf)/K
+  return(Linf*(1-exp(-K*(age-t0))))
+}
+
+Pauly= function(age, t0, Linf, K, D){
+  return(Linf/10*((1-exp(-K*D*(age-t0)))^(1/D)))
+}
+
+vonBertalaffyIncrement = function(nStep, L0, Linf, K, deltaT, sigma, withTempEffect=FALSE){
+  tempEffect = temperatureEffect( c(7.753891, 14.979708, 19.782974, 11.108207) , 3, 17, 26)
+  L=matrix(nrow=nStep+1)
+  L[1]=L0
+  for (i in 1:nStep) {
+    mu = log((Linf-L[i])*(1-exp(-K*deltaT))) - sigma*sigma/2
+    increment = exp(rnorm(1, mu, sigma))
+    if (withTempEffect){
+      increment = increment * tempEffect[((i-1) %% 4)+1]
+    }
+    L[i+1]=L[i]+increment
+  }
+  return(L)
+}
+
+vonBertalaffyIncrement(6/.25, 0, 60, 0.3900707, .25, .2)
+
+age=seq(0,6,.25)
+plot(age,vonBertalaffyGrowth(age, 2, 60, 0.3900707), type="l")
+for (i in 1:100) {
+  lines(age, vonBertalaffyIncrement(6/.25, 2, 60, 0.3900707, .25, .2), col='red')
+}
+lines(age, vonBertalaffyGrowth(age, 2, 60, 0.3900707), lwd=3, col='black')
+abline(h=40)
+for (i in 1:100) {
+  lines(age, vonBertalaffyIncrement(6/.25, 2, 60, 0.3900707, .25, .2, withTempEffect = TRUE), col='green')
+}
+ lines(age, vonBertalaffyGrowth(age, 2, 60, 0.3900707*mean(temperatureEffect(Tref, 3, 17, 26))), lty=2, lwd = 2)
+
+abline(h=40)
+
+
+
+nbRep=1000
+res=matrix(nrow=nbRep)
+for (i in 1:nbRep) {
+ prov = vonBertalaffyIncrement(24, 2, 60, 0.3900707, .25, .2)
+ res[i] = prov[max(which(prov < 40))+4]
+}
+mean(res)
+hist(res,20)
+abline(v=mean(res), col='red')
+
+res2=matrix(nrow=nbRep)
+for (i in 1:nbRep) {
+  res2[i] = vonBertalaffyIncrement(4, 40, 60, 0.3900707, .25, .2)[5]
+}
+mean(res2)
+hist(res2,20)
+abline(v=mean(res2), col='red')
+
+
+# ======================================================================
+# exploration of growth for male and female
+# ======================================================================
+correction=mean(temperatureEffect(Tref, 3, 17, 26))
+
+age=seq(0,10,.25)
+present = vonBertalaffyGrowth(age, 2, 60, 0.3900707 * correction)
+plot(age, present, type='l', lwd=3, ylim =c(0,80))
+present[age == 5]
+abline(v=5)
+
+male =  vonBertalaffyGrowth(age, 2, 65, 0.3900707 * correction)
+lines(age, male, type='l', lwd=3, ylim =c(0,80), col='blue')
+abline(h=40, col='blue', lwd=2, lty=2)
+male[age == 5]
+
+female =  vonBertalaffyGrowth(age, 2, 75, 0.3900707*55/40 * correction)
+lines(age, female, lwd=3, col='red')
+abline(h=55, col='red', lwd=2, lty=2)
+female[age == 5]
+
+## a partir d'individus en mer donc à croissance  de plus lente à mesure qu'ils sont agées
+(taverny = Pauly(age,t0=-0.7294, Linf=701.59, K=0.4491, D=.5912))
+lines (age, taverny, lwd=2, col ='green')
+taverny[age == 5]
+
+# ===================================================================
+# GR3D outputs
+# =====================================================================
+simData=read.csv("/home/patrick.lambert/Documents/workspace/GR3D/data/output/lengthAgeDistribution_1-RCP85.csv", sep=";", row.names = NULL)
+
+simGaronne= simData[simData$basin =="Garonne",]
+sel=simGaronne$nbSpawn == 0
+tapply(simGaronne$length[sel], simGaronne[sel,c('year')],quantile, na.rm = TRUE, probs=.5)
+
+
+# masse des gonades avant
+sel = (dataBruch$LOT =='Tuilières' | dataBruch$LOT =='Golfech') & !is.na(dataBruch$`M.gonades.(g)`) & dataBruch$Sexe =='F'
+mean(dataBruch$`M.gonades.(g)`[sel]/dataBruch$`M.tot.(g)`[sel])
+
+
+
+sel = (dataBruch$LOT =='Tuilières' | dataBruch$LOT =='Golfech') & dataBruch$Sexe =='F'
+sum(sel)
+Wpre = mean(dataBruch$`M.tot.(g)`[sel])
+Wgonad =mean(dataBruch$`M.gonades.(g)`[sel], na.rm = TRUE)
+
+sel = (! (dataBruch$LOT =='Tuilières' | dataBruch$LOT =='Golfech')) & dataBruch$Sexe =='F'
+Wpost= mean(dataBruch$`M.tot.(g)`[sel])
+WgonadSpent =mean(dataBruch$`M.gonades.(g)`[sel], na.rm = TRUE)
+(Wloss=(Wpre - Wpost)/Wpre)

src/main/java/environment/Time.java

@@ -20,9 +20,15 @@ public final class Time {
 				((BasinNetwork) pilot.getAquaticWorld().getEnvironment()).getYearOfTheBegin();		
 	}
 	
+	public static long getNbYearFromBegin(Pilot pilot){
+		return getNbYearFromBegin(pilot.getCurrentTime());
+
+	}
+	
 	public static long getNbYearFromBegin(long time){
 		return (long) Math.floor(time / Season.values().length);		
 	}
+	
 
 	/**
 	 * @return the duration of season (time step) 

src/main/java/species/Calibrate.java

@@ -164,13 +164,13 @@ class GR3DObjeciveFunction  implements   IObjectiveFunction {
 		try {
 			pilot.load();
 
-			ReflectUtils.setFieldValueFromPath(pilot, "aquaticWorld.aquaNismsGroupsList.0.processes.processesEachStep.6.tempMinRep", par[0]);
+			ReflectUtils.setFieldValueFromPath(pilot, "aquaticWorld.aquaNismsGroupsList.0.processes.processesEachStep.7.tempMinRep", par[0]);
 			//System.out.println("tempMinRep: " + (double)  ReflectUtils.getValueFromPath(pilot, "aquaticWorld.aquaNismsGroupsList.0.processes.processesEachStep.6.getTempMinRep"));
 
-			ReflectUtils.setFieldValueFromPath(pilot, "aquaticWorld.aquaNismsGroupsList.0.processes.processesEachStep.3.kOptForFemale", par[1]);
+			ReflectUtils.setFieldValueFromPath(pilot, "aquaticWorld.aquaNismsGroupsList.0.processes.processesEachStep.4.kOptForFemale", par[1]);
 			//System.out.println("KOptFemale: " + (double)  ReflectUtils.getValueFromPath(pilot, "aquaticWorld.aquaNismsGroupsList.0.processes.processesEachStep.3.getkOptForFemale"));
 
-			ReflectUtils.setFieldValueFromPath(pilot, "aquaticWorld.aquaNismsGroupsList.0.processes.processesEachStep.3.kOptForMale",par[2]);
+			ReflectUtils.setFieldValueFromPath(pilot, "aquaticWorld.aquaNismsGroupsList.0.processes.processesEachStep.4.kOptForMale",par[2]);
 			//System.out.println("KOptMale: " + (double)  ReflectUtils.getValueFromPath(pilot, "aquaticWorld.aquaNismsGroupsList.0.processes.processesEachStep.3.getkOptForMale"));
 
 		} catch (Exception   e1) {

src/main/java/species/DiadromousFishGroup.java

@@ -33,6 +33,7 @@ import java.util.regex.Pattern;
 import miscellaneous.Duo;
 import miscellaneous.TreeMapForCentile;
 import species.DiadromousFish.Gender;
+import species.DiadromousFish.Stage;
 
 import org.openide.util.lookup.ServiceProvider;
 
@@ -915,6 +916,56 @@ public class DiadromousFishGroup extends AquaNismsGroup< DiadromousFish, BasinNe
 		return eff;
 	}
 
+	/**
+	 * @return sum of spawner effectives in all the river basins
+	 */
+	@Observable(description = "Number of spawners in river basins")
+	public double getSpawnerEffective() {
+		long eff = 0;
+		for (RiverBasin basin : this.getEnvironment().getRiverBasins()){
+			if (basin.getFishs(this) != null) {
+				for (DiadromousFish fish : basin.getFishs(this)) {
+					if (fish.getStage() == Stage.MATURE)
+						eff += fish.getAmount();
+				}
+			}
+		}
+		return eff;
+	}
+	
+	/**
+	 * @return sum of male spawner effectives in all the river basins
+	 */
+	@Observable(description = "Number of male spawners in all basins")
+	public double getMaleSpawnerEffective() {
+		long eff = 0;
+		for (Basin basin : this.getEnvironment().getBasins()){
+			if (basin.getFishs(this) != null) {
+				for (DiadromousFish fish : basin.getFishs(this)) {
+					if (fish.getStage() == Stage.MATURE && fish.getGender() == Gender. MALE)
+						eff += fish.getAmount();
+				}
+			}
+		}
+		return eff;
+	}
+
+	/**
+	 * @return sum of female spawner effectives in all the river basins
+	 */
+	@Observable(description = "Number of female spawners in all basins")
+	public double getFemaleSpawnerEffective() {
+		long eff = 0;
+		for (Basin basin : this.getEnvironment().getBasins()){
+			if (basin.getFishs(this) != null) {
+				for (DiadromousFish fish : basin.getFishs(this)) {
+					if (fish.getStage() == Stage.MATURE && fish.getGender() == Gender. FEMALE)
+						eff += fish.getAmount();
+				}
+			}
+		}
+		return eff;
+	}
 
 	@Override
 	public void addAquaNism(DiadromousFish fish) {

src/main/java/species/PopulateBasinNetworkSeveralTimes.java

+package species;
+
+import com.thoughtworks.xstream.XStream;
+import com.thoughtworks.xstream.io.xml.DomDriver;
+
+import environment.Basin;
+import environment.Time;
+import environment.Time.Season;
+import fr.cemagref.simaqualife.kernel.processes.AquaNismsGroupProcess;
+import species.DiadromousFish.Gender;
+
+import org.openide.util.lookup.ServiceProvider;
+
+@ServiceProvider(service = AquaNismsGroupProcess.class)
+public class PopulateBasinNetworkSeveralTimes extends AquaNismsGroupProcess<DiadromousFish, DiadromousFishGroup> {
+
+	private int nbSIPerBasin=200;
+	private int nbFishPerSI=2500;
+	private double initialLength =20.;
+
+	private int timesOfPopulate = 10;
+	private Season populateSeason = Season.SPRING;
+
+	public static void main(String[] args) { System.out.println((new
+			XStream(new DomDriver())) .toXML(new PopulateBasinNetworkSeveralTimes())); }
+
+	@Override
+	public void doProcess(DiadromousFishGroup group) {
+
+		if (Time.getNbYearFromBegin(group.getPilot()) <= timesOfPopulate && Time.getSeason(group.getPilot()) == populateSeason){
+
+			int nbFemaleSIPerBasin = nbSIPerBasin / 2;
+			int nbMaleSIPerBasin = nbSIPerBasin - nbFemaleSIPerBasin;
+
+			for (Basin basin : group.getEnvironment().getRiverBasins()){
+				for (int i=0; i < nbFemaleSIPerBasin; i++){
+					group.addAquaNism(new DiadromousFish(group.getPilot(), basin, initialLength, nbFishPerSI, Gender.FEMALE));
+				}
+				for (int i=0; i < nbMaleSIPerBasin; i++){
+					group.addAquaNism(new DiadromousFish(group.getPilot(), basin, initialLength, nbFishPerSI, Gender.MALE));
+				}
+			}
+		}
+		
+		System.out.println("fem: "+ group. getFemaleSpawnerEffective()+ " \tmal"+group.getMaleSpawnerEffective());
+	}
+}

src/main/java/species/ReproduceAndSurviveAfterReproductionWithDiagnose.java

@@ -152,6 +152,7 @@ public class ReproduceAndSurviveAfterReproductionWithDiagnose extends AquaNismsG
 
 					// effective temperature for reproduction (priority to the ANG value) 
 					double tempEffectRep;
+					double currentTemp = riverBasin.getCurrentTemperature(group.getPilot());
 					if (Double.isNaN(group.getTempMinRep())){
 						tempEffectRep = Miscellaneous.temperatureEffect(riverBasin.getCurrentTemperature(group.getPilot()), tempMinRep, tempOptRep, tempMaxRep);
 					}

src/main/java/species/Survive.java

@@ -24,7 +24,7 @@ public class Survive extends AquaNismsGroupProcess<DiadromousFish, DiadromousFis
 	public double tempOptSurvivalSpawnerInRiv = 20.;
 	public double survivalProbOptSpawnerInRiv = 1.;
 
-	public double mortalityRateInRiver = 0.4;
+	//public double mortalityRateInRiver = 0.4;
 	public double mortalityRateInSea = 0.4;
 	public double mortalityRateInOffshore = 0.4;