Index: /issm/trunk-jpl/src/c/classes/Elements/Penta.cpp
===================================================================
--- /issm/trunk-jpl/src/c/classes/Elements/Penta.cpp	(revision 19136)
+++ /issm/trunk-jpl/src/c/classes/Elements/Penta.cpp	(revision 19137)
@@ -2166,8 +2166,10 @@
    IssmDouble agd[NUMVERTICES];             // surface mass balance
    IssmDouble monthlytemperatures[NUMVERTICES][12],monthlyprec[NUMVERTICES][12];
+   IssmDouble yearlytemperatures[NUMVERTICES]={0.};
    IssmDouble tmp[NUMVERTICES];
    IssmDouble h[NUMVERTICES],s[NUMVERTICES];
    IssmDouble rho_water,rho_ice,desfac,s0p,s0t,rlaps,rlapslgm;
    IssmDouble PfacTime,TdiffTime,sealevTime;
+   IssmDouble mavg=1./12.; //factor for monthly average
 
    /*Get material parameters :*/
@@ -2194,5 +2196,6 @@
        gauss->GaussVertex(iv);
        input->GetInputValue(&monthlytemperatures[iv][month],gauss,time+month/12.*yts);
-       monthlytemperatures[iv][month]=monthlytemperatures[iv][month]-273.15; // conversion from Kelvin to celcius
+       yearlytemperatures[iv]=yearlytemperatures[iv]+monthlytemperatures[iv][month]*mavg; // Has to be in Kelvin 
+       monthlytemperatures[iv][month]=monthlytemperatures[iv][month]-273.15; // conversion from Kelvin to celcius for PDD module
        input2->GetInputValue(&monthlyprec[iv][month],gauss,time+month/12.*yts);
      }
@@ -2241,4 +2244,5 @@
 
 
+   this->inputs->AddInput(new PentaInput(TemperatureEnum,&yearlytemperatures[0],P1Enum));
    this->inputs->AddInput(new PentaInput(SurfaceforcingsMassBalanceEnum,&agd[0],P1Enum));
    // this->inputs->AddInput(NewTemperatureInput);
@@ -2246,4 +2250,5 @@
    // //this->inputs->AddInput(new PentaVertexInput(ThermalSpcTemperatureEnum,&Tsurf[0]));
     this->InputExtrude(SurfaceforcingsMassBalanceEnum,-1);
+    this->InputExtrude(TemperatureEnum,-1);
    // this->InputExtrude(SurfaceforcingsMonthlytemperaturesEnum,-1);
    // this->InputExtrude(SurfaceforcingsPrecipitationEnum,-1);
Index: /issm/trunk-jpl/src/c/classes/Elements/Tria.cpp
===================================================================
--- /issm/trunk-jpl/src/c/classes/Elements/Tria.cpp	(revision 19136)
+++ /issm/trunk-jpl/src/c/classes/Elements/Tria.cpp	(revision 19137)
@@ -2574,8 +2574,10 @@
    IssmDouble agd[NUMVERTICES];             // surface mass balance
    IssmDouble monthlytemperatures[NUMVERTICES][12],monthlyprec[NUMVERTICES][12];
+   IssmDouble yearlytemperatures[NUMVERTICES]={0.};
    IssmDouble tmp[NUMVERTICES];
    IssmDouble h[NUMVERTICES],s[NUMVERTICES];
    IssmDouble rho_water,rho_ice,desfac,s0p,s0t,rlaps,rlapslgm;
    IssmDouble PfacTime,TdiffTime,sealevTime;
+   IssmDouble mavg=1./12.; //factor for monthly average
    
    /*Get material parameters :*/
@@ -2590,5 +2592,5 @@
   rlapslgm=matpar->GetMaterialParameter(SurfaceforcingsRlapslgmEnum);
 
-   /*Recover monthly temperatures and precipitation and Present day and LGm ones*/
+   /*Recover monthly temperatures and precipitation and compute the yearly mean temperatures*/
    Input*     input=inputs->GetInput(SurfaceforcingsMonthlytemperaturesEnum); _assert_(input);
    Input*     input2=inputs->GetInput(SurfaceforcingsPrecipitationEnum); _assert_(input2);
@@ -2598,13 +2600,13 @@
    this->parameters->FindParam(&yts,ConstantsYtsEnum);
 
-
    for(int month=0;month<12;month++) {
      for(int iv=0;iv<NUMVERTICES;iv++) {
        gauss->GaussVertex(iv);
        input->GetInputValue(&monthlytemperatures[iv][month],gauss,time+month/12.*yts);
-       monthlytemperatures[iv][month]=monthlytemperatures[iv][month]-273.15; // conversion from Kelvin to celcius
+       yearlytemperatures[iv]=yearlytemperatures[iv]+monthlytemperatures[iv][month]*mavg; // Has to be in Kelvin 
+       monthlytemperatures[iv][month]=monthlytemperatures[iv][month]-273.15; // conversion from Kelvin to celcius for PDD module
        input2->GetInputValue(&monthlyprec[iv][month],gauss,time+month/12.*yts);
      }
-   }
+   }     
 
   /*Recover Pfac, Tdiff and sealev at time t:
@@ -2626,9 +2628,9 @@
    /*measure the surface mass balance*/
   for (int iv = 0; iv<NUMVERTICES; iv++){
-     agd[iv]=PddSurfaceMassBalance(&monthlytemperatures[iv][0], &monthlyprec[iv][0],
+    agd[iv]=PddSurfaceMassBalance(&monthlytemperatures[iv][0], &monthlyprec[iv][0],
 				  pdds, pds, signorm, yts, h[iv], s[iv],
 				  desfac, s0t, s0p,rlaps,rlapslgm,TdiffTime,sealevTime,
 				  rho_water,rho_ice);
-   }
+  }
 
    /*Update inputs*/    
@@ -2648,5 +2650,6 @@
 
 
-   this->inputs->AddInput(new TriaInput(SurfaceforcingsMassBalanceEnum,&agd[0],P1Enum));
+  this->inputs->AddInput(new TriaInput(TemperatureEnum,&yearlytemperatures[0],P1Enum));
+  this->inputs->AddInput(new TriaInput(SurfaceforcingsMassBalanceEnum,&agd[0],P1Enum));
    // this->inputs->AddInput(NewTemperatureInput);
    // this->inputs->AddInput(NewPrecipitationInput);
