Index: /issm/trunk-jpl/src/c/analyses/SmbAnalysis.cpp =================================================================== --- /issm/trunk-jpl/src/c/analyses/SmbAnalysis.cpp (revision 24625) +++ /issm/trunk-jpl/src/c/analyses/SmbAnalysis.cpp (revision 24626) @@ -114,8 +114,38 @@ iomodel->FetchDataToDatasetInput(inputs2,elements,"md.smb.precipitations_presentday",SmbPrecipitationsPresentdayEnum); if(!istemperaturescaled){ - iomodel->FetchDataToInput(inputs2,elements,"md.smb.temperatures_reconstructed",SmbTemperaturesReconstructedEnum); + /*Fetch array*/ + IssmDouble* doublearray = NULL; + int M,N; + iomodel->FetchData(&doublearray,&M,&N,"md.smb.temperatures_reconstructed"); + if(M!=iomodel->numberofvertices+1) _error_("md.smb.temperatures_reconstructed should have nbv+1 rows"); + if(N%12!=0) _error_("md.smb.temperatures_reconstructed should have a multiple of 12 columns (since it is monthly)"); + + /*Build indices*/ + int* ids = xNew(N); for(int i=0;iSize();i++){ + Element* element=xDynamicCast(elements->GetObjectByOffset(i)); + element->DatasetInputCreate(doublearray,M-1,N,ids,N,inputs2,iomodel,SmbTemperaturesReconstructedEnum); + } + xDelete(ids); + iomodel->DeleteData(doublearray,"md.smb.temperatures_reconstructed"); } if(!isprecipscaled){ - iomodel->FetchDataToInput(inputs2,elements,"md.smb.precipitations_reconstructed",SmbPrecipitationsReconstructedEnum); + /*Fetch array*/ + IssmDouble* doublearray = NULL; + int M,N; + iomodel->FetchData(&doublearray,&M,&N,"md.smb.precipitations_reconstructed"); + if(M!=iomodel->numberofvertices+1) _error_("md.smb.precipitations_reconstructed should have nbv+1 rows"); + if(N%12!=0) _error_("md.smb.precipitations_reconstructed should have a multiple of 12 columns (since it is monthly)"); + + /*Build indices*/ + int* ids = xNew(N); for(int i=0;iSize();i++){ + Element* element=xDynamicCast(elements->GetObjectByOffset(i)); + element->DatasetInputCreate(doublearray,M-1,N,ids,N,inputs2,iomodel,SmbPrecipitationsReconstructedEnum); + } + xDelete(ids); + iomodel->DeleteData(doublearray,"md.smb.precipitations_reconstructed"); } } @@ -248,4 +278,5 @@ iomodel->DeleteData(temp,"md.smb.delta18o_surface"); } + break; case SMBpddSicopolisEnum: @@ -266,5 +297,62 @@ parameters->AddObject(new TransientParam(SmbDelta18oEnum,&temp[0],&temp[M],interp,M)); iomodel->DeleteData(temp,"md.smb.delta18o"); - } + + IssmDouble yts; + bool istemperaturescaled,isprecipscaled; + iomodel->FindConstant(&yts,"md.constants.yts"); + iomodel->FindConstant(&istemperaturescaled,"md.smb.istemperaturescaled"); + iomodel->FindConstant(&isprecipscaled,"md.smb.isprecipscaled"); + if(!istemperaturescaled){ + /*Fetch array*/ + IssmDouble* doublearray = NULL; + int M,N; + iomodel->FetchData(&doublearray,&M,&N,"md.smb.temperatures_reconstructed"); + if(M!=iomodel->numberofvertices+1) _error_("md.smb.temperatures_reconstructed should have nbv+1 rows"); + if(N%12!=0) _error_("md.smb.temperatures_reconstructed should have a multiple of 12 columns (since it is monthly)"); + int numyears = N/12; _assert_(numyears*12==N); + + /*Check times*/ + #ifdef _ISSM_DEBUG_ + for(int i=0;idoublearray[(M-1)*N+i*12+j-1]); + } + } + #endif + + /*Build time*/ + IssmDouble* times = xNew(numyears); for(int i=0;iAddObject(new DoubleVecParam(SmbTemperaturesReconstructedYearsEnum,times,numyears)); + xDelete(times); + iomodel->DeleteData(doublearray,"md.smb.temperatures_reconstructed"); + } + if(!isprecipscaled){ + /*Fetch array*/ + IssmDouble* doublearray = NULL; + int M,N; + iomodel->FetchData(&doublearray,&M,&N,"md.smb.precipitations_reconstructed"); + if(M!=iomodel->numberofvertices+1) _error_("md.smb.precipitations_reconstructed should have nbv+1 rows"); + if(N%12!=0) _error_("md.smb.precipitations_reconstructed should have a multiple of 12 columns (since it is monthly)"); + int numyears = N/12; _assert_(numyears*12==N); + + /*Check times*/ + #ifdef _ISSM_DEBUG_ + for(int i=0;idoublearray[(M-1)*N+i*12+j-1]); + } + } + #endif + + /*Build time*/ + IssmDouble* times = xNew(numyears); for(int i=0;iAddObject(new DoubleVecParam(SmbPrecipitationsReconstructedYearsEnum,times,numyears)); + xDelete(times); + iomodel->DeleteData(doublearray,"md.smb.precipitations_reconstructed"); + } + } + break; case SMBgradientsEnum: Index: /issm/trunk-jpl/src/c/classes/Elements/Element.cpp =================================================================== --- /issm/trunk-jpl/src/c/classes/Elements/Element.cpp (revision 24625) +++ /issm/trunk-jpl/src/c/classes/Elements/Element.cpp (revision 24626) @@ -563,15 +563,20 @@ /*Recover present day temperature and precipitation*/ - TransientInput2 *tinput3 = NULL; - TransientInput2 *tinput4 = NULL; + DatasetInput2 *dinput3 = NULL; + DatasetInput2 *dinput4 = NULL; + int offset_t,offset_p,N; if(!isTemperatureScaled){ - tinput3=this->inputs2->GetTransientInput(SmbTemperaturesReconstructedEnum); _assert_(tinput3); - offset = tinput3->GetTimeInputOffset(time_climt); - time_climt = tinput3->GetTimeByOffset(offset-11-fmod(offset-11,12.)); + IssmDouble* time_temp_scaled = NULL; + parameters->FindParam(&time_temp_scaled,&N,SmbTemperaturesReconstructedYearsEnum); + if(!binary_search(&offset_t,time_climt,time_temp_scaled,N)) _error_("time not sorted?"); + xDelete(time_temp_scaled); + dinput3=this->GetDatasetInput2(SmbTemperaturesReconstructedEnum); _assert_(dinput3); } if(!isPrecipScaled){ - tinput4=this->inputs2->GetTransientInput(SmbPrecipitationsReconstructedEnum); _assert_(tinput4); - offset = tinput4->GetTimeInputOffset(time_climp); - time_climp = tinput4->GetTimeByOffset(offset-11-fmod(offset-11,12.)); + IssmDouble* time_precip_scaled = NULL; + parameters->FindParam(&time_precip_scaled,&N,SmbPrecipitationsReconstructedYearsEnum); + if(!binary_search(&offset_t,time_climt,time_precip_scaled,N)) _error_("time not sorted?"); + xDelete(time_precip_scaled); + dinput4=this->GetDatasetInput2(SmbPrecipitationsReconstructedEnum); _assert_(dinput4); } @@ -583,13 +588,4 @@ Gauss* gauss=this->NewGauss(); for(int month=0;month<12;month++) { - - Input2 *input3 = NULL; - Input2 *input4 = NULL; - if(!isTemperatureScaled){ - input3=this->GetInput2(SmbTemperaturesReconstructedEnum,time_climt+month/12.*yts); _assert_(input3); - } - if(!isPrecipScaled){ - input4=this->GetInput2(SmbPrecipitationsReconstructedEnum,time_climp+month/12.*yts); - } for(int iv=0;ivGaussVertex(iv); @@ -599,8 +595,8 @@ if(!isTemperatureScaled){ - input3->GetInputValue(&TemperaturesReconstructed[iv*12+month],gauss); + dinput3->GetInputValue(&TemperaturesReconstructed[iv*12+month],gauss,offset_t*12+month); } if(!isPrecipScaled){ - input4->GetInputValue(&PrecipitationsReconstructed[iv*12+month],gauss); + dinput4->GetInputValue(&PrecipitationsReconstructed[iv*12+month],gauss,offset_p*12+month); PrecipitationsReconstructed[iv*12+month]=PrecipitationsReconstructed[iv*12+month]*yts; } Index: /issm/trunk-jpl/src/c/shared/Enum/Enum.vim =================================================================== --- /issm/trunk-jpl/src/c/shared/Enum/Enum.vim (revision 24625) +++ /issm/trunk-jpl/src/c/shared/Enum/Enum.vim (revision 24626) @@ -383,4 +383,6 @@ syn keyword cConstant SmbTdiffEnum syn keyword cConstant SmbThermoDeltaTScalingEnum +syn keyword cConstant SmbTemperaturesReconstructedYearsEnum +syn keyword cConstant SmbPrecipitationsReconstructedYearsEnum syn keyword cConstant SmoothThicknessMultiplierEnum syn keyword cConstant SolutionTypeEnum Index: /issm/trunk-jpl/src/c/shared/Enum/EnumDefinitions.h =================================================================== --- /issm/trunk-jpl/src/c/shared/Enum/EnumDefinitions.h (revision 24625) +++ /issm/trunk-jpl/src/c/shared/Enum/EnumDefinitions.h (revision 24626) @@ -377,4 +377,6 @@ SmbTdiffEnum, SmbThermoDeltaTScalingEnum, + SmbTemperaturesReconstructedYearsEnum, + SmbPrecipitationsReconstructedYearsEnum, SmoothThicknessMultiplierEnum, SolutionTypeEnum, Index: /issm/trunk-jpl/src/c/shared/Enum/EnumToStringx.cpp =================================================================== --- /issm/trunk-jpl/src/c/shared/Enum/EnumToStringx.cpp (revision 24625) +++ /issm/trunk-jpl/src/c/shared/Enum/EnumToStringx.cpp (revision 24626) @@ -385,4 +385,6 @@ case SmbTdiffEnum : return "SmbTdiff"; case SmbThermoDeltaTScalingEnum : return "SmbThermoDeltaTScaling"; + case SmbTemperaturesReconstructedYearsEnum : return "SmbTemperaturesReconstructedYears"; + case SmbPrecipitationsReconstructedYearsEnum : return "SmbPrecipitationsReconstructedYears"; case SmoothThicknessMultiplierEnum : return "SmoothThicknessMultiplier"; case SolutionTypeEnum : return "SolutionType"; Index: /issm/trunk-jpl/src/c/shared/Enum/StringToEnumx.cpp =================================================================== --- /issm/trunk-jpl/src/c/shared/Enum/StringToEnumx.cpp (revision 24625) +++ /issm/trunk-jpl/src/c/shared/Enum/StringToEnumx.cpp (revision 24626) @@ -394,4 +394,6 @@ else if (strcmp(name,"SmbTdiff")==0) return SmbTdiffEnum; else if (strcmp(name,"SmbThermoDeltaTScaling")==0) return SmbThermoDeltaTScalingEnum; + else if (strcmp(name,"SmbTemperaturesReconstructedYears")==0) return SmbTemperaturesReconstructedYearsEnum; + else if (strcmp(name,"SmbPrecipitationsReconstructedYears")==0) return SmbPrecipitationsReconstructedYearsEnum; else if (strcmp(name,"SmoothThicknessMultiplier")==0) return SmoothThicknessMultiplierEnum; else if (strcmp(name,"SolutionType")==0) return SolutionTypeEnum; @@ -504,10 +506,10 @@ else if (strcmp(name,"CalvinglevermannCoeff")==0) return CalvinglevermannCoeffEnum; else if (strcmp(name,"CalvingratexAverage")==0) return CalvingratexAverageEnum; - else if (strcmp(name,"Calvingratex")==0) return CalvingratexEnum; - else if (strcmp(name,"CalvingrateyAverage")==0) return CalvingrateyAverageEnum; else stage=5; } if(stage==5){ - if (strcmp(name,"Calvingratey")==0) return CalvingrateyEnum; + if (strcmp(name,"Calvingratex")==0) return CalvingratexEnum; + else if (strcmp(name,"CalvingrateyAverage")==0) return CalvingrateyAverageEnum; + else if (strcmp(name,"Calvingratey")==0) return CalvingrateyEnum; else if (strcmp(name,"CalvingFluxLevelset")==0) return CalvingFluxLevelsetEnum; else if (strcmp(name,"CalvingMeltingFluxLevelset")==0) return CalvingMeltingFluxLevelsetEnum; @@ -627,10 +629,10 @@ else if (strcmp(name,"InversionVelObs")==0) return InversionVelObsEnum; else if (strcmp(name,"InversionVxObs")==0) return InversionVxObsEnum; - else if (strcmp(name,"InversionVyObs")==0) return InversionVyObsEnum; - else if (strcmp(name,"LevelsetfunctionSlopeX")==0) return LevelsetfunctionSlopeXEnum; else stage=6; } if(stage==6){ - if (strcmp(name,"LevelsetfunctionSlopeY")==0) return LevelsetfunctionSlopeYEnum; + if (strcmp(name,"InversionVyObs")==0) return InversionVyObsEnum; + else if (strcmp(name,"LevelsetfunctionSlopeX")==0) return LevelsetfunctionSlopeXEnum; + else if (strcmp(name,"LevelsetfunctionSlopeY")==0) return LevelsetfunctionSlopeYEnum; else if (strcmp(name,"LoadingforceX")==0) return LoadingforceXEnum; else if (strcmp(name,"LoadingforceY")==0) return LoadingforceYEnum; @@ -750,10 +752,10 @@ else if (strcmp(name,"SmbMeanSHF")==0) return SmbMeanSHFEnum; else if (strcmp(name,"SmbMeanULW")==0) return SmbMeanULWEnum; - else if (strcmp(name,"SmbMelt")==0) return SmbMeltEnum; - else if (strcmp(name,"SmbMonthlytemperatures")==0) return SmbMonthlytemperaturesEnum; else stage=7; } if(stage==7){ - if (strcmp(name,"SmbNetLW")==0) return SmbNetLWEnum; + if (strcmp(name,"SmbMelt")==0) return SmbMeltEnum; + else if (strcmp(name,"SmbMonthlytemperatures")==0) return SmbMonthlytemperaturesEnum; + else if (strcmp(name,"SmbNetLW")==0) return SmbNetLWEnum; else if (strcmp(name,"SmbNetSW")==0) return SmbNetSWEnum; else if (strcmp(name,"SmbPAir")==0) return SmbPAirEnum; @@ -873,10 +875,10 @@ else if (strcmp(name,"Outputdefinition20")==0) return Outputdefinition20Enum; else if (strcmp(name,"Outputdefinition21")==0) return Outputdefinition21Enum; - else if (strcmp(name,"Outputdefinition22")==0) return Outputdefinition22Enum; - else if (strcmp(name,"Outputdefinition23")==0) return Outputdefinition23Enum; else stage=8; } if(stage==8){ - if (strcmp(name,"Outputdefinition24")==0) return Outputdefinition24Enum; + if (strcmp(name,"Outputdefinition22")==0) return Outputdefinition22Enum; + else if (strcmp(name,"Outputdefinition23")==0) return Outputdefinition23Enum; + else if (strcmp(name,"Outputdefinition24")==0) return Outputdefinition24Enum; else if (strcmp(name,"Outputdefinition25")==0) return Outputdefinition25Enum; else if (strcmp(name,"Outputdefinition26")==0) return Outputdefinition26Enum; @@ -996,10 +998,10 @@ else if (strcmp(name,"Boundary")==0) return BoundaryEnum; else if (strcmp(name,"BuddJacka")==0) return BuddJackaEnum; - else if (strcmp(name,"CalvingDev2")==0) return CalvingDev2Enum; - else if (strcmp(name,"CalvingHab")==0) return CalvingHabEnum; else stage=9; } if(stage==9){ - if (strcmp(name,"CalvingLevermann")==0) return CalvingLevermannEnum; + if (strcmp(name,"CalvingDev2")==0) return CalvingDev2Enum; + else if (strcmp(name,"CalvingHab")==0) return CalvingHabEnum; + else if (strcmp(name,"CalvingLevermann")==0) return CalvingLevermannEnum; else if (strcmp(name,"CalvingVonmises")==0) return CalvingVonmisesEnum; else if (strcmp(name,"Cfdragcoeffabsgrad")==0) return CfdragcoeffabsgradEnum; @@ -1119,10 +1121,10 @@ else if (strcmp(name,"Indexed")==0) return IndexedEnum; else if (strcmp(name,"IntExternalResult")==0) return IntExternalResultEnum; - else if (strcmp(name,"IntInput")==0) return IntInputEnum; - else if (strcmp(name,"ElementInput2")==0) return ElementInput2Enum; else stage=10; } if(stage==10){ - if (strcmp(name,"SegInput2")==0) return SegInput2Enum; + if (strcmp(name,"IntInput")==0) return IntInputEnum; + else if (strcmp(name,"ElementInput2")==0) return ElementInput2Enum; + else if (strcmp(name,"SegInput2")==0) return SegInput2Enum; else if (strcmp(name,"TriaInput2")==0) return TriaInput2Enum; else if (strcmp(name,"PentaInput2")==0) return PentaInput2Enum; @@ -1242,10 +1244,10 @@ else if (strcmp(name,"Riftfront")==0) return RiftfrontEnum; else if (strcmp(name,"SIAApproximation")==0) return SIAApproximationEnum; - else if (strcmp(name,"SMBcomponents")==0) return SMBcomponentsEnum; - else if (strcmp(name,"SMBd18opdd")==0) return SMBd18opddEnum; else stage=11; } if(stage==11){ - if (strcmp(name,"SMBforcing")==0) return SMBforcingEnum; + if (strcmp(name,"SMBcomponents")==0) return SMBcomponentsEnum; + else if (strcmp(name,"SMBd18opdd")==0) return SMBd18opddEnum; + else if (strcmp(name,"SMBforcing")==0) return SMBforcingEnum; else if (strcmp(name,"SMBgcm")==0) return SMBgcmEnum; else if (strcmp(name,"SMBgemb")==0) return SMBgembEnum;