C#集成Aunspline气象数据插值软件
关注公众号,分享GIS知识、ArcGIS教程、SCI论文与科研日常等
C#集成Aunspline气象数据插值软件
//*******************这一部分只是用来读协变量栅格数据的属性值********************
IWorkspaceFactory workspaceFactory = new RasterWorkspaceFactory();
IWorkspace workspace;
workspace = workspaceFactory.OpenFromFile(@"D:\interpolation\B", 0);
IRasterWorkspace rastWork = (IRasterWorkspace)workspace;
IRasterDataset rasterdataset = rastWork.OpenRasterDataset("NeiM_dem_49n.tif");
IRaster raster = rasterdataset.CreateDefaultRaster();
IRasterProps rasterprops = (IRasterProps)raster;
int width = rasterprops.Width;
int height = rasterprops.Height;
double xcell = rasterprops.MeanCellSize().X;
double ycell = rasterprops.MeanCellSize().Y;
double xmin = rasterprops.Extent.XMin;
double ymin = rasterprops.Extent.YMin;
//*******************这一部分只是用来读协变量栅格数据的属性值********************
//******************也可以在arcgis种打开直接读取数据*******************************
//******************接下来批量根据储存好的dat文件编写对应的cmd文件**************
StreamWriter sw;
string[] findfiles = Directory.GetFiles(@"D:\GISer\GuoDa\Temp_dat", "*.dat");
for (int i = 0; i < findfiles.Count(); i++)
{
string name = System.IO.Path.GetFileNameWithoutExtension(findfiles[i]);
sw = new StreamWriter("D:\\GISer\\GuoDa\\Temp_dat\\" + name + ".cmd", false, Encoding.Default);
// StreamWriter函数编写cmd文件
//文件名
sw.WriteLine(name);
sw.WriteLine(5);
//摄氏度。一般有9个类型,分别为:0-未定义;1-米;2-英尺;3-千米;4-英里;5-度;6-弧度;7-毫米;8-兆焦耳;在此做气温所以用5,降水时用8。
sw.WriteLine(2);//两个自变量
sw.WriteLine(1);//一个协变量
sw.WriteLine(0);//曲面样条变量个数,每个曲面采用不同的样条变量时使用
sw.WriteLine(0);//曲面协变量个数,每个曲面采用不同的协变量时使用
sw.WriteLine((xmin - xcell) + " " + (xmin + xcell * width + xcell) + " " + " " + 0 + " " + 1);
//x(一般为经度)的最小范围与最大范围(一般比dem里的最左端坐标读值与最右端坐标读值要大一点);是否转换数据,一般不用转换,为0;所用单位(在此要注意该文件里的范围要比dem数据范围大,一般可以左右多加一列栅格计算范围,且范围要满足公式:Xmax=Xmin+cellsize*rows)。此处读取协变量文件根据属性值直接计算比较准确方便。
sw.WriteLine((ymin - ycell) + " " + (ymin + ycell * height + ycell) + " " + " " + 0 + " " + 1);
//y(一般为纬度)的最小范围与最大范围(一般比dem里的最下端坐标读值与最上端坐标读值大一点);是否转换数据,一般不用转换,为0;所用单位(在此要注意该文件里的范围要比dem数据范围大,一般可以上下多加一行栅格计算范围,且范围要满足公式:Ymax=Ymin+cellsize*cols)
sw.WriteLine(-200 + " " + 8000 + " " + 1 + " " + 1);
//高程范围:最低高程 最高高程;是否需要转换,在此高程需要转换,因此选1,一般定义为0-不转换,1-转换;转换方式,在此有8个转换方式,从简单数学计算到复杂函数公式都由包括,在此不一一描述,本例中使用了1,其转化方式为x/a,即除法公式,x为此的未知数高程,a为自己设置的常数。
sw.WriteLine(1000);
//转换时常数a的数值,在此为1000,表示我对原数据缩小1000倍,即将单位为米的高程数据转化为单位为千米
sw.WriteLine(0);//因变量转换。具有的转换方法有:0-不转换;1-对数转换;2平方根转换;5-类似二值化转化,正数设为1 ,负数忽略。
sw.WriteLine(3);//样条次数
sw.WriteLine(1);//输出表面个数,与下面所给的插值数据个数对应。本例只做一个变量一个时间的插值。
sw.WriteLine(0);//相对误差方差个数,一般为0表示采用同意权重。
sw.WriteLine(1);//优化参数指标
sw.WriteLine(1);//平滑参数选择方式。其中:1-GCV;2-MSE;3-固定值;4-GML
sw.WriteLine(name + ".dat");//输入的插值数据文件
sw.WriteLine(250);//插值样本总数,比实际站点稍多且为整数
sw.WriteLine(5);//站标字符数
sw.WriteLine("(a5,2f14.4,f7.1,f9.4)");//输入文件数据格式定义
sw.WriteLine(name + ".res");
sw.WriteLine(name + ".opt");
sw.WriteLine(name + ".sur");
sw.WriteLine(name + ".lis");
sw.WriteLine(name + ".cov");
sw.WriteLine();//空几行,这个比较重要,负责可能出错。
sw.WriteLine();
sw.WriteLine();
sw.WriteLine();
sw.WriteLine();
sw.WriteLine();
sw.WriteLine();
sw.WriteLine();
sw.Close();
sw = new StreamWriter("D:\\GISer\\GuoDa\\Temp_dat\\" + name + "grd.cmd", false, Encoding.Default); StreamWriter函数编写lapgrd.cmd文件
sw.WriteLine(name + ".sur");//曲面文件
sw.WriteLine(1);//输出表面个数
sw.WriteLine(1);//计算类型。0-只计算统计信息;1-拟合曲面
sw.WriteLine(name + ".cov");//误差协方差文件
sw.WriteLine(2);//误差类型,1-模型标准误差;2-预测标准误差;3-95%模型置信区间;4-95%预测置信区间
sw.WriteLine();//最大标准误差填写,可以不填,但要空出来
sw.WriteLine(1);
sw.WriteLine(1);
sw.WriteLine(xmin + " " + (xmin + width * xcell) + " " + xcell);//X范围,出图时使用的dem范围,如果不使用dem插值出图,则按照.sur文件填写;1000为分辨率
sw.WriteLine(2);
sw.WriteLine(ymin + " " + (ymin + height * ycell) + " " + ycell);//Y范围,出图时使用的dem范围,如果不使用dem插值出图,则按照.sur文件填写;1000为分辨率
sw.WriteLine(0);//掩模方式。0-不掩模
sw.WriteLine(2);//独立协变量数据格式,1-普通栅格;2-Arc/Info格式栅格;3-Idrisi格式栅格
sw.WriteLine("neim_dem1k.txt");//协变量文件,dem数据,用arcgis转换为ASCII格式
sw.WriteLine(2);//输出变量格式
sw.WriteLine(-9999);//指示空格数据
sw.WriteLine(name+".grd");//插值结果,有几个写几行,本例只有一个。
sw.WriteLine("(100f10.2)");//输出文件数据格式。
sw.WriteLine(2);//输出数据格式,下同。
sw.WriteLine(-9999);
sw.WriteLine(name + "_cov.grd");
sw.WriteLine("(100f10.2)");
sw.WriteLine();
sw.WriteLine();
sw.WriteLine();
sw.WriteLine();
sw.WriteLine();
sw.WriteLine();
sw.WriteLine();
sw.WriteLine();
sw.Close();
}
//****************集成运行cmd的函数*********************
static bool RunCmd(string command,string path, string name)
{
bool result = false;
try
{
Process p = new Process();
p.StartInfo.FileName = "cmd.exe";
p.StartInfo.UseShellExecute = false;//是否使用操作系统shell启动
p.StartInfo.RedirectStandardInput = true;//接受来自调用程序的输入信息
p.StartInfo.RedirectStandardOutput = true;//由调用程序获取输出信息
p.StartInfo.RedirectStandardError = true;//重定向标准错误输出
p.StartInfo.CreateNoWindow = true;//不显示程序窗口
p.Start();//启动程序
//string strCMD2="G:"+
//string strCMD = "\"" + @"C:\Program Files (x86)\Microsoft SDKs\Windows\v8.0A\bin\NETFX 4.0 Tools\SvcUtil.exe" + "\" " + ""
// + " /r:" + "\"" + @"C:\Program Files (x86)\Reference Assemblies\Microsoft\Framework\.NETFramework\v4.0\System.Data.dll" + "\"" + " /syncOnly";
向cmd窗口发送输入信息
//p.StandardInput.WriteLine(strCMD + "&exit");
p.StandardInput.WriteLine("G:");
p.StandardInput.WriteLine("cd " + path);
p.StandardInput.WriteLine(command+" <" + name +".cmd> " + name + ".log" + "&exit");
p.StandardInput.AutoFlush = true;
//获取cmd窗口的输出信息
string output = p.StandardOutput.ReadToEnd();
//等待程序执行完退出进程
p.WaitForExit();
p.Close();
result=true;
}
catch (Exception ex)
{
}
return result;
}
//**************************批量运行cmd(调用cmd函数)*********************
string path = @"D:\GISer\GuoDa\Temp_dat";
string[] findfiles = Directory.GetFiles(path, "*.cmd");
string[] findgrdfiles = Directory.GetFiles(path, "*grd.cmd");
List<string> find_nongrd_files = new List<string>(); ;
for (int i = 0; i < findfiles.Count(); i++)
{
bool b=true;
for (int j = 0; j < findgrdfiles.Count();j++)
{
if (findfiles[i].Equals(findgrdfiles[j]))
{
b = false;
break;
}
}
if (b == true)
find_nongrd_files.Add(findfiles[i]);
}
//string command1 = "splina";
string command2 = "lapgrd";
for (int i = 0; i < find_nongrd_files.Count(); i++)
{
//string command = "";
//command = command2;
bool bo;
//switch (command)
//{
//bo = RunCmd(command1, path, System.IO.Path.GetFileNameWithoutExtension(find_nongrd_files[i]));
//break;
bo = RunCmd(command2, path, System.IO.Path.GetFileNameWithoutExtension(findgrdfiles[i]));
//break;
//}
//bool bo = RunCmd(command2,path, System.IO.Path.GetFileNameWithoutExtension(find_nongrd_files[i]));
}