Index: /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/c/classes/objects/Elements/Tria.cpp =================================================================== --- /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/c/classes/objects/Elements/Tria.cpp (revision 12865) +++ /u/astrid-r1b/morlighe/issmuci/trunk-jpl/../trunk-jpl/src/c/classes/objects/Elements/Tria.cpp (revision 12866) @@ -930,17 +930,17 @@ input->GetInputValue(&TemperaturesPresentday[iv][month],gauss,month/12.*yts); input2->GetInputValue(&TemperaturesLgm[iv][month],gauss,month/12.*yts); input3->GetInputValue(&PrecipitationsPresentday[iv][month],gauss,month/12.*yts); - monthlyprec[iv][month]=monthlyprec[iv][month]*yts; // convertion to m/yr + PrecipitationsPresentday[iv][month]=PrecipitationsPresentday[iv][month]/yts; // converion in m/sec } } /*Recover delta18o and Delta18oSurface at present day, lgm and at time t*/ - this->parameters->FindParam(&Delta18oPresent,SurfaceforcingsDelta18oEnum,finaltime*yts); - this->parameters->FindParam(&Delta18oLgm,SurfaceforcingsDelta18oEnum,(finaltime-21000)*yts); - this->parameters->FindParam(&Delta18oTime,SurfaceforcingsDelta18oEnum,time*yts); - this->parameters->FindParam(&Delta18oSurfacePresent,SurfaceforcingsDelta18oSurfaceEnum,finaltime*yts); - this->parameters->FindParam(&Delta18oSurfaceLgm,SurfaceforcingsDelta18oSurfaceEnum,(finaltime-21000)*yts); - this->parameters->FindParam(&Delta18oSurfaceTime,SurfaceforcingsDelta18oSurfaceEnum,time*yts); - + this->parameters->FindParam(&Delta18oPresent,SurfaceforcingsDelta18oEnum,finaltime); + this->parameters->FindParam(&Delta18oLgm,SurfaceforcingsDelta18oEnum,(finaltime-(21000*yts))); + this->parameters->FindParam(&Delta18oTime,SurfaceforcingsDelta18oEnum,time); + this->parameters->FindParam(&Delta18oSurfacePresent,SurfaceforcingsDelta18oSurfaceEnum,finaltime); + this->parameters->FindParam(&Delta18oSurfaceLgm,SurfaceforcingsDelta18oSurfaceEnum,(finaltime-(21000*yts))); + this->parameters->FindParam(&Delta18oSurfaceTime,SurfaceforcingsDelta18oSurfaceEnum,time); + /*Compute the temperature and precipitation*/ for(int iv=0;ivAddTimeInput(newmonthinput1,imonth/12.*yts); - TriaP1Input* newmonthinput2 = new TriaP1Input(SurfaceforcingsPrecipitationEnum,&monthlyprec[iv][imonth]); - NewPrecipitationInput->AddTimeInput(newmonthinput2,imonth/12.*yts); - } + TriaP1Input* newmonthinput1 = new TriaP1Input(SurfaceforcingsMonthlytemperaturesEnum,&TemperaturesPresentday[0][imonth]); + //TriaP1Input* newmonthinput1 = new TriaP1Input(SurfaceforcingsMonthlytemperaturesEnum,&monthlytemperatures[0][imonth]); + NewTemperatureInput->AddTimeInput(newmonthinput1,time+imonth/12.*yts); + TriaP1Input* newmonthinput2 = new TriaP1Input(SurfaceforcingsPrecipitationEnum,&monthlyprec[0][imonth]); + NewPrecipitationInput->AddTimeInput(newmonthinput2,time+imonth/12.*yts); } this->inputs->AddInput(NewTemperatureInput); this->inputs->AddInput(NewPrecipitationInput); @@ -1384,7 +1385,7 @@ void Tria::InputToResult(int enum_type,int step,IssmDouble time){ int i; - Input *input = NULL; + Input *input = NULL; /*Go through all the input objects, and find the one corresponding to enum_type, if it exists: */ if (enum_type==MaterialsRheologyBbarEnum) input=this->matice->inputs->GetInput(enum_type); @@ -2153,7 +2154,7 @@ IssmDouble agd[NUMVERTICES]; // surface mass balance IssmDouble monthlytemperatures[NUMVERTICES][12],monthlyprec[NUMVERTICES][12]; - IssmDouble h[NUMVERTICES],s[NUMVERTICES]; // ,b + IssmDouble h[NUMVERTICES],s[NUMVERTICES]; IssmDouble rho_water,rho_ice; /*Recover monthly temperatures and precipitation*/ @@ -2169,7 +2170,7 @@ input->GetInputValue(&monthlytemperatures[iv][month],gauss,time+month/12.*yts); monthlytemperatures[iv][month]=monthlytemperatures[iv][month]-273.15; // conversion from Kelvin to celcius input2->GetInputValue(&monthlyprec[iv][month],gauss,time+month/12.*yts); - monthlyprec[iv][month]=monthlyprec[iv][month]*yts; // convertion to m/yr + monthlyprec[iv][month]=monthlyprec[iv][month]*yts; // convertion in m/y } } @@ -2181,6 +2182,9 @@ rho_ice=matpar->GetRhoIce(); rho_water=matpar->GetRhoFreshwater(); + //if(id==1) printf(" monthlytemperatures %f\n",monthlytemperatures[1][0]); + //if(id==1) printf("monthlyprec %f\n",monthlyprec[1][0]); + /*measure the surface mass balance*/ for (int iv = 0; iv