/////////////////////////////////////////////////////////////////////////////// 
// Copyright (C) 2002-2025, Open Design Alliance (the "Alliance"). 
// All rights reserved. 
// 
// This software and its documentation and related materials are owned by 
// the Alliance. The software may only be incorporated into application 
// programs owned by members of the Alliance, subject to a signed 
// Membership Agreement and Supplemental Software License Agreement with the
// Alliance. The structure and organization of this software are the valuable  
// trade secrets of the Alliance and its suppliers. The software is also 
// protected by copyright law and international treaty provisions. Application  
// programs incorporating this software must include the following statement 
// with their copyright notices:
//   
//   This application incorporates Open Design Alliance software pursuant to a license 
//   agreement with Open Design Alliance.
//   Open Design Alliance Copyright (C) 2002-2025 by Open Design Alliance. 
//   All rights reserved.
//
// By use of this software, its documentation or related materials, you 
// acknowledge and accept the above terms.
///////////////////////////////////////////////////////////////////////////////

#include "OdaCommon.h"
#include "PdfImportGeo.h"
#include "DbViewportTable.h"
#include "DbDatabase.h"
#include "DbDictionary.h"
#include "DbDictionaryVar.h"
#include "OdGeoMapType.h"


const double BasePdfImportGeoDataTransformer::cPdfDPI = 72.;

std::vector<OdPdfImportPathData> BasePdfImportGeoDataTransformer::Transform(const OdPdfImportPathData* path_data, const unsigned long path_point_count)
{
  std::vector<OdPdfImportPathData> ret(path_point_count);
  for (size_t i = 0; i < path_point_count; ++i)
  {
    ret[i].m_Point = Transform(&path_data[i].m_Point);
    ret[i].m_CloseFigure = path_data[i].m_CloseFigure;
    ret[i].m_Type = path_data[i].m_Type;
  }
  return ret;
}

void PdfImportNewGeoDataTransformer::Transform(const double& x, const double& y, double& new_x, double& new_y)
{
  applyHomography(x * cPdfDPI / m_Width, y * cPdfDPI / m_Height, new_x, new_y);
}

PdfImportNewGeoDataTransformer::PdfImportNewGeoDataTransformer(OdDbDatabasePtr& pDb, const std::vector<FPDF_PageGeoMeasureDictInfo>& geo_mes,
  const float width, const float height)
  :BasePdfImportGeoDataTransformer(geo_mes, width, height)
{
  m_pGeoData = OdDbGeoData::createObject();
  OdDbObjectId active_space = pDb->getActiveLayoutBTRId();
  m_pGeoData->setBlockTableRecordId(active_space);
  OdDbObjectId objId;
  m_pGeoData->postToDb(objId); 
  
  m_pGeoData->setCoordinateSystem(geo_mes[0].WKT.c_str());

  OdGePoint3d ptBuf;
  OdGePoint2dArray LPTS_array, GPTS_array;
  LPTS_array.resize((int)geo_mes[0].GPTS.size() / 2);
  GPTS_array.resize((int)geo_mes[0].GPTS.size() / 2);
  for (OdUInt32 i = 0; i < geo_mes[0].GPTS.size(); i += 2)
  {
    LPTS_array[i / 2] = OdGePoint2d(geo_mes[0].LPTS[i], geo_mes[0].LPTS[i + 1]);
    m_pGeoData->transformFromLonLatAlt(OdGePoint3d(geo_mes[0].GPTS[i + 1], geo_mes[0].GPTS[i], 0.), ptBuf);
    GPTS_array[i / 2] = OdGePoint2d(ptBuf.x, ptBuf.y);
  }

  computeHomography(LPTS_array, GPTS_array);

  {
    OdDbViewportTablePtr pVt = pDb->getViewportTableId().openObject();
    OdDbObjectId viewportId = pVt->getActiveViewportId();

    OdDbObjectPtr pVp = viewportId.openObject(OdDb::kForWrite);
    OdDbObjectId dictId = pVp->extensionDictionary();
    if (dictId.isNull())
    {
      pVp->createExtensionDictionary();
      dictId = pVp->extensionDictionary();
    }

    OdDbDictionaryPtr pExtDictionary = dictId.openObject(OdDb::kForWrite);
    OdDbObjectId idVariableDictionary = pExtDictionary->getAt(L"AcDbVariableDictionary");
    if (idVariableDictionary.isNull())
    {
      pExtDictionary->setAt(L"AcDbVariableDictionary", OdDbDictionary::createObject());
      idVariableDictionary = pExtDictionary->getAt(L"AcDbVariableDictionary");
    }

    OdDbDictionaryPtr pExtVariableDictionary = idVariableDictionary.openObject(OdDb::kForWrite);

    OdDbObjectId idGEOMAPMODE = pExtVariableDictionary->getAt(L"GEOMAPMODE");
    if (idGEOMAPMODE.isNull())
    {
      pExtVariableDictionary->setAt(L"GEOMAPMODE", OdDbDictionaryVar::createObject());
      idGEOMAPMODE = pExtVariableDictionary->getAt(L"GEOMAPMODE");
    }
    OdDbDictionaryVarPtr pDictionaryVar = pExtVariableDictionary->getAt(L"GEOMAPMODE").openObject(OdDb::kForWrite);
    OdGeoMapType iGeoMapType = kNoMap;
    pDictionaryVar->setValue(OdString().format(L"%d", iGeoMapType));
  }
}

void PdfImportExistingGeoDataTransformer::Transform(const double& x, const double& y, double& new_x, double& new_y)
{
  OdGePoint3d ptBuf;
  applyHomography(x * cPdfDPI / m_Width, y * cPdfDPI / m_Height, new_x, new_y);
  // toLLA
  m_pTransformer->transformPoint(OdGePoint3d(new_x, new_y, 0.), ptBuf);
  // to Existing GEO
  m_pGeoData->transformFromLonLatAlt(ptBuf, ptBuf);
  new_x = ptBuf.x;
  new_y = ptBuf.y;
}

PdfImportExistingGeoDataTransformer::PdfImportExistingGeoDataTransformer(OdDbObjectId& idGeoData, const std::vector<FPDF_PageGeoMeasureDictInfo>& geo_mes,
  const float width, const float height)
  :BasePdfImportGeoDataTransformer(geo_mes, width, height)
{
  if (idGeoData)
  {
    m_pGeoData = idGeoData.openObject();
  }

  OdDbGeoCoordinateSystemTransformerPtr pFromLLTransformer;
  OdDbGeoCoordinateSystemTransformer::create(L"LL84", geo_mes[0].WKT.c_str(), pFromLLTransformer);
  OdGePoint3d ptBuf;
  OdGePoint2dArray LPTS_array, GPTS_array;
  LPTS_array.resize((int)geo_mes[0].GPTS.size() / 2);
  GPTS_array.resize((int)geo_mes[0].GPTS.size() / 2);
  for (OdUInt32 i = 0; i < geo_mes[0].GPTS.size(); i += 2)
  {
    LPTS_array[i / 2] = OdGePoint2d(geo_mes[0].LPTS[i], geo_mes[0].LPTS[i + 1]);
    pFromLLTransformer->transformPoint(OdGePoint3d(geo_mes[0].GPTS[i + 1], geo_mes[0].GPTS[i], 0.), ptBuf);
    GPTS_array[i / 2] = OdGePoint2d(ptBuf.x, ptBuf.y);
  }

  computeHomography(LPTS_array, GPTS_array);

  OdDbGeoCoordinateSystemTransformer::create(geo_mes[0].WKT.c_str(), L"LL84", m_pTransformer);
}

std::unique_ptr<BasePdfImportGeoDataTransformer> createTransformer2PdfImport(OdDbDatabasePtr& pDb, const std::vector<FPDF_PageGeoMeasureDictInfo>& geo_mes, 
  const float width, const float height)
{
  if (0 != geo_mes.size())
  {
    if(geo_mes[0].LPTS.size() < 4 || geo_mes[0].LPTS.size() != geo_mes[0].GPTS.size() )
      return nullptr;

    OdDbObjectId idGeoData;
    oddbGetGeoDataObjId(pDb, idGeoData);
    if (idGeoData)
    {
      return std::unique_ptr<PdfImportExistingGeoDataTransformer>(new PdfImportExistingGeoDataTransformer(idGeoData, geo_mes, width, height));
    }
    else
    {
      return std::unique_ptr<PdfImportNewGeoDataTransformer>(new PdfImportNewGeoDataTransformer(pDb, geo_mes, width, height));
    }
  }
  return nullptr;
}

BasePdfImportGeoDataTransformer::BasePdfImportGeoDataTransformer(const std::vector<FPDF_PageGeoMeasureDictInfo>& geo_mes, const float width, const float height)
  :m_Width(width)
  ,m_Height(std::fabs(height))
{

}

OdPdfImportPoint BasePdfImportGeoDataTransformer::Transform(const OdPdfImportPoint* pt)
{
  OdPdfImportPoint ret;
  Transform(pt->x, pt->y, ret.x, ret.y);
  return ret;
}

void BasePdfImportGeoDataTransformer::computeHomography(const OdGePoint2dArray& src, const OdGePoint2dArray& dst)
{
  OdGePoint2dArray::size_type n = src.size();

  const OdGePoint2dArray::size_type N = 2 * n;
  const OdGePoint2dArray::size_type C = 9;
  std::vector<std::array<double, C>> A(N, { 0. });

  for (OdGePoint2dArray::size_type i = 0; i < n; ++i)
  {
    double x = src[i].x, y = src[i].y;
    double X = dst[i].x, Y = dst[i].y;

    A[2 * i][0] = x;
    A[2 * i][1] = y;
    A[2 * i][2] = 1.0;
    A[2 * i][6] = -x * X;
    A[2 * i][7] = -y * X;
    A[2 * i][8] = -X;

    A[2 * i + 1][3] = x;
    A[2 * i + 1][4] = y;
    A[2 * i + 1][5] = 1.0;
    A[2 * i + 1][6] = -x * Y;
    A[2 * i + 1][7] = -y * Y;
    A[2 * i + 1][8] = -Y;
  }

  double AtA[C][C] = { 0 };
  for (OdGePoint2dArray::size_type i = 0; i < C; ++i) 
  {
    for (OdGePoint2dArray::size_type j = 0; j < C; ++j) 
    {
      for (OdGePoint2dArray::size_type k = 0; k < N; ++k)
        AtA[i][j] += A[k][i] * A[k][j];
    }
  }

  double b[C] = { 0 };
  b[C - 1] = 1.0;

  double M[C][C + 1];
  for (OdGePoint2dArray::size_type i = 0; i < C; ++i) 
  {
    for (OdGePoint2dArray::size_type j = 0; j < C; ++j)
      M[i][j] = AtA[i][j];
    M[i][C] = b[i];
  }

  for (OdGePoint2dArray::size_type i = 0; i < C; ++i) 
  {
    double maxEl = std::fabs(M[i][i]);
    OdGePoint2dArray::size_type maxRow = i;
    for (OdGePoint2dArray::size_type k = i + 1; k < C; ++k) 
    {
      if (std::fabs(M[k][i]) > maxEl) {
        maxEl = std::fabs(M[k][i]);
        maxRow = k;
      }
    }

    for (OdGePoint2dArray::size_type k = i; k <= C; ++k)
      std::swap(M[maxRow][k], M[i][k]);

    for (OdGePoint2dArray::size_type k = i + 1; k < C; ++k) 
    {
      double c = -M[k][i] / M[i][i];
      for (OdGePoint2dArray::size_type j = i; j <= C; ++j) 
      {
        if (i == j)
          M[k][j] = 0;
        else
          M[k][j] += c * M[i][j];
      }
    }
  }

  m_Homography.resize(C);
  for (OdGePoint2dArray::size_type i = C; i > 0; --i)
  {
    const OdGePoint2dArray::size_type idx = i - 1;
    m_Homography[idx] = M[idx][C];
    for (OdGePoint2dArray::size_type j = idx + 1; j < C; ++j)
      m_Homography[idx] -= M[idx][j] * m_Homography[j];
    m_Homography[idx] /= M[idx][idx];
  }

  for (OdGePoint2dArray::size_type i = 0; i < 9; ++i)
    m_Homography[i] /= m_Homography[8];

}

void BasePdfImportGeoDataTransformer::applyHomography(const double& xx, const double& yy, double& new_x, double& new_y)
{
  if (m_Homography.size() < 9)
  {
    new_x = xx;
    new_y = yy;
    return;
  }

  double x = xx;
  double y = yy;
  double denom = m_Homography[6] * x + m_Homography[7] * y + m_Homography[8];
  new_x = (m_Homography[0] * x + m_Homography[1] * y + m_Homography[2]) / denom;
  new_y = (m_Homography[3] * x + m_Homography[4] * y + m_Homography[5]) / denom;
  return;
}