WSL/SLF GitLab Repository

Commit 4f51c789 authored by Mathias Bavay's avatar Mathias Bavay
Browse files

More size_t fixes (ie trying to be consistent) and a bug was found when...

More size_t fixes (ie trying to be consistent) and a bug was found when spatially interpolating using a few specific algorithms that would not reset their internal vectors between two calls (thanks Matteo for helping uncover this bug!)
parent 6ecfe7bc
......@@ -35,7 +35,7 @@ IOPlugin& IOPlugin::operator=(const IOPlugin& source)
const std::string IOPlugin::toString() const {
std::stringstream os;
const unsigned int pt_w=8;
const int pt_w=8;
os << "<IOPlugin>" << std::setw(10) << classname;
os << "," << std::showbase << std::setw(pt_w+2) << io;
os << "," << std::showbase << std::setw(pt_w+2) << creator_func << "</IOPlugin>\n";
......
......@@ -788,7 +788,7 @@ std::string printFractionalDay(const double& fractional) {
return tmp.str();
}
void getArraySliceParams(const unsigned int& dimx, const unsigned int& nbworkers, const unsigned int &wk, unsigned int& startx, unsigned int& nx)
void getArraySliceParams(const size_t& dimx, const unsigned int& nbworkers, const unsigned int &wk, size_t& startx, size_t& nx)
{
if(nbworkers>dimx) {
std::stringstream ss;
......@@ -796,8 +796,8 @@ void getArraySliceParams(const unsigned int& dimx, const unsigned int& nbworkers
throw InvalidArgumentException(ss.str(), AT);
}
const unsigned int nx_avg = dimx / nbworkers;
const unsigned int remainder = dimx % nbworkers;
const size_t nx_avg = dimx / nbworkers;
const size_t remainder = dimx % nbworkers;
if(wk<=remainder) { //distribute remainder, 1 extra column per worker, on first workers
nx = nx_avg+1;
......
......@@ -374,7 +374,7 @@ namespace IOUtils {
* @param[out] startx calculated start index for the current slice
* @param[out] nx calculated number of cells (in the desired dimension) of the current slice
*/
void getArraySliceParams(const unsigned int& dimx, const unsigned int& nbworkers, const unsigned int &wk, unsigned int& startx, unsigned int& nx);
void getArraySliceParams(const size_t& dimx, const unsigned int& nbworkers, const unsigned int &wk, size_t& startx, size_t& nx);
/**
* @class file_indexer
......
......@@ -382,8 +382,8 @@ double RHAlgorithm::getQualityRating(const Date& i_date, const MeteoData::Parame
date = i_date;
param = in_param;
vecData.clear();
vecMeta.clear();
vecData.clear(); vecMeta.clear();
vecDataTA.clear(); vecDataRH.clear();
nrOfMeasurments = 0;
iomanager.getMeteoData(date, vecMeteo);
......@@ -396,7 +396,7 @@ double RHAlgorithm::getQualityRating(const Date& i_date, const MeteoData::Parame
}
}
if (vecDataTA.empty() || nrOfMeasurments==0)
if (nrOfMeasurments==0)
return 0.0;
if( (nrOfMeasurments<vecDataRH.size()/2) || ( nrOfMeasurments<2 ) )
return 0.6;
......@@ -414,10 +414,9 @@ void RHAlgorithm::calculate(const DEMObject& dem, Grid2DObject& grid)
Grid2DObject ta;
mi.interpolate(date, dem, MeteoData::TA, ta); //get TA interpolation from call back to Meteo2DInterpolator
//here, RH->Td, interpolations, Td->RH
std::vector<double> vecTd(vecDataRH.size(), 0.0); // init to 0.0
//RH->Td, interpolations, Td->RH
//Compute dew point temperatures at stations
std::vector<double> vecTd(vecDataRH.size());
for (size_t ii=0; ii<vecDataRH.size(); ii++){
vecTd[ii] = Atmosphere::RhtoDewPoint(vecDataRH[ii], vecDataTA[ii], 1);
}
......@@ -430,8 +429,8 @@ void RHAlgorithm::calculate(const DEMObject& dem, Grid2DObject& grid)
retrend(dem, trend, grid);
//Recompute Rh from the interpolated td
for (unsigned int jj=0; jj<grid.nrows; jj++) {
for (unsigned int ii=0; ii<grid.ncols; ii++) {
for (size_t jj=0; jj<grid.nrows; jj++) {
for (size_t ii=0; ii<grid.ncols; ii++) {
double &value = grid.grid2D(ii,jj);
if(value!=IOUtils::nodata)
value = Atmosphere::DewPointtoRh(value, ta.grid2D(ii,jj), 1);
......@@ -447,8 +446,8 @@ double ILWRAlgorithm::getQualityRating(const Date& i_date, const MeteoData::Para
date = i_date;
param = in_param;
vecData.clear();
vecMeta.clear();
vecData.clear(); vecMeta.clear();
vecDataEA.clear();
nrOfMeasurments = 0;
for (size_t ii=0; ii<vecMeteo.size(); ii++){
......@@ -459,7 +458,7 @@ double ILWRAlgorithm::getQualityRating(const Date& i_date, const MeteoData::Para
}
}
if (vecDataEA.empty() || nrOfMeasurments==0)
if (nrOfMeasurments==0)
return 0.0;
return 0.9;
......@@ -483,8 +482,8 @@ void ILWRAlgorithm::calculate(const DEMObject& dem, Grid2DObject& grid)
retrend(dem, trend, grid);
//Recompute Rh from the interpolated td
for (unsigned int jj=0; jj<grid.nrows; jj++) {
for (unsigned int ii=0; ii<grid.ncols; ii++) {
for (size_t jj=0; jj<grid.nrows; jj++) {
for (size_t ii=0; ii<grid.ncols; ii++) {
double &value = grid.grid2D(ii,jj);
value = Atmosphere::blkBody_Radiation(value, ta.grid2D(ii,jj));
}
......@@ -499,8 +498,8 @@ double SimpleWindInterpolationAlgorithm::getQualityRating(const Date& i_date, co
date = i_date;
param = in_param;
vecData.clear();
vecMeta.clear();
vecData.clear(); vecMeta.clear();
vecDataVW.clear(); vecDataDW.clear();
nrOfMeasurments = 0;
iomanager.getMeteoData(date, vecMeteo);
......@@ -513,7 +512,7 @@ double SimpleWindInterpolationAlgorithm::getQualityRating(const Date& i_date, co
}
}
if (vecDataVW.empty() || vecDataDW.empty() || nrOfMeasurments==0)
if (nrOfMeasurments==0)
return 0.0;
if( (nrOfMeasurments<vecDataVW.size()/2) || ( nrOfMeasurments<2 ) )
return 0.6;
......@@ -543,18 +542,18 @@ void SimpleWindInterpolationAlgorithm::calculate(const DEMObject& dem, Grid2DObj
std::string USERInterpolation::getGridFileName() const
{
//HACK: use read2DGrid(grid, MeteoGrid::Parameters, Date) instead?
const std::string ext=std::string(".asc");
if (vecArgs.size() != 1){
throw InvalidArgumentException("Please provide the path to the grids for the USER interpolation algorithm", AT);
}
const std::string ext(".asc");
const std::string& grid_path = vecArgs[0];
std::string gridname = grid_path + std::string("/");
std::string gridname = grid_path + "/";
if(!vecMeteo.empty()) {
const Date& timestep = vecMeteo.at(0).date;
gridname = gridname + timestep.toString(Date::NUM) + std::string("_") + MeteoData::getParameterName(param) + ext;
gridname = gridname + timestep.toString(Date::NUM) + "_" + MeteoData::getParameterName(param) + ext;
} else {
gridname = gridname + std::string("Default") + std::string("_") + MeteoData::getParameterName(param) + ext;
gridname = gridname + "Default" + "_" + MeteoData::getParameterName(param) + ext;
}
return gridname;
......@@ -564,8 +563,7 @@ double USERInterpolation::getQualityRating(const Date& i_date, const MeteoData::
{
date = i_date;
param = in_param;
nrOfMeasurments = 0;
const std::string filename = getGridFileName();
filename = getGridFileName();
if (!IOUtils::validFileName(filename)) {
cerr << "[E] Invalid grid filename for USER interpolation algorithm: " << filename << "\n";
......@@ -580,8 +578,6 @@ double USERInterpolation::getQualityRating(const Date& i_date, const MeteoData::
void USERInterpolation::calculate(const DEMObject& dem, Grid2DObject& grid)
{
const std::string filename = getGridFileName();
iomanager.read2DGrid(grid, filename);
if(!grid.isSameGeolocalization(dem)) {
throw InvalidArgumentException("[E] trying to load a grid(" + filename + ") that has not the same georeferencing as the DEM!", AT);
......@@ -640,13 +636,13 @@ void OrdinaryKrigingAlgorithm::getDataForVariogram(std::vector<double> &distData
distData.clear();
variData.clear();
for(unsigned int j=0; j<nrOfMeasurments; j++) {
for(size_t j=0; j<nrOfMeasurments; j++) {
const Coords& st1 = vecMeta[j].position;
const double x1 = st1.getEasting();
const double y1 = st1.getNorthing();
const double val1 = vecData[j];
for(unsigned int i=0; i<j; i++) {
for(size_t i=0; i<j; i++) {
//compute distance between stations
const Coords& st2 = vecMeta[i].position;
const double val2 = vecData[i];
......
......@@ -372,11 +372,12 @@ class USERInterpolation : public InterpolationAlgorithm {
USERInterpolation(Meteo2DInterpolator& i_mi,
const std::vector<std::string>& i_vecArgs,
const std::string& i_algo, IOManager& iom)
: InterpolationAlgorithm(i_mi, i_vecArgs, i_algo, iom) {}
: InterpolationAlgorithm(i_mi, i_vecArgs, i_algo, iom), filename() {nrOfMeasurments=0;}
virtual double getQualityRating(const Date& i_date, const MeteoData::Parameters& in_param);
virtual void calculate(const DEMObject& dem, Grid2DObject& grid);
private:
std::string getGridFileName() const;
std::string filename;
};
/**
......
......@@ -22,7 +22,7 @@ using namespace std;
namespace mio {
Meteo1DInterpolator::Meteo1DInterpolator(const Config& in_cfg)
: cfg(in_cfg), window_size(2.*86400.),
: cfg(in_cfg), window_size(86400.),
mapAlgorithms()
{
//default window_size is 2 julian days
......@@ -33,7 +33,7 @@ Meteo1DInterpolator::Meteo1DInterpolator(const Config& in_cfg)
//read the Config object to create the resampling algorithms for each
//MeteoData::Parameters parameter (i.e. each member variable like ta, p, hnw, ...)
for (unsigned int ii=0; ii<MeteoData::nrOfParameters; ii++){ //loop over all MeteoData member variables
for (size_t ii=0; ii<MeteoData::nrOfParameters; ii++){ //loop over all MeteoData member variables
const std::string parname = MeteoData::getParameterName(ii); //Current parameter name
vector<string> vecArgs;
const string algo_name = getInterpolationForParameter(parname, vecArgs);
......
......@@ -236,7 +236,7 @@ const std::string Meteo2DInterpolator::toString() const {
std::map<std::string, std::vector<InterpolationAlgorithm*> >::const_iterator iter;
for (iter = mapAlgorithms.begin(); iter != mapAlgorithms.end(); ++iter) {
os << setw(10) << iter->first << "::";
for(unsigned int jj=0; jj<iter->second.size(); jj++) {
for(size_t jj=0; jj<iter->second.size(); jj++) {
os << iter->second[jj]->algo << " ";
}
os << "\n";
......
......@@ -103,50 +103,50 @@ void fulfillDoubleArray(const Grid2DObject& p,
dest[5] = 1.; //reserved ...
if(cellOrder == "llur"){
unsigned int knx = 6;
unsigned int iknx = 0;
for (unsigned int kk = 0; kk < p.nrows; kk++){
for (unsigned int ll=0; ll < p.ncols; ll++)
size_t knx = 6;
size_t iknx = 0;
for (size_t kk = 0; kk < p.nrows; kk++){
for (size_t ll=0; ll < p.ncols; ll++)
dest[knx + ll] = p.grid2D(ll+iknx);
knx+=p.ncols;
iknx+=p.ncols;
}
}
else if(cellOrder == "urll" ){
unsigned int knx = (p.nrows-1)*p.ncols+6;
unsigned int iknx = 0;
for (unsigned int kk = p.nrows-1; kk >=0; kk--){
for (unsigned int ll=p.ncols -1; ll >=0; ll--)
size_t knx = (p.nrows-1)*p.ncols+6;
size_t iknx = 0;
for (size_t kk = p.nrows-1; kk >=0; kk--){
for (size_t ll=p.ncols -1; ll >=0; ll--)
dest[knx + ll] = p.grid2D(p.ncols -1-ll+iknx);
knx-=p.ncols;
iknx+=p.ncols;
}
}
else if(cellOrder == "lrul" ){
unsigned int knx = 6;
unsigned int iknx = 0;
for (unsigned int kk = 0; kk < p.nrows; kk++){
for (unsigned int ll=p.ncols -1; ll >=0; ll--)
size_t knx = 6;
size_t iknx = 0;
for (size_t kk = 0; kk < p.nrows; kk++){
for (size_t ll=p.ncols -1; ll >=0; ll--)
dest[knx+ ll] = p.grid2D(p.ncols -1-ll+iknx);
knx+=p.ncols;
iknx+=p.ncols;
}
}
else if(cellOrder == "ullr"){
unsigned int knx = (p.nrows-1)*p.ncols+6;
unsigned int iknx = 0;
for (unsigned int kk = (signed)p.nrows-1; kk >=0; kk--){
for (unsigned int ll=0; ll < p.ncols; ll++)
size_t knx = (p.nrows-1)*p.ncols+6;
size_t iknx = 0;
for (size_t kk = (signed)p.nrows-1; kk >=0; kk--){
for (size_t ll=0; ll < p.ncols; ll++)
dest[knx + ll] = p.grid2D(ll+iknx);
knx-=p.ncols;
iknx+=p.ncols;
}
}
else{
unsigned int knx = 6;
unsigned int iknx = 0;
for (unsigned int kk = 0; kk < p.nrows; kk++){
for (unsigned int ll=0; ll < p.ncols; ll++)
size_t knx = 6;
size_t iknx = 0;
for (size_t kk = 0; kk < p.nrows; kk++){
for (size_t ll=0; ll < p.ncols; ll++)
dest[knx + ll] = p.grid2D(ll+iknx);
knx+=p.ncols;
iknx+=p.ncols;
......
......@@ -392,18 +392,18 @@ void marshal_DOUBLE2D(POPBuffer &buf, DOUBLE2D &data, int maxsize, int flag, POP
(void)maxsize;
(void)*temp;
if (flag & FLAG_MARSHAL) {
unsigned int nx,ny;
size_t nx,ny;
data.size(nx,ny);
buf.Pack(&nx,1);
buf.Pack(&ny,1);
bool keep_nodata = data.getKeepNodata();
buf.Pack(&keep_nodata,1);
if (nx>0 && ny>0) {
const unsigned int nxy=nx*ny;
const size_t nxy=nx*ny;
buf.Pack(&data(0,0), nxy);
}
} else {
unsigned int nx,ny;
size_t nx,ny;
buf.UnPack(&nx,1);
buf.UnPack(&ny,1);
bool keep_nodata;
......@@ -411,7 +411,7 @@ void marshal_DOUBLE2D(POPBuffer &buf, DOUBLE2D &data, int maxsize, int flag, POP
data.setKeepNodata(keep_nodata);
if (nx>0 && ny>0) {
data.resize(nx,ny);
const unsigned int nxy=nx*ny;
const size_t nxy=nx*ny;
buf.UnPack(&data(0,0),nxy);
} else
data.clear();
......@@ -423,7 +423,7 @@ void marshal_DOUBLE3D(POPBuffer &buf, DOUBLE3D &data, int maxsize, int flag, POP
(void)maxsize;
(void)*temp;
if (flag & FLAG_MARSHAL) {
unsigned int nx,ny,nz;
size_t nx,ny,nz;
data.size(nx,ny,nz);
buf.Pack(&nx,1);
buf.Pack(&ny,1);
......@@ -431,11 +431,11 @@ void marshal_DOUBLE3D(POPBuffer &buf, DOUBLE3D &data, int maxsize, int flag, POP
bool keep_nodata = data.getKeepNodata();
buf.Pack(&keep_nodata,1);
if (nx>0 && ny>0 && nz>0) {
const unsigned int nxyz=nx*ny*nz;
const size_t nxyz=nx*ny*nz;
buf.Pack(&data(0,0,0),nxyz);
}
} else {
unsigned int nx,ny,nz;
size_t nx,ny,nz;
buf.UnPack(&nx,1);
buf.UnPack(&ny,1);
buf.UnPack(&nz,1);
......@@ -444,7 +444,7 @@ void marshal_DOUBLE3D(POPBuffer &buf, DOUBLE3D &data, int maxsize, int flag, POP
data.setKeepNodata(keep_nodata);
if (nx>0 && ny>0 && nz>0) {
data.resize(nx,ny,nz);
const unsigned int nxyz=nx*ny*nz;
const size_t nxyz=nx*ny*nz;
buf.UnPack(&data(0,0,0),nxyz);
} else
data.clear();
......@@ -456,18 +456,18 @@ void marshal_INT2D(POPBuffer &buf, INT2D &data, int maxsize, int flag, POPMemspo
(void)maxsize;
(void)*temp;
if (flag & FLAG_MARSHAL) {
unsigned int nx, ny;
size_t nx, ny;
data.size(nx,ny);
buf.Pack(&nx,1);
buf.Pack(&ny,1);
bool keep_nodata = data.getKeepNodata();
buf.Pack(&keep_nodata,1);
if (nx>0 && ny>0) {
const unsigned int nxy=nx*ny;
const size_t nxy=nx*ny;
buf.Pack(&data(0,0), nxy);
}
} else {
unsigned int nx,ny;
size_t nx,ny;
buf.UnPack(&nx,1);
buf.UnPack(&ny,1);
bool keep_nodata;
......@@ -475,7 +475,7 @@ void marshal_INT2D(POPBuffer &buf, INT2D &data, int maxsize, int flag, POPMemspo
data.setKeepNodata(keep_nodata);
if (nx>0 && ny>0) {
data.resize(nx,ny);
const unsigned int nxy=nx*ny;
const size_t nxy=nx*ny;
buf.UnPack(&data(0,0), nxy);
} else
data.clear();
......@@ -487,18 +487,18 @@ void marshal_CHAR2D(POPBuffer &buf, CHAR2D &data, int maxsize, int flag, POPMems
(void)maxsize;
(void)*temp;
if (flag & FLAG_MARSHAL) {
unsigned int nx, ny;
size_t nx, ny;
data.size(nx,ny);
buf.Pack(&nx,1);
buf.Pack(&ny,1);
bool keep_nodata = data.getKeepNodata();
buf.Pack(&keep_nodata,1);
if (nx>0 && ny>0) {
const unsigned int nxy=nx*ny;
const size_t nxy=nx*ny;
buf.Pack(&data(0,0), nxy);
}
} else {
unsigned int nx,ny;
size_t nx,ny;
buf.UnPack(&nx,1);
buf.UnPack(&ny,1);
bool keep_nodata;
......@@ -506,7 +506,7 @@ void marshal_CHAR2D(POPBuffer &buf, CHAR2D &data, int maxsize, int flag, POPMems
data.setKeepNodata(keep_nodata);
if (nx>0 && ny>0) {
data.resize(nx,ny);
const unsigned int nxy=nx*ny;
const size_t nxy=nx*ny;
buf.UnPack(&data(0,0), nxy);
} else
data.clear();
......
......@@ -325,7 +325,7 @@ const std::string SunObject::toString() const
os << "Lat/Long/Alt\t" << std::setw(7) << latitude << "° " << std::setw(7) << longitude << "° " << std::setprecision(0) << std::setw(4) << altitude << "\n";
os << "Elev. thresh.\t" << std::setprecision(1) << elevation_threshold << \n";
const unsigned int colw=10;
const int colw=10;
os << std::setw(colw) << "" << std::fixed << std::setw(colw) << std::setprecision(1) << "toa";
os << std::fixed << std::setw(colw) << std::setprecision(1) << "direct";
os << std::fixed << std::setw(colw) << std::setprecision(1) << "diffuse";
......
......@@ -70,7 +70,7 @@ std::vector<double> Interpol1D::quantiles(const std::vector<double>& X, const st
else if(q>=1.) vecResults[ii] = vecTemp.back();
else {
const double pos = static_cast<double>(vecSize - 1) * q;
const unsigned int ind = static_cast<unsigned int>(pos);
const size_t ind = static_cast<size_t>(pos);
const double delta = pos - static_cast<double>(ind);
const double i1 = vecTemp[ind];
const double i2 = vecTemp[ind+1];
......@@ -505,8 +505,8 @@ int Interpol1D::NoisyLinRegression(const std::vector<double>& in_X, const std::v
const double r_thres = 0.7;
const size_t nb_pts = in_X.size();
//we want at least 4 points AND 85% of the initial data set kept in the regression
const unsigned int min_dataset = (unsigned int)Optim::floor( 0.85*(double)nb_pts );
const unsigned int min_pts = (min_dataset>4)? min_dataset : 4;
const size_t min_dataset = (size_t)Optim::floor( 0.85*(double)nb_pts );
const size_t min_pts = (min_dataset>4)? min_dataset : 4;
double a=A,b,r; //a needs to be initiallized to A in case of fixed_rate
LinRegression(in_X, in_Y, A, B, R, mesg, fixed_rate);
......
......@@ -145,8 +145,8 @@ void Interpol2D::stdPressure(const DEMObject& dem, Grid2DObject& grid) {
grid.set(dem.ncols, dem.nrows, dem.cellsize, dem.llcorner);
//provide each point with an altitude dependant pressure... it is worth what it is...
for (unsigned int j=0; j<grid.nrows; j++) {
for (unsigned int i=0; i<grid.ncols; i++) {
for (size_t j=0; j<grid.nrows; j++) {
for (size_t i=0; i<grid.ncols; i++) {
const double& cell_altitude=dem.grid2D(i,j);
if (cell_altitude!=IOUtils::nodata) {
grid.grid2D(i,j) = Atmosphere::stdAirPressure(cell_altitude);
......@@ -169,8 +169,8 @@ void Interpol2D::constant(const double& value, const DEMObject& dem, Grid2DObjec
grid.set(dem.ncols, dem.nrows, dem.cellsize, dem.llcorner);
//fills a data table with constant values
for (unsigned int j=0; j<grid.nrows; j++) {
for (unsigned int i=0; i<grid.ncols; i++) {
for (size_t j=0; j<grid.nrows; j++) {
for (size_t i=0; i<grid.ncols; i++) {
if (dem.grid2D(i,j)!=IOUtils::nodata) {
grid.grid2D(i,j) = value;
} else {
......@@ -218,8 +218,8 @@ double Interpol2D::IDWCore(const double& x, const double& y, const std::vector<d
grid.set(dem.ncols, dem.nrows, dem.cellsize, dem.llcorner);
//run algorithm
for (unsigned int j=0; j<grid.nrows; j++) {
for (unsigned int i=0; i<grid.ncols; i++) {
for (size_t j=0; j<grid.nrows; j++) {
for (size_t i=0; i<grid.ncols; i++) {
//LL_IDW_pixel returns nodata when appropriate
double r;
const double value = LLIDW_pixel(i,j,vecData_in, vecStations_in, dem, nrOfNeighbors, r); //TODO: precompute in vectors
......@@ -237,7 +237,7 @@ double Interpol2D::IDWCore(const double& x, const double& y, const std::vector<d
}
//calculate a local pixel for LocalLapseIDW
double Interpol2D::LLIDW_pixel(const unsigned int& i, const unsigned int& j,
double Interpol2D::LLIDW_pixel(const size_t& i, const size_t& j,
const std::vector<double>& vecData_in, const std::vector<StationData>& vecStations_in,
const DEMObject& dem, const size_t& nrOfNeighbors, double& r2)
{
......@@ -315,8 +315,8 @@ void Interpol2D::IDW(const std::vector<double>& vecData_in, const std::vector<St
const double xllcorner = dem.llcorner.getEasting();
const double yllcorner = dem.llcorner.getNorthing();
const double cellsize = dem.cellsize;
for (unsigned int j=0; j<grid.nrows; j++) {
for (unsigned int i=0; i<grid.ncols; i++) {
for (size_t j=0; j<grid.nrows; j++) {
for (size_t i=0; i<grid.ncols; i++) {
if (dem.grid2D(i,j)!=IOUtils::nodata) {
grid.grid2D(i,j) = IDWCore((xllcorner+i*cellsize), (yllcorner+j*cellsize),
vecData_in, vecEastings, vecNorthings);
......@@ -370,8 +370,8 @@ void Interpol2D::SimpleDEMWindInterpolate(const DEMObject& i_dem, Grid2DObject&
if(dem_range_slope==0.) dem_range_slope = 1.; //to avoid division by zero below
if(dem_range_curvature==0.) dem_range_curvature = 1.; //to avoid division by zero below
for (unsigned int j=0;j<VW.nrows;j++) {
for (unsigned int i=0;i<VW.ncols;i++){
for (size_t j=0;j<VW.nrows;j++) {
for (size_t i=0;i<VW.ncols;i++){
speed = VW.grid2D(i,j);
if(speed==0.) continue; //we can not apply any correction factor!
dir = DW.grid2D(i,j);
......@@ -442,8 +442,8 @@ void Interpol2D::PrecipSnow(const DEMObject& dem, const Grid2DObject& ta, Grid2D
}
const double dem_max_curvature=dem.max_curvature, dem_range_curvature=(dem.max_curvature-dem.min_curvature);
for (unsigned int j=0;j<grid.nrows;j++) {
for (unsigned int i=0;i<grid.ncols;i++) {
for (size_t j=0;j<grid.nrows;j++) {
for (size_t i=0;i<grid.ncols;i++) {
// Get input data
const double slope = dem.slope(i, j);
const double curvature = dem.curvature(i, j);
......
......@@ -79,7 +79,7 @@ class Interpol2D {
static double IDWCore(const double& x, const double& y,
const std::vector<double>& vecData_in,
const std::vector<double>& vecEastings, const std::vector<double>& vecNorthings);
static double LLIDW_pixel(const unsigned int& i, const unsigned int& j,
static double LLIDW_pixel(const size_t& i, const size_t& j,
const std::vector<double>& vecData_in,
const std::vector<StationData>& vecStations_in,
const DEMObject& dem, const size_t& nrOfNeighbors, double& r2);
......
......@@ -164,7 +164,7 @@ void ARCIO::cleanup() throw()
void ARCIO::read2DGrid_internal(Grid2DObject& grid_out, const std::string& full_name)
{
int i_ncols, i_nrows;
unsigned int ncols, nrows;
size_t ncols, nrows;
double xllcorner, yllcorner, cellsize, plugin_nodata;
double tmp;
std::string line;
......@@ -209,8 +209,8 @@ void ARCIO::read2DGrid_internal(Grid2DObject& grid_out, const std::string& full_
if((i_ncols<0) || (i_nrows<0)) {
throw IOException("Number of rows or columns in 2D Grid read as \"nodata\", in file: " + full_name, AT);
}
ncols = (unsigned int)i_ncols;
nrows = (unsigned int)i_nrows;
ncols = (size_t)i_ncols;
nrows = (size_t)i_nrows;
//compute/check WGS coordinates (considered as the true reference) according to the projection as defined in cfg
Coords location(coordin, coordinparam);
......@@ -219,9 +219,9 @@ void ARCIO::read2DGrid_internal(Grid2DObject& grid_out, const std::string& full_
//Initialize the 2D grid
grid_out.set(ncols, nrows, cellsize, location);
unsigned int nr_empty=0;
size_t nr_empty=0;
//Read one line after the other and parse values into Grid2DObject
for (unsigned int kk=nrows-1; (kk < nrows); kk--) {
for (size_t kk=nrows-1; (kk < nrows); kk--) {
getline(fin, line, eoln);
if(line.empty()) { //so we can tolerate empty lines
kk++; //to keep the same kk at the next iteration
......@@ -232,7 +232,7 @@ void ARCIO::read2DGrid_internal(Grid2DObject& grid_out, const std::string& full_
iss.setf(std::ios::fixed);
iss.precision(std::numeric_limits<double>::digits10);
for (unsigned int ll=0; ll < ncols; ll++){
for (size_t