Как провести геометрическую калибровку в Java с помощью GDAL? - PullRequest
0 голосов
/ 22 марта 2019

Я хочу создать некоторый gcp и использовать калибровку TPS, чтобы выполнить геометрическую калибровку. Я пробую метод: 'AutoCreateWarpedVRT', но не смог, когда я пытаюсь вывести данные в свой драйвер в формате tif. Я успешно выполнил калибровку TPS в Python, но я хочу сделать калибровку в Java. JAVA:

        gdal.AllRegister();
        String cal_out = "E:\\FY\\GUI_design\\process_java_file\\test_cal.tif";
        String raw = "E:\\FY\\GUI_design\\process_java_file\\4000.tif";
        Dataset cal_dataset = gdal.Open(cal_out,gdalconstConstants.GA_ReadOnly);
        Dataset ref_dataset = gdal.Open(raw,gdalconstConstants.GA_ReadOnly);
        int Xsize = ref_dataset.GetRasterXSize();
        int Ysize = ref_dataset.GetRasterYSize();
        float[][] lon = ReadBand.readband(ref_dataset,1);
        float[][] lat = ReadBand.readband(ref_dataset,2);
        ArrayList<GCP> gcps = new ArrayList<>();
        for (int i=0;i<Xsize;i++){
            for (int j=0;j<Ysize;j++){
                if (lon[i][j]<360){
                    gcps.add(new GCP((double) lat[i][j],(double) lon[i][j],(double)j,(double)i ));
                }
            }
        }
        GCP[] gcps2 = new GCP[gcps.size()/10000 + 1];

        int j =0;
        for (int i=0;i<gcps.size();i=i+10000){
            gcps2[j] = gcps.get(i);
            j=j+1;
        }

        Dataset ds = gdal.Open(cal_out,gdalconstConstants.GA_Update);
        Vector<String> options = new Vector<>();

        SpatialReference sr = new SpatialReference();
        sr.SetWellKnownGeogCS("WGS84");
        ds.SetGCPs(gcps2,sr.ExportToWkt());
        WarpOptions warpOptions = new WarpOptions(options)
        Dataset warped = gdal.AutoCreateWarpedVRT(ds,null,sr.ExportToWkt(),gdalconstConstants.GRA_NearestNeighbour,0.5);
        System.out.println(warped.GetRasterXSize()+"   "+warped.GetRasterYSize());
        Driver wDriver =  warped.GetDriver();
        Dataset dstdataset = wDriver.Create("E:\\FY\\GUI_design\\test.tif",warped.GetRasterXSize(),warped.GetRasterYSize(),warped.getRasterCount(),gdalconstConstants.GDT_Float32);
        System.out.println("lat,lon:"+warped.GetGeoTransform()+" \n "+warped.GetProjection());
        Band testband = warped.GetRasterBand(1);
        float[] buf = new float[warped.GetRasterXSize()*warped.getRasterYSize()];
        testband.ReadRaster(0,0,warped.GetRasterXSize(),warped.getRasterYSize(),buf);
        int [] bandlist = new int[]{1};
        dstdataset.WriteRaster(0,0,warped.GetRasterXSize(),warped.getRasterYSize(),warped.GetRasterXSize(),warped.getRasterYSize(),gdalconstConstants.GDT_Float32,buf,bandlist);

Python:

gcps = []
for _i in range(0, fy1.mask.shape[0]):
    for _j in range(0, 2748):
        if fy1.mask[_i, _j] == 1 and _refer_data[_i, _j, 1] < 360:
            gcps.append(gdal.GCP(_refer_data[_i, _j, 1], _refer_data[_i, _j, 0], 0, _j, _i))
gcps_list = gcps[::10000] 

print("GCP count:", len(gcps_list))

ds: gdal.Dataset = gdal.Open(fn_output, gdal.GA_Update) ##wtf is this
sr = osr.SpatialReference()
sr.SetWellKnownGeogCS('WGS84')
ds.SetGCPs(gcps_list, sr.ExportToWkt())

print("Translating...")
dstDS = gdal.Warp(r'fytest_tps3.tif', ds, format='GTiff', tps=True, xRes=0.05, yRes=0.05, dstNodata=65535,
                  srcNodata=65535, resampleAlg=gdal.GRIORA_NearestNeighbour)

print("Done.")
...