Ignore:
Timestamp:
07/26/12 10:08:15 (13 years ago)
Author:
lemorzad
Message:

addind delta18o temperature and precipitation methode

File:
1 edited

Legend:

Unmodified
Added
Removed
  • issm/trunk-jpl/src/c/objects/Elements/Tria.cpp

    r12704 r12748  
    906906        this->results=new Results();
    907907
     908}
     909/*}}}*/
     910/*FUNCTION Tria::Delta18oParameterization{{{*/
     911void  Tria::Delta18oParameterization(void){
     912
     913        IssmDouble monthlytemperatures[NUMVERTICES][12],monthlyprec[NUMVERTICES][12];
     914        IssmDouble Delta18oTemperaturesPresentday[NUMVERTICES][12],Delta18oTemperaturesLgm[NUMVERTICES][12];
     915        IssmDouble PrecipitationsPresentday[NUMVERTICES][12];
     916        IssmDouble Delta18oPresent,Delta18oLgm,Delta18oTime;
     917        IssmDouble Delta18oSurfacePresent,Delta18oSurfaceLgm,Delta18oSurfaceTime;
     918        IssmDouble time,yts,finaltime;
     919        this->parameters->FindParam(&time,TimeEnum);
     920        this->parameters->FindParam(&yts,ConstantsYtsEnum);
     921        this->parameters->FindParam(&finaltime,TimesteppingFinalTimeEnum);
     922
     923        /*Recover present day temperature and precipitation*/
     924        Input*     input=inputs->GetInput(SurfaceforcingsDelta18oTemperaturesPresentdayEnum); _assert_(input);
     925        Input*     input2=inputs->GetInput(SurfaceforcingsDelta18oTemperaturesLgmEnum); _assert_(input2);
     926        Input*     input3=inputs->GetInput(SurfaceforcingsPrecipitationsPresentdayEnum); _assert_(input3);
     927        GaussTria* gauss=new GaussTria();
     928        for(int month=0;month<12;month++) {
     929                for(int iv=0;iv<NUMVERTICES;iv++) {
     930                        gauss->GaussVertex(iv);
     931                        input->GetInputValue(&Delta18oTemperaturesPresentday[iv][month],gauss,month/12.*yts);
     932                        input2->GetInputValue(&Delta18oTemperaturesLgm[iv][month],gauss,month/12.*yts);
     933                        input3->GetInputValue(&PrecipitationsPresentday[iv][month],gauss,month/12.*yts);
     934                        monthlyprec[iv][month]=monthlyprec[iv][month]*yts; // convertion to m/yr
     935                }
     936        }
     937
     938        /*Recover delta18o and Delta18oSurface at present day, lgm and at time t*/
     939        this->parameters->FindParam(&Delta18oPresent,SurfaceforcingsDelta18oEnum,finaltime*yts);
     940        this->parameters->FindParam(&Delta18oLgm,SurfaceforcingsDelta18oEnum,(finaltime-21000)*yts);
     941        this->parameters->FindParam(&Delta18oTime,SurfaceforcingsDelta18oEnum,time*yts);
     942        this->parameters->FindParam(&Delta18oSurfacePresent,SurfaceforcingsDelta18oSurfaceEnum,finaltime*yts);
     943        this->parameters->FindParam(&Delta18oSurfaceLgm,SurfaceforcingsDelta18oSurfaceEnum,(finaltime-21000)*yts);
     944        this->parameters->FindParam(&Delta18oSurfaceTime,SurfaceforcingsDelta18oSurfaceEnum,time*yts);
     945   
     946        /*Compute the temperature and precipitation*/
     947        for(int iv=0;iv<NUMVERTICES;iv++){
     948          ComputeDelta18oTemperaturePrecipitation(Delta18oSurfacePresent, Delta18oSurfaceLgm, Delta18oSurfaceTime,
     949                                              Delta18oPresent, Delta18oLgm, Delta18oTime,
     950                                              &PrecipitationsPresentday[iv][0],
     951                                              &Delta18oTemperaturesLgm[iv][0], &Delta18oTemperaturesPresentday[iv][0],
     952                                              &monthlytemperatures[iv][0], &monthlyprec[iv][0]);
     953        }
     954 
     955        /*Update inputs*/
     956        TransientInput* NewTemperatureInput = new TransientInput(SurfaceforcingsMonthlytemperaturesEnum);
     957        TransientInput* NewPrecipitationInput = new TransientInput(SurfaceforcingsPrecipitationEnum);
     958        for (int imonth=0;imonth<12;imonth++) {
     959                for(int iv=0;iv<NUMVERTICES;iv++) {
     960                        TriaP1Input* newmonthinput1 = new TriaP1Input(SurfaceforcingsMonthlytemperaturesEnum,&monthlytemperatures[iv][imonth]);
     961                        NewTemperatureInput->AddTimeInput(newmonthinput1,imonth/12.*yts);
     962                        TriaP1Input* newmonthinput2 = new TriaP1Input(SurfaceforcingsPrecipitationEnum,&monthlyprec[iv][imonth]);
     963                        NewPrecipitationInput->AddTimeInput(newmonthinput2,imonth/12.*yts);
     964                }
     965        }
     966        this->inputs->AddInput(NewTemperatureInput);
     967        this->inputs->AddInput(NewPrecipitationInput);
    908968}
    909969/*}}}*/
Note: See TracChangeset for help on using the changeset viewer.