Получить строка/столбец из Ш/Д
-
- Участник
- Сообщения: 93
- Зарегистрирован: 04 дек 2013, 02:14
- Репутация: 1
Получить строка/столбец из Ш/Д
Есть GDALDataset. Как получить координаты в пикселях (причем в вещественных числах) из широты/долготы? Как обратно сделать я знаю, но вот так - нет. Кто-то мне давал ссылку на функцию int GDALInvGeoTransform(double * gt_in, double * gt_out), но это явно не оно.
- Дмитрий Барышников
- Гуру
- Сообщения: 2572
- Зарегистрирован: 17 ноя 2009, 19:17
- Репутация: 261
- Откуда: Москва
Re: Получить строка/столбец из Ш/Д
Прежде чем заводить тему стоит посмотреть, нет ли похожих.
viewtopic.php?f=30&t=16017&p=108576
viewtopic.php?f=30&t=16017&p=108576
-
- Участник
- Сообщения: 93
- Зарегистрирован: 04 дек 2013, 02:14
- Репутация: 1
Re: Получить строка/столбец из Ш/Д
Прежде чем говорить о похожих темах стоит удостовериться, что в них есть решение проблемы.
И да, я смотрел ту тему. И да, там есть такой вопрос. Но нет, там нет решения, а лишь ваше высказывание:
Поэтому вопрос остался.
И да, я смотрел ту тему. И да, там есть такой вопрос. Но нет, там нет решения, а лишь ваше высказывание:
которое ничем помочь конечно же не может.Дмитрий Барышников писал(а):1) Для преобразования из географических координат в пиксельные OGRSpatialReference не подходит. Он вообще не может ничего преобразовывать. Скорее OGRCoordinateTransofrmation - да и то, только для преобразования между различными системами координат.
Поэтому вопрос остался.
-
- Участник
- Сообщения: 93
- Зарегистрирован: 04 дек 2013, 02:14
- Репутация: 1
Re: Получить строка/столбец из Ш/Д
Ну как же черт возьми получить X/Y в пикселях по широте/долготе???? Меня уже несколько недель мучает эта проблема. Отчаянно прошу помощи.... И нету на этом форуме темы с решением!!! Дайте простой код получения нужной информации из GDALDataSetH!!!!!!!!!!!!!!
-
- Участник
- Сообщения: 93
- Зарегистрирован: 04 дек 2013, 02:14
- Репутация: 1
Re: Получить строка/столбец из Ш/Д
up. Вопрос встал остро! хелпми!!!
-
- Гуру
- Сообщения: 4069
- Зарегистрирован: 15 окт 2010, 08:33
- Репутация: 1064
- Ваше звание: программист
- Откуда: Казань
Re: Получить строка/столбец из Ш/Д
код за вас никто писать, конечно, не будет, но схема действий такая (если нет извратов типа повернутого изображения)wowka1319 писал(а):up. Вопрос встал остро! хелпми!!!
1) достаете строку +proj и переводите долготу/широту в координаты проекции снимка
2) достаете bbox (охватывающий прямоугольник) Xmin,Ymin,Xmsx,Ymax, вычисляете dx,dy (я не уверен, но вроде bbox считается по внешним краям пикселей), dx=(Xmax-Xmin)/Nx, Nx x Ny - размеры снимка, и вроде dx, dy там готовые есть
3) вычисляете номер пикселя Xpix=(X-Xmin)/dx, нумерация с 0
4) profit
P.S. я все это делаю в R, там обращение в gdal через обертку, но принцип тот же.
-
- Участник
- Сообщения: 93
- Зарегистрирован: 04 дек 2013, 02:14
- Репутация: 1
Re: Получить строка/столбец из Ш/Д
как? обратный процесс мне пол силу, но вот так - нет.gamm писал(а): 1) достаете строку +proj и переводите долготу/широту в координаты проекции снимка
А дальше вроде бы понятно, но все-таки опишу чтоб вы если что поправили:
2) вычисляем размер пикселя (dx и dy) в условных единицах проекции снимка по X и Y (они разве не равны?)
3) отталкиваясь от Xmin и Ymax (не очень ясно понимаю где их брать, но думаю разберусь) считаю (X-Xmin)/dx и (Y-Ymin)/dy
4)
верно? ну прошу разъяснить первый пункт - именно в нем загвоздка!!!
- Дмитрий Барышников
- Гуру
- Сообщения: 2572
- Зарегистрирован: 17 ноя 2009, 19:17
- Репутация: 261
- Откуда: Москва
Re: Получить строка/столбец из Ш/Д
Из пикселей в координаты
Из координат в пиксели
Если координаты снимка отличаются от координат на входе - их надо перепроецировать в координаты снимка при помощи OGRCoordinateTransformation
Код: Выделить всё
GDALDataset *pDataset = (GDALDataset *)GDALOpen("test.tiff", GA_ReadOnly);
double adfGeoTransform[6] = {0, 0, 0, 0, 0, 0};
pDataset->GetGeoTransform(adfGeoTransform);
double dfRow = 1.5, dfCol = 2.3;
double dfX, dfY;
GDALApplyGeoTransform(adfInvGeoTransform, dfRow, dfCol, &dfX, &dfY);
Код: Выделить всё
GDALDataset *pDataset = (GDALDataset *)GDALOpen("test.tiff", GA_ReadOnly);
double adfGeoTransform[6] = {0, 0, 0, 0, 0, 0};
pDataset->GetGeoTransform(adfGeoTransform);
double adfInvGeoTransform[6] = {0, 0, 0, 0, 0, 0};
GDALInvGeoTransform(adfGeoTransform, adfInvGeoTransform);
double dfRow, dfCol;
double dfX = 37.55, dfY = 55.37;
GDALApplyGeoTransform(adfInvGeoTransform, dfX, dfY, &dfRow, &dfCol);
-
- Участник
- Сообщения: 93
- Зарегистрирован: 04 дек 2013, 02:14
- Репутация: 1
Re: Получить строка/столбец из Ш/Д
увы но ваш код работает неверно:
Во-первых:
Из пикселей в координаты будет так: (это по мануалу GDAL)
не понимаю зачем вы вообще приводили первый код, ведь этого не требовалось
Во-вторых: ваш пример Из координат в пиксели работает НЕВЕРНО.
может вы не поняли задачу? повторю - из координат долгота/широта в градусах получить координаты изображения X/Y в пикселях, т.е. нужен код работающий ровно наоборот коду, приведенному выше (в этом сообщении).
Вопрос остался. Пожалуйста помогите - проблема серьезно встала, не дает работать дальше!!!!!!!
Во-первых:
Из пикселей в координаты будет так: (это по мануалу GDAL)
Код: Выделить всё
double x,y;
double geoTransform[6];
dataset_->GetGeoTransform(geoTransform);
GDALApplyGeoTransform(geoTransform, x, y, &x, &y); //использую те же переменные
const char* projection = dataset_->GetProjectionRef();
OGRSpatialReference* spatialReference = new OGRSpatialReference(projection);
if (spatialReference->IsProjected())
{
OGRSpatialReference* spatialReferenceGeogCS = spatialReference->CloneGeogCS();
OGRCoordinateTransformation* coordinateTransformation = OGRCreateCoordinateTransformation(spatialReference, spatialReferenceGeogCS);
if (coordinateTransformation) coordinateTransformation->Transform(1, &x, &y); //опять использую те же переменные
}
// теперь в x - долгота, а в y - широта в градусах!!!
// проверено - РАБОТАЕТ!!!
Во-вторых: ваш пример Из координат в пиксели работает НЕВЕРНО.
может вы не поняли задачу? повторю - из координат долгота/широта в градусах получить координаты изображения X/Y в пикселях, т.е. нужен код работающий ровно наоборот коду, приведенному выше (в этом сообщении).
Вопрос остался. Пожалуйста помогите - проблема серьезно встала, не дает работать дальше!!!!!!!
-
- Гуру
- Сообщения: 4069
- Зарегистрирован: 15 окт 2010, 08:33
- Репутация: 1064
- Ваше звание: программист
- Откуда: Казань
Re: Получить строка/столбец из Ш/Д
работает верно, вы просто не поняли, что там написано - это код для перевода пикселей в координаты проекции тифа и обратно линейным преобразованием.wowka1319 писал(а):увы но ваш код работает неверно
Если проекция долгота/широта, то получатся искомые пиксели. Если нет, то сначала нужно преобразовать градусы в проекцию тифа. У вас все части есть, в чем вопрос - непонятно.
-
- Участник
- Сообщения: 93
- Зарегистрирован: 04 дек 2013, 02:14
- Репутация: 1
Re: Получить строка/столбец из Ш/Д
Вот именно, что мне нужен был код (сотый раз повторяю) из долготы/широты (НЕ из проекции тифа в метрах, а именно градусов) в пиксели (наоборот не надо - уже давноооо есть). Так никто рабочий код и не привел(((
Но взглянув на обратный код (который я привел постом выше), я по аналогии все-таки добился своего.
Итак, код для тех возможно сталкнется с этой проблемой (думаю Qt типы никого не собьют с толку):
Из пикселей X/Y в градусы Д/Ш (то, что уже приводил постом выше)
dataset_ - это GDALDataset*
И, наконец, искомый код - Из градус Д/Ш в пиксели X/Y:
dataset_ - это GDALDataset*
Но взглянув на обратный код (который я привел постом выше), я по аналогии все-таки добился своего.
Итак, код для тех возможно сталкнется с этой проблемой (думаю Qt типы никого не собьют с толку):
Из пикселей X/Y в градусы Д/Ш (то, что уже приводил постом выше)
Код: Выделить всё
QPointF MMGeoImage::geoCoords(const QPointF &pixelCoords)
{
if (!dataset_) return QPointF();
double point[2];
point[0] = pixelCoords.x();
point[1] = pixelCoords.y();
double geoTransform[6];
dataset_->GetGeoTransform(geoTransform);
GDALApplyGeoTransform(geoTransform, point[0], point[1], &point[0], &point[1]);
const char* projection = dataset_->GetProjectionRef();
OGRSpatialReference* spatialReference = new OGRSpatialReference(projection);
if (spatialReference->IsProjected())
{
OGRSpatialReference* spatialReferenceGeogCS = spatialReference->CloneGeogCS();
OGRCoordinateTransformation* coordinateTransformation = OGRCreateCoordinateTransformation(spatialReference, spatialReferenceGeogCS);
if (coordinateTransformation)
{
coordinateTransformation->Transform(1, &point[0], &point[1]);
return QPointF(point[0], point[1]); //OK
}
else return QPointF();
}
else return QPointF();
}
И, наконец, искомый код - Из градус Д/Ш в пиксели X/Y:
Код: Выделить всё
QPointF MMGeoImage::pixelCoords(const QPointF &geoCoords)
{
if (!dataset_) return QPointF();
double point[2];
point[0] = geoCoords.x();
point[1] = geoCoords.y();
const char* projection = dataset_->GetProjectionRef();
OGRSpatialReference* spatialReference = new OGRSpatialReference(projection);
if (spatialReference->IsProjected())
{
OGRSpatialReference* spatialReferenceGeogCS = spatialReference->CloneGeogCS();
OGRCoordinateTransformation* coordinateTransformation = OGRCreateCoordinateTransformation(spatialReferenceGeogCS, spatialReference);
if (coordinateTransformation)
{
coordinateTransformation->Transform(1, &point[0], &point[1]);
double geoTransform[6];
dataset_->GetGeoTransform(geoTransform);
double invGeoTransform[6];
GDALInvGeoTransform(geoTransform, invGeoTransform);
GDALApplyGeoTransform(invGeoTransform, point[0], point[1], &point[0], &point[1]);
return QPointF(point[0], point[1]); //OK
}
else return QPointF();
}
else return QPointF();
}
Кто сейчас на конференции
Сейчас этот форум просматривают: нет зарегистрированных пользователей и 2 гостя