Changeset 12748 for issm/trunk-jpl/src/c/objects/Elements/Tria.cpp

Timestamp:

07/26/12 10:08:15 (13 years ago)

Author:

lemorzad

Message:

addind delta18o temperature and precipitation methode

File:

• issm/trunk-jpl/src/c/objects/Elements/Tria.cpp (modified) (1 diff)

issm/trunk-jpl/src/c/objects/Elements/Tria.cpp

@@ -906 +906 @@
         this->results=new Results();
 
+}
+/*}}}*/
+/*FUNCTION Tria::Delta18oParameterization{{{*/
+void  Tria::Delta18oParameterization(void){
+
+        IssmDouble monthlytemperatures[NUMVERTICES][12],monthlyprec[NUMVERTICES][12];
+        IssmDouble Delta18oTemperaturesPresentday[NUMVERTICES][12],Delta18oTemperaturesLgm[NUMVERTICES][12];
+        IssmDouble PrecipitationsPresentday[NUMVERTICES][12];
+        IssmDouble Delta18oPresent,Delta18oLgm,Delta18oTime;
+        IssmDouble Delta18oSurfacePresent,Delta18oSurfaceLgm,Delta18oSurfaceTime;
+        IssmDouble time,yts,finaltime;
+        this->parameters->FindParam(&time,TimeEnum);
+        this->parameters->FindParam(&yts,ConstantsYtsEnum);
+        this->parameters->FindParam(&finaltime,TimesteppingFinalTimeEnum);
+
+        /*Recover present day temperature and precipitation*/
+        Input*     input=inputs->GetInput(SurfaceforcingsDelta18oTemperaturesPresentdayEnum); _assert_(input);
+        Input*     input2=inputs->GetInput(SurfaceforcingsDelta18oTemperaturesLgmEnum); _assert_(input2);
+        Input*     input3=inputs->GetInput(SurfaceforcingsPrecipitationsPresentdayEnum); _assert_(input3);
+        GaussTria* gauss=new GaussTria();
+        for(int month=0;month<12;month++) {
+                for(int iv=0;iv<NUMVERTICES;iv++) {
+                        gauss->GaussVertex(iv);
+                        input->GetInputValue(&Delta18oTemperaturesPresentday[iv][month],gauss,month/12.*yts);
+                        input2->GetInputValue(&Delta18oTemperaturesLgm[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
+                }
+        }
+
+        /*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);
+   
+        /*Compute the temperature and precipitation*/
+        for(int iv=0;iv<NUMVERTICES;iv++){
+          ComputeDelta18oTemperaturePrecipitation(Delta18oSurfacePresent, Delta18oSurfaceLgm, Delta18oSurfaceTime, 
+                                              Delta18oPresent, Delta18oLgm, Delta18oTime,
+                                              &PrecipitationsPresentday[iv][0], 
+                                              &Delta18oTemperaturesLgm[iv][0], &Delta18oTemperaturesPresentday[iv][0], 
+                                              &monthlytemperatures[iv][0], &monthlyprec[iv][0]);
+        }
+  
+        /*Update inputs*/ 
+        TransientInput* NewTemperatureInput = new TransientInput(SurfaceforcingsMonthlytemperaturesEnum);
+        TransientInput* NewPrecipitationInput = new TransientInput(SurfaceforcingsPrecipitationEnum);
+        for (int imonth=0;imonth<12;imonth++) {
+                for(int iv=0;iv<NUMVERTICES;iv++) {
+                        TriaP1Input* newmonthinput1 = new TriaP1Input(SurfaceforcingsMonthlytemperaturesEnum,&monthlytemperatures[iv][imonth]);
+                        NewTemperatureInput->AddTimeInput(newmonthinput1,imonth/12.*yts);
+                        TriaP1Input* newmonthinput2 = new TriaP1Input(SurfaceforcingsPrecipitationEnum,&monthlyprec[iv][imonth]);
+                        NewPrecipitationInput->AddTimeInput(newmonthinput2,imonth/12.*yts);
+                }
+        }
+        this->inputs->AddInput(NewTemperatureInput);
+        this->inputs->AddInput(NewPrecipitationInput);
 }
 /*}}}*/