45using std::normal_distribution;
58bool bCurvaturePairCompareDescending(
const pair<int, double>& prLeft,
const pair<int, double>& prRight)
61 return prLeft.second > prRight.second;
73 for (
unsigned int nCoast = 0; nCoast <
m_VCoast.size(); nCoast++)
76 int const nCoastSize =
m_VCoast[nCoast].nGetCoastlineSize();
79 vector<bool> bVCoastPointDone(nCoastSize,
false);
82 vector<pair<int, double>> prVCurvature;
84 for (
int nCoastPoint = 0; nCoastPoint < nCoastSize; nCoastPoint++)
91 dCurvature =
m_VCoast[nCoast].dGetSmoothCurvature(nCoastPoint);
97 dCurvature =
m_VCoast[nCoast].dGetDetailedCurvature(nCoastPoint);
100 prVCurvature.push_back(make_pair(nCoastPoint, dCurvature));
104 sort(prVCurvature.begin(), prVCurvature.end(), bCurvaturePairCompareDescending);
118 bVCoastPointDone[n] =
true;
120 int const m = nCoastSize - n - 1;
123 bVCoastPointDone[m] =
true;
135 strErr +=
". Check the SWL";
158 m_VCoast[nCoast].InsertProfilesInProfileCoastPointIndex();
179 for (
int nCoastPoint = 0; nCoastPoint < nCoastSize; nCoastPoint++)
181 if (
m_VCoast[nCoast].bIsProfileAtCoastPoint(nCoastPoint))
183 pThisProfile =
m_VCoast[nCoast].pGetProfileAtCoastPoint(nCoastPoint);
186 if (nCoastPoint == 0)
192 pLastProfile = pThisProfile;
201 if (nCoastPoint == nCoastSize - 1)
206 pLastProfile = pThisProfile;
212 m_VCoast[nCoast].CreateProfileDownCoastIndex();
279 int const nCoastSize =
m_VCoast[nCoast].nGetCoastlineSize();
282 for (
int n = nCoastSize - 1; n >= 0;
286 int nStillToSearch = 0;
288 for (
int m = 0; m < nCoastSize; m++)
289 if (! pbVCoastPointDone->at(m))
292 if (nStillToSearch == 0)
297 int const nNormalPoint = prVCurvature->at(n).first;
300 if ((nNormalPoint == 0) || (nNormalPoint == nCoastSize - 1))
305 if (! pbVCoastPointDone->at(nNormalPoint))
308 bool bIntervention =
false;
313 bIntervention =
true;
319 int const nRet =
nCreateProfile(nCoast, nCoastSize, nNormalPoint, nProfile, bIntervention, &PtiThis);
337 pbVCoastPointDone->at(nNormalPoint) =
true;
427 for (
int m = 1; m < dNumToMark; m++)
429 int nTmpPoint = nNormalPoint + m;
431 if (nTmpPoint < nCoastSize)
432 pbVCoastPointDone->at(nTmpPoint) =
true;
434 nTmpPoint = nNormalPoint - m;
437 pbVCoastPointDone->at(nTmpPoint) =
true;
463 int const nXEnd = PtiEnd.
nGetX();
464 int const nYEnd = PtiEnd.
nGetY();
467 if (!
m_pRasterGrid->m_Cell[nXEnd][nYEnd].bIsInContiguousSea())
488 vector<CGeom2DPoint> VNormal;
489 VNormal.push_back(PtStart);
490 VNormal.push_back(PtEnd);
500 m_VCoast[nCoast].AppendProfile(pProfile);
519 LogStream <<
m_ulIter <<
": coast " << nCoast <<
" profile " << nProfile <<
" (nCoastID = " << pProfile->
nGetProfileID() <<
") created at coast point " << nProfileStartPoint <<
" from [" << pPtiStart->
nGetX() <<
"][" << pPtiStart->
nGetY() <<
"] = {" << PtStart.
dGetX() <<
", " << PtStart.
dGetY() <<
"} to [" << PtiEnd.
nGetX() <<
"][" << PtiEnd.
nGetY() <<
"] = {" << PtEnd.
dGetX() <<
", " << PtEnd.
dGetY() <<
"}" << (pProfile->
bIsIntervention() ?
", from intervention" :
"") << endl;
529 int const nCoastSize =
m_VCoast[nCoast].nGetCoastlineSize();
530 int const nHandedness =
m_VCoast[nCoast].nGetSeaHandedness();
532 int nProfileStartEdge;
535 vector<CGeom2DIPoint> VPtiNormalPoints;
540 PtiProfileStart = *
m_VCoast[nCoast].pPtiGetCellMarkedAsCoastline(0);
541 nProfileStartEdge =
m_VCoast[nCoast].nGetStartEdge();
546 PtiProfileStart = *
m_VCoast[nCoast].pPtiGetCellMarkedAsCoastline(nCoastSize - 1);
547 nProfileStartEdge =
m_VCoast[nCoast].nGetEndEdge();
550 VPtiNormalPoints.push_back(PtiProfileStart);
565 int nPos =
static_cast<int>(it -
m_VEdgeCell.begin());
568 for (
int n = 0; n < nProfileLen; n++)
633 int nProfileStartPoint;
639 nProfileStartPoint = 0;
642 pProfile =
new CGeomProfile(nCoast, nProfileStartPoint, nProfile, &PtiProfileStart, &PtiDummy,
false);
649 nProfileStartPoint = nCoastSize - 1;
652 pProfile =
new CGeomProfile(nCoast, nProfileStartPoint, nProfile, &PtiProfileStart, &PtiDummy,
false);
659 for (
unsigned int n = 0; n < VPtiNormalPoints.size(); n++)
661 int const nX = VPtiNormalPoints[n].nGetX();
662 int const nY = VPtiNormalPoints[n].nGetY();
665 m_pRasterGrid->m_Cell[nX][nY].SetCoastAndProfileID(nCoast, nProfile);
676 int const nEndX = VPtiNormalPoints.back().nGetX();
677 int const nEndY = VPtiNormalPoints.back().nGetY();
682 double const dDeepWaterWaveHeight =
m_pRasterGrid->m_Cell[nEndX][nEndY].dGetCellDeepWaterWaveHeight();
683 double const dDeepWaterWaveAngle =
m_pRasterGrid->m_Cell[nEndX][nEndY].dGetCellDeepWaterWaveAngle();
684 double const dDeepWaterWavePeriod =
m_pRasterGrid->m_Cell[nEndX][nEndY].dGetCellDeepWaterWavePeriod();
696 m_VCoast[nCoast].AppendProfile(pProfile);
697 m_VCoast[nCoast].SetProfileAtCoastPoint(nProfileStartPoint, pProfile);
700 LogStream <<
m_ulIter <<
": coast " << nCoast <<
" grid-edge profile " << nProfile <<
" (nCoastID = " << pProfile->
nGetProfileID() <<
") created at coast " << (bCoastStart ?
"start" :
"end") <<
" point " << (bCoastStart ? 0 : nCoastSize - 1) <<
", from [" << PtiProfileStart.
nGetX() <<
"][" << PtiProfileStart.
nGetY() <<
"] = {" <<
dGridCentroidXToExtCRSX(PtiProfileStart.
nGetX()) <<
", " <<
dGridCentroidYToExtCRSY(PtiProfileStart.
nGetY()) <<
"} to [" << VPtiNormalPoints.back().nGetX() <<
"][" << VPtiNormalPoints.back().nGetY() <<
"] = {" <<
dGridCentroidXToExtCRSX(VPtiNormalPoints.back().nGetX()) <<
", " <<
dGridCentroidYToExtCRSY(VPtiNormalPoints.back().nGetY()) <<
"} profile length = " << VPtiNormalPoints.size() << endl;
712 int const AVGSIZE = 21;
720 PtBefore = *
m_VCoast[nCoast].pPtGetCoastlinePointExtCRS(nStartCoastPoint - 1);
721 PtAfter = *
m_VCoast[nCoast].pPtGetCoastlinePointExtCRS(nStartCoastPoint + 1);
727 vector<CGeom2DPoint> PtBeforeToAverage;
729 for (
int n = 1; n <= AVGSIZE; n++)
731 int const nPoint = nStartCoastPoint - n;
736 PtBeforeToAverage.push_back(*
m_VCoast[nCoast].pPtGetCoastlinePointExtCRS(nPoint));
740 vector<CGeom2DPoint> PtAfterToAverage;
742 for (
int n = 1; n <= AVGSIZE; n++)
744 int const nPoint = nStartCoastPoint + n;
746 if (nPoint > nCoastSize - 1)
749 PtAfterToAverage.push_back(*
m_VCoast[nCoast].pPtGetCoastlinePointExtCRS(nPoint));
753 PtBefore =
PtAverage(&PtBeforeToAverage);
758 double const dYDiff = PtAfter.
dGetY() - PtBefore.
dGetY();
759 double const dXDiff = PtAfter.
dGetX() - PtBefore.
dGetX();
761 double dXEnd1 = 0, dXEnd2 = 0, dYEnd1 = 0, dYEnd2 = 0;
766 dXEnd1 = dXEnd2 = pPtStart->
dGetX();
767 dYEnd1 = pPtStart->
dGetY() + dLineLength;
768 dYEnd2 = pPtStart->
dGetY() - dLineLength;
774 dYEnd1 = dYEnd2 = pPtStart->
dGetY();
775 dXEnd1 = pPtStart->
dGetX() + dLineLength;
776 dXEnd2 = pPtStart->
dGetX() - dLineLength;
782 double const dA = dYDiff / dXDiff;
785 double const dAPerp = -1 / dA;
786 double const dBPerp = pPtStart->
dGetY() - (dAPerp * pPtStart->
dGetX());
789 double const dQuadA = 1 + (dAPerp * dAPerp);
790 double const dQuadB = 2 * ((dBPerp * dAPerp) - (dAPerp * pPtStart->
dGetY()) - pPtStart->
dGetX());
791 double const dQuadC = ((pPtStart->
dGetX() * pPtStart->
dGetX()) + (pPtStart->
dGetY() * pPtStart->
dGetY()) + (dBPerp * dBPerp) - (2 * pPtStart->
dGetY() * dBPerp) - (dLineLength * dLineLength));
794 double const dDiscriminant = (dQuadB * dQuadB) - (4 * dQuadA * dQuadC);
796 if (dDiscriminant < 0)
798 LogStream <<
ERR <<
"timestep " <<
m_ulIter <<
": discriminant < 0 when finding profile end point on coastline " << nCoast <<
", from coastline point " << nStartCoastPoint <<
"), ignored" << endl;
802 dXEnd1 = (-dQuadB + sqrt(dDiscriminant)) / (2 * dQuadA);
803 dYEnd1 = (dAPerp * dXEnd1) + dBPerp;
804 dXEnd2 = (-dQuadB - sqrt(dDiscriminant)) / (2 * dQuadA);
805 dYEnd2 = (dAPerp * dXEnd2) + dBPerp;
809 int const nSeaHand =
m_VCoast[nCoast].nGetSeaHandedness();
810 *pPtEnd =
PtChooseEndPoint(nSeaHand, &PtBefore, &PtAfter, dXEnd1, dYEnd1, dXEnd2, dYEnd2);
850 PtChosen.
SetX(dXEnd1);
851 PtChosen.
SetY(dYEnd1);
856 PtChosen.
SetX(dXEnd2);
857 PtChosen.
SetY(dYEnd2);
861 else if (PtAfter->
dGetY() < PtBefore->
dGetY())
866 PtChosen.
SetX(dXEnd1);
867 PtChosen.
SetY(dYEnd1);
872 PtChosen.
SetX(dXEnd2);
873 PtChosen.
SetY(dYEnd2);
885 PtChosen.
SetX(dXEnd1);
886 PtChosen.
SetY(dYEnd1);
891 PtChosen.
SetX(dXEnd2);
892 PtChosen.
SetY(dYEnd2);
901 PtChosen.
SetX(dXEnd1);
902 PtChosen.
SetY(dYEnd1);
907 PtChosen.
SetX(dXEnd2);
908 PtChosen.
SetY(dYEnd2);
922 PtChosen.
SetX(dXEnd1);
923 PtChosen.
SetY(dYEnd1);
928 PtChosen.
SetX(dXEnd2);
929 PtChosen.
SetY(dYEnd2);
933 else if (PtAfter->
dGetY() < PtBefore->
dGetY())
938 PtChosen.
SetX(dXEnd1);
939 PtChosen.
SetY(dYEnd1);
944 PtChosen.
SetX(dXEnd2);
945 PtChosen.
SetY(dYEnd2);
957 PtChosen.
SetX(dXEnd1);
958 PtChosen.
SetY(dYEnd1);
963 PtChosen.
SetX(dXEnd2);
964 PtChosen.
SetY(dYEnd2);
973 PtChosen.
SetX(dXEnd1);
974 PtChosen.
SetY(dYEnd1);
979 PtChosen.
SetX(dXEnd2);
980 PtChosen.
SetY(dYEnd2);
995 int const nCoastLines =
static_cast<int>(
m_VCoast.size());
997 for (
int nCoast = 0; nCoast < nCoastLines; nCoast++)
999 int const nCoastSize =
m_VCoast[nCoast].nGetCoastlineSize();
1002 for (
int nCoastPoint = 0; nCoastPoint < nCoastSize; nCoastPoint++)
1004 if (!
m_VCoast[nCoast].bIsProfileAtCoastPoint(nCoastPoint))
1024 nStartPoint = nCoastPoint + 1;
1026 nStartPoint = nCoastPoint - 1;
1028 for (
int nSecondCoastPoint = nStartPoint; (nDirection ==
DIRECTION_DOWNCOAST) ? nSecondCoastPoint < nCoastSize : nSecondCoastPoint >= 0; (nDirection ==
DIRECTION_DOWNCOAST) ? nSecondCoastPoint++ : nSecondCoastPoint--)
1038 if (
m_VCoast[nCoast].bIsProfileAtCoastPoint(nSecondCoastPoint))
1061 int nProf1LineSeg = 0;
1062 int nProf2LineSeg = 0;
1063 double dIntersectX = 0;
1064 double dIntersectY = 0;
1065 double dAvgEndX = 0;
1066 double dAvgEndY = 0;
1068 if (
bCheckForIntersection(pFirstProfile, pSecondProfile, nProf1LineSeg, nProf2LineSeg, dIntersectX, dIntersectY, dAvgEndX, dAvgEndY))
1075 LogStream <<
m_ulIter <<
": profiles " << nFirstProfile <<
" and " << nSecondProfile <<
" intersect, truncate " << nFirstProfile <<
" since it is an intervention profile" << endl;
1083 LogStream <<
m_ulIter <<
": profiles " << nFirstProfile <<
" and " << nSecondProfile <<
" intersect, truncate " << nSecondProfile <<
" since it is an intervention profile" << endl;
1092 LogStream <<
m_ulIter <<
": profiles " << nFirstProfile <<
" and " << nSecondProfile <<
" intersect, but point {" << dIntersectX <<
", " << dIntersectY <<
"} is already present in profile " << nFirstProfile <<
" as point " << nPoint << endl;
1101 LogStream <<
m_ulIter <<
": profiles " << nFirstProfile <<
" and " << nSecondProfile <<
" intersect, but point {" << dIntersectX <<
", " << dIntersectY <<
"} is already present in profile " << nSecondProfile <<
" as point " << nPoint << endl;
1117 if ((nProf1LineSeg == (nFirstProfileLineSegments - 1)) && (nProf2LineSeg == (nSecondProfileLineSegments - 1)))
1120 MergeProfilesAtFinalLineSegments(nCoast, pFirstProfile, pSecondProfile, nFirstProfileLineSegments, nSecondProfileLineSegments, dIntersectX, dIntersectY, dAvgEndX, dAvgEndY);
1159 else if (nFirstProfileLineSegments < nSecondProfileLineSegments)
1167 else if (nFirstProfileLineSegments > nSecondProfileLineSegments)
1213 for (
unsigned int nCoast = 0; nCoast <
m_VCoast.size(); nCoast++)
1215 for (
int n = 0; n <
m_VCoast[nCoast].nGetNumProfiles(); n++)
1228 m_VCoast[nCoast].pGetProfile(nProfile)->SetTooShort(
true);
1229 LogStream <<
"Profile " << nProfile <<
" is too short, size = " << nSize << endl;
1235 int nXEnd = PtiEnd.
nGetX();
1236 int nYEnd = PtiEnd.
nGetY();
1241 nXEnd =
tMax(nXEnd, 0);
1242 nYEnd =
tMax(nYEnd, 0);
1245 m_VCoast[nCoast].pGetProfile(nProfile)->SetEndPoint(&PtiEnd);
1250 LogStream <<
m_ulIter <<
": coast " << nCoast <<
", profile " << nProfile <<
" is invalid, is too short for depth of closure " <<
m_dDepthOfClosure <<
" at end point [" << nXEnd <<
"][" << nYEnd <<
"] = {" << pPtEnd->
dGetX() <<
", " << pPtEnd->
dGetY() <<
"}, flagging as too short" << endl;
1253 m_VCoast[nCoast].pGetProfile(nProfile)->SetTooShort(
true);
1259 int nValidProfiles = 0;
1262 if (nValidProfiles == 0)
1265 cerr <<
m_ulIter <<
": " <<
ERR <<
"no coastline-normal profiles created" << endl;
1331 for (
int i = 0; i < nProfile1NumSegments; i++)
1333 for (
int j = 0; j < nProfile2NumSegments; j++)
1336 double const dX1 = pVProfile1->
pPtVGetPoints()->at(i).dGetX();
1337 double const dY1 = pVProfile1->
pPtVGetPoints()->at(i).dGetY();
1338 double const dX2 = pVProfile1->
pPtVGetPoints()->at(i + 1).dGetX();
1339 double const dY2 = pVProfile1->
pPtVGetPoints()->at(i + 1).dGetY();
1341 double const dX3 = pVProfile2->
pPtVGetPoints()->at(j).dGetX();
1342 double const dY3 = pVProfile2->
pPtVGetPoints()->at(j).dGetY();
1343 double const dX4 = pVProfile2->
pPtVGetPoints()->at(j + 1).dGetX();
1344 double const dY4 = pVProfile2->
pPtVGetPoints()->at(j + 1).dGetY();
1347 double const dDiffX1 = dX2 - dX1;
1348 double const dDiffY1 = dY2 - dY1;
1349 double const dDiffX2 = dX4 - dX3;
1350 double const dDiffY2 = dY4 - dY3;
1356 dTmp = -dDiffX2 * dDiffY1 + dDiffX1 * dDiffY2;
1359 dS = (-dDiffY1 * (dX1 - dX3) + dDiffX1 * (dY1 - dY3)) / dTmp;
1361 dTmp = -dDiffX2 * dDiffY1 + dDiffX1 * dDiffY2;
1364 dT = (dDiffX2 * (dY1 - dY3) - dDiffY2 * (dX1 - dX3)) / dTmp;
1366 if (dS >= 0 && dS <= 1 && dT >= 0 && dT <= 1)
1369 dXIntersect = dX1 + (dT * dDiffX1);
1370 dYIntersect = dY1 + (dT * dDiffY1);
1373 dXAvgEnd = (dX2 + dX4) / 2;
1374 dYAvgEnd = (dY2 + dY4) / 2;
1377 nProfile1LineSegment = i;
1378 nProfile2LineSegment = j;
1396 int const nProfiles =
m_VCoast[nCoast].nGetNumProfiles();
1407 static bool bDownCoast =
true;
1410 for (
int n = 0; n < nProfiles; n++)
1415 pProfile =
m_VCoast[nCoast].pGetProfileWithDownCoastSeq(n);
1417 pProfile =
m_VCoast[nCoast].pGetProfileWithUpCoastSeq(n);
1437 m_VCoast[nCoast].pGetProfile(nProfile)->SetTooShort(
true);
1440 LogStream <<
m_ulIter <<
": coast " << nCoast <<
", profile " << nProfile <<
" is invalid, has only " << nPoints <<
" points" << endl;
1446 vector<CGeom2DIPoint> VCellsToMark;
1447 vector<bool> bVShared;
1448 bool bTooShort =
false;
1449 bool bTruncatedSameCoast =
false;
1450 bool bHitCoast =
false;
1451 bool bHitLand =
false;
1452 bool bHitIntervention =
false;
1453 bool bHitAnotherProfile =
false;
1455 CreateRasterizedProfile(nCoast, pProfile, &VCellsToMark, &bVShared, bTooShort, bTruncatedSameCoast, bHitCoast, bHitLand, bHitIntervention, bHitAnotherProfile);
1457 if ((bTruncatedSameCoast && (!
ACCEPT_TRUNCATED_PROFILES)) || bTooShort || bHitCoast || bHitLand || bHitIntervention || bHitAnotherProfile || VCellsToMark.size() == 0)
1463 for (
unsigned int k = 0; k < VCellsToMark.size(); k++)
1466 if ((k > 0) && (VCellsToMark[k] ==
m_VCoast[nCoast].pGetProfile(nProfile)->pPtiGetLastCellInProfile()))
1470 int nXTmp = VCellsToMark[k].nGetX();
1471 int nYTmp = VCellsToMark[k].nGetY();
1472 m_pRasterGrid->m_Cell[nXTmp][nYTmp].SetCoastAndProfileID(nCoast, nProfile);
1475 m_VCoast[nCoast].pGetProfile(nProfile)->AppendCellInProfile(nXTmp, nYTmp);
1491 double const dDeepWaterWaveHeight =
m_pRasterGrid->m_Cell[VCellsToMark.back().nGetX()][VCellsToMark.back().nGetY()].dGetCellDeepWaterWaveHeight();
1492 double const dDeepWaterWaveAngle =
m_pRasterGrid->m_Cell[VCellsToMark.back().nGetX()][VCellsToMark.back().nGetY()].dGetCellDeepWaterWaveAngle();
1493 double const dDeepWaterWavePeriod =
m_pRasterGrid->m_Cell[VCellsToMark.back().nGetX()][VCellsToMark.back().nGetY()].dGetCellDeepWaterWavePeriod();
1496 m_VCoast[nCoast].pGetProfile(nProfile)->SetProfileDeepWaterWaveHeight(dDeepWaterWaveHeight);
1497 m_VCoast[nCoast].pGetProfile(nProfile)->SetProfileDeepWaterWaveAngle(dDeepWaterWaveAngle);
1498 m_VCoast[nCoast].pGetProfile(nProfile)->SetProfileDeepWaterWavePeriod(dDeepWaterWavePeriod);
1501 bDownCoast = ! bDownCoast;
1507void CSimulation::CreateRasterizedProfile(
int const nCoast,
CGeomProfile* pProfile, vector<CGeom2DIPoint>* pVIPointsOut, vector<bool>* pbVShared,
bool& bTooShort,
bool& bTruncatedSameCoast,
bool& bHitCoast,
bool& bHitLand,
bool& bHitIntervention,
bool& bHitAnotherProfile)
1513 pVIPointsOut->clear();
1522 for (nSeg = 0; nSeg < nNumSegments; nSeg++)
1546 if (PtiSegStart == PtiSegEnd)
1549 int const nXStart = PtiSegStart.
nGetX();
1550 int const nYStart = PtiSegStart.
nGetY();
1551 int const nXEnd = PtiSegEnd.
nGetX();
1552 int const nYEnd = PtiSegEnd.
nGetY();
1554 bool bShared =
false;
1561 if ((nXStart == nXStartLast) && (nYStart == nYStartLast) && (nXEnd == nXEndLast) && (nYEnd == nYEndLast))
1566 double dXInc = nXEnd - nXStart;
1567 double dYInc = nYEnd - nYStart;
1573 double dX = nXStart;
1574 double dY = nYStart;
1577 for (
int m = 0; m <=
nRound(dLength); m++)
1579 int const nX =
nRound(dX);
1580 int const nY =
nRound(dY);
1589 bTruncatedSameCoast =
true;
1604 int nHitProfileCoast1 =
m_pRasterGrid->m_Cell[nX][nY].nGetProfileCoastID();
1605 int nHitProfileCoast2 =
m_pRasterGrid->m_Cell[nX][nY+1].nGetProfileCoastID();
1607 if ((nHitProfileCoast1 == nCoast) || (nHitProfileCoast2 == nCoast))
1610 bHitAnotherProfile =
true;
1625 int const nHitCoast =
m_pRasterGrid->m_Cell[nX][nY].nGetCoastline();
1647 bHitIntervention =
true;
1660 int const nHitProfile =
m_pRasterGrid->m_Cell[nX][nY].nGetProfileID();
1661 int const nHitProfileCoast =
m_pRasterGrid->m_Cell[nX][nY].nGetProfileCoastID();
1664 if (nCoast == nHitProfileCoast)
1671 bHitAnotherProfile =
true;
1681 pbVShared->push_back(bShared);
1688 nXStartLast = nXStart;
1689 nYStartLast = nYStart;
1693 if (bTruncatedSameCoast)
1697 if (bTruncatedSameCoast)
1699 if (nSeg < (nNumSegments - 1))
1704 int const nLastX = pVIPointsOut->at(pVIPointsOut->size() - 1).nGetX();
1705 int const nLastY = pVIPointsOut->at(pVIPointsOut->size() - 1).nGetY();
1719 if (pVIPointsOut->size() < 3)
1728 LogStream <<
m_ulIter <<
": profile " << nProfile <<
" is invalid, is too short, only " << pVIPointsOut->size() <<
" points, HitLand?" << bHitLand <<
". From [" << pVIPointsOut->at(0).nGetX() <<
"][" << pVIPointsOut->at(0).nGetY() <<
"] = {" <<
dGridCentroidXToExtCRSX(pVIPointsOut->at(0).nGetX()) <<
", " <<
dGridCentroidYToExtCRSY(pVIPointsOut->at(0).nGetY()) <<
"} to [" << pVIPointsOut->at(pVIPointsOut->size() - 1).nGetX() <<
"][" << pVIPointsOut->at(pVIPointsOut->size() - 1).nGetY() <<
"] = {" <<
dGridCentroidXToExtCRSX(pVIPointsOut->at(pVIPointsOut->size() - 1).nGetX()) <<
", " <<
dGridCentroidYToExtCRSY(pVIPointsOut->at(pVIPointsOut->size() - 1).nGetY()) <<
"}" << endl;
1739 int nCombinedLastSeg = 0;
1740 vector<pair<int, int>> prVCombinedProfilesCoincidentProfilesLastSeg;
1744 pair<int, int> prTmp;
1748 bool bFound =
false;
1750 for (
unsigned int m = 0; m < prVCombinedProfilesCoincidentProfilesLastSeg.size(); m++)
1752 if (prVCombinedProfilesCoincidentProfilesLastSeg[m].first == prTmp.first)
1761 prVCombinedProfilesCoincidentProfilesLastSeg.push_back(prTmp);
1768 pair<int, int> prTmp;
1772 bool bFound =
false;
1774 for (
unsigned int m = 0; m < prVCombinedProfilesCoincidentProfilesLastSeg.size(); m++)
1776 if (prVCombinedProfilesCoincidentProfilesLastSeg[m].first == prTmp.first)
1785 prVCombinedProfilesCoincidentProfilesLastSeg.push_back(prTmp);
1791 for (
int m = 0; m < nCombinedLastSeg; m++)
1792 prVCombinedProfilesCoincidentProfilesLastSeg[m].second++;
1796 int const nNumFirstProfileCoincidentProfilesLastSeg =
static_cast<int>(prVFirstProfileCoincidentProfilesLastSeg.size());
1797 int const nNumSecondProfileCoincidentProfilesLastSeg =
static_cast<int>(prVSecondProfileCoincidentProfilesLastSeg.size());
1802 for (
int n = 0; n < nNumFirstProfileCoincidentProfilesLastSeg; n++)
1804 int const nThisProfile = prVFirstProfileCoincidentProfilesLastSeg[n].first;
1813 for (
int n = 0; n < nNumSecondProfileCoincidentProfilesLastSeg; n++)
1815 int const nThisProfile = prVSecondProfileCoincidentProfilesLastSeg[n].first;
1824 for (
int nThisLineSeg = 0; nThisLineSeg < nNumFirstProfileCoincidentProfilesLastSeg; nThisLineSeg++)
1826 int const nThisProfile = prVFirstProfileCoincidentProfilesLastSeg[nThisLineSeg].first;
1835 for (
int m = 0; m < nCombinedLastSeg; m++)
1840 for (
int nThisLineSeg = 0; nThisLineSeg < nNumSecondProfileCoincidentProfilesLastSeg; nThisLineSeg++)
1842 int const nThisProfile = prVSecondProfileCoincidentProfilesLastSeg[nThisLineSeg].first;
1851 for (
int m = 0; m < nCombinedLastSeg; m++)
1922 int const nRet =
nInsertPointIntoProfilesIfNeededThenUpdate(nCoast, pProfileToRetain, dIntersectX, dIntersectY, nProfileToRetainIntersectLineSeg, pProfileToTruncate, nProfileToTruncateIntersectLineSeg, bAlreadyPresent);
1931 vector<CGeom2DPoint> PtVProfileLastPart;
1932 vector<vector<pair<int, int>>> prVLineSegLastPart;
1934 if (bAlreadyPresent)
1950 TruncateProfileAndAppendNew(nCoast, pProfileToTruncate, nProfileToTruncateIntersectLineSeg, &PtVProfileLastPart, &prVLineSegLastPart);
2026 int const nProfileToRetain = pProfileToRetain->
nGetProfileID();
2030 int const nNumCoincident =
static_cast<int>(prVCoincidentProfiles.size());
2031 vector<int> nLineSegAfterIntersect(nNumCoincident, -1);
2034 for (
int nn = 0; nn < nNumCoincident; nn++)
2036 int const nThisProfile = prVCoincidentProfiles[nn].first;
2037 int const nThisLineSeg = prVCoincidentProfiles[nn].second;
2041 if (! bAlreadyPresent)
2047 LogStream <<
WARN <<
m_ulIter <<
": cannot insert a line segment after the final line segment (" << nThisLineSeg <<
") for " << (nThisProfile == nProfileToRetain ?
"main" :
"co-incident") <<
" profile (" << nThisProfile <<
"), abandoning" << endl;
2056 nLineSegAfterIntersect[nn] = nThisLineSeg + 1;
2064 int const nNumToTruncateCoincident =
static_cast<int>(prVToTruncateCoincidentProfiles.size());
2067 for (
int nn = 0; nn < nNumCoincident; nn++)
2069 int const nThisProfile = prVCoincidentProfiles[nn].first;
2076 for (
int nLineSeg = nLineSegAfterIntersect[nn], nIncr = 0; nLineSeg < nNumLineSegs; nLineSeg++, nIncr++)
2087 for (
int m = 0; m < nNumToTruncateCoincident; m++)
2089 int const nProfileToAdd = prVToTruncateCoincidentProfiles[m].first;
2090 int const nProfileToAddLineSeg = prVToTruncateCoincidentProfiles[m].second;
2191 int const nNumCoincident =
static_cast<int>(prVCoincidentProfiles.size());
2193 for (
int nn = 0; nn < nNumCoincident; nn++)
2196 int const nThisProfile = prVCoincidentProfiles[nn].first;
2197 int const nThisProfileLineSeg = prVCoincidentProfiles[nn].second;
2211 for (
unsigned int mm = 0; mm < pPtVProfileLastPart->size(); mm++)
2218 for (
unsigned int mm = 0; mm < pprVLineSegLastPart->size(); mm++)
2220 vector<pair<int, int>> prVTmp = pprVLineSegLastPart->at(mm);
2227 vector<int> nVProfsLineSeg;
2233 int const nProf = pThisProfile->
nGetProf(nSeg, nCoinc);
2236 auto it = find(nVProf.begin(), nVProf.end(), nProf);
2238 if (it == nVProf.end())
2241 nVProf.push_back(nProf);
2242 nVProfsLineSeg.push_back(nProfsLineSeg);
2248 int const nPos =
static_cast<int>(it - nVProf.begin());
2249 int nNewProfsLineSeg = nVProfsLineSeg[nPos];
2252 nVProfsLineSeg[nPos] = nNewProfsLineSeg;
Contains CGeom2DPoint definitions.
Contains CGeom2DIPoint definitions.
vector< CGeom2DPoint > * pPtVGetPoints(void)
Returns the address of the vector which represents this 2D shape.
Geometry class used to represent 2D point objects with integer coordinates.
int nGetY(void) const
Returns the CGeom2DIPoint object's integer Y coordinate.
void SetXY(int const, int const)
The two integer parameters set values for the CGeom2DIPoint object's X and Y coordinates.
int nGetX(void) const
Returns the CGeom2DIPoint object's integer X coordinate.
Geometry class used to represent 2D point objects with floating-point coordinates.
void SetY(double const)
The double parameter sets a value for the CGeom2DIPoint object's Y coordinate.
double dGetY(void) const
Returns the CGeom2DPoint object's double Y coordinate.
double dGetX(void) const
Returns the CGeom2DPoint object's double X coordinate.
void SetX(double const)
The double parameter sets a value for the CGeom2DIPoint object's X coordinate.
int nGetProfsLineSeg(int const, int const) const
Returns the profile's own line segment, given a line segment and the index of the co-incident profile...
bool bFindProfileInCoincidentProfilesOfLastLineSegment(int const)
Returns true if the given profile number is amongst the coincident profiles of the CGeomMultiLine obj...
int nGetNumCoincidentProfilesInLineSegment(int const)
Returns the count of coincident profiles in a specified line segment, or -1 if the line segment does ...
void AppendLineSegment(void)
Appends a new empty line segment.
int nGetProf(int const, int const) const
Returns the profile number, given a line segment and the index of the co-incident profile for that li...
void TruncateLineSegments(int const)
Cuts short the number of line segments.
vector< vector< pair< int, int > > > prVVGetAllLineSegAfter(int const)
Returns a vector of the line segments which succeed the specified line segment number.
void AppendCoincidentProfileToLineSegments(pair< int, int > const)
Appends a coincident profile pair to the CGeomMultiLine object's final line segment.
void SetProfsLineSeg(int const, int const, int const)
Sets a profile's own line segment number, given a line segment and the index of the co-incident profi...
vector< pair< int, int > > * pprVGetPairedCoincidentProfilesForLineSegment(int const)
Returns a vector of pairs (a line segment)
void AddCoincidentProfileToExistingLineSegment(int const, int const, int const)
Adds a coincident profile to a pre-existing line segment of the CGeomMultiLine object.
int nGetNumLineSegments(void) const
Appends a line segment which then inherits from the preceding line segments.
Geometry class used to represent coast profile objects.
void TruncateProfile(int const)
Truncates the profile's CGeomLine (external CRS points)
void SetEndOfCoast(bool const)
Sets a switch to indicate whether this is an end-of-coast profile.
void AppendPointInProfile(double const, double const)
Appends a point (external CRS) to the profile.
int nGetProfileID(void) const
Returns the profile's this-coast ID.
void SetProfileDeepWaterWavePeriod(double const)
Sets the deep-water wave period for this profile.
bool bIsPointInProfile(double const, double const)
Removes a line segment from the profile.
void SetHitCoast(bool const)
Sets a switch which indicates whether this profile has hit a coast.
bool bProfileOKIncTruncated(void) const
Returns true if this is a problem-free profile, and is not a start-of-coast and is not an end-of-coas...
void SetStartOfCoast(bool const)
Sets a switch to indicate whether this is a start-of-coast profile.
void SetUpCoastAdjacentProfile(CGeomProfile *)
Sets the up-coast adjacent profile.
void SetHitAnotherProfile(bool const)
Sets a switch which indicates whether this profile hits another profile badly.
vector< CGeom2DPoint > PtVGetThisPointAndAllAfter(int const)
Returns a given external CRS point from the profile, and all points after this.
bool bInsertIntersection(double const, double const, int const)
Inserts an intersection (at a point specified in external CRS, with a line segment) into the profile.
void SetTruncatedSameCoast(bool const)
Sets a switch which indicates whether this profile is truncated, due to hitting another profile from ...
bool bProfileOK(void) const
Returns true if this is a problem-free profile, and is not a start-of-coast and is not an end-of-coas...
void SetProfileDeepWaterWaveHeight(double const)
Sets the deep-water wave height for this profile.
void SetDownCoastAdjacentProfile(CGeomProfile *)
Sets the down-coast adjacent profile.
bool bIsIntervention(void) const
Returns true if this is an intervention profile.
void SetHitIntervention(bool const)
Sets a switch which indicates whether this profile has hit an intervention.
CGeom2DPoint * pPtGetPointInProfile(int const)
Returns a single point (external CRS) from the profile.
void AppendCellInProfile(CGeom2DIPoint const *)
Appends a cell (grid CRS) to the profile.
void SetPointInProfile(int const, double const, double const)
Sets a single point (external CRS) in the profile.
int nGetProfileSize(void) const
Returns the number of external CRS points in the profile (only two, initally; and always just two for...
CGeom2DIPoint * pPtiGetStartPoint(void)
Returns a pointer to the location of the cell (grid CRS) on which the profile starts.
bool bStartOfCoast(void) const
Returns the switch to indicate whether this is a start-of-coast profile.
void SetEndPoint(CGeom2DIPoint const *)
Sets the the location of the cell (grid CRS) on which the profile ends.
void SetHitLand(bool const)
Sets a switch which indicates whether this profile has hit land.
void SetPointsInProfile(vector< CGeom2DPoint > const *)
Sets points (external CRS) in the profile. Note that only two points, the start and end point,...
bool bEndOfCoast(void) const
Returns the switch to indicate whether this is an end-of-coast profile.
void SetProfileDeepWaterWaveAngle(double const)
Sets the deep-water wave orientation for this profile.
void SetTooShort(bool const)
Sets a switch which indicates whether this profile is too short to be useful.
int m_nLogFileDetail
The level of detail in the log file output. Can be LOG_FILE_LOW_DETAIL, LOG_FILE_MIDDLE_DETAIL,...
double m_dThisIterSWL
The still water level for this timestep (this includes tidal changes and any long-term SWL change)
CGeomRasterGrid * m_pRasterGrid
Pointer to the raster grid object.
int m_nXGridSize
The size of the grid in the x direction.
void LocateAndCreateProfiles(int const, int &, vector< bool > *, vector< pair< int, double > > const *)
For a single coastline, locate the start points for all coastline-normal profiles (except the grid-ed...
CGeom2DIPoint PtiExtCRSToGridRound(CGeom2DPoint const *) const
Transforms a pointer to a CGeom2DPoint in the external CRS to the equivalent CGeom2DIPoint in the ras...
vector< CRWCoast > m_VCoast
The coastline objects.
double m_dCoastNormalLength
Length of the coastline-normal profiles, in m.
void MarkProfilesOnGrid(int const, int &)
For this coastline, marks all coastline-normal profiles (apart from the two 'special' ones at the sta...
static CGeom2DPoint PtChooseEndPoint(int const, CGeom2DPoint const *, CGeom2DPoint const *, double const, double const, double const, double const)
Choose which end point to use for the coastline-normal profile.
int m_nYGridSize
The size of the grid in the y direction.
normal_distribution< double > m_dGetFromUnitNormalDist
c++11 unit normal distribution (mean = 0, stdev = 1)
void CreateRasterizedProfile(int const, CGeomProfile *, vector< CGeom2DIPoint > *, vector< bool > *, bool &, bool &, bool &, bool &, bool &, bool &)
Given a pointer to a coastline-normal profile, returns an output vector of cells which are 'under' ev...
static CGeom2DPoint PtAverage(CGeom2DPoint const *, CGeom2DPoint const *)
Returns a point (external CRS) which is the average of (i.e. is midway between) two other external CR...
void KeepWithinValidGrid(int &, int &) const
int m_nCoastNormalSpacing
Average spacing between coastline normals, measured in cells.
void MergeProfilesAtFinalLineSegments(int const, CGeomProfile *, CGeomProfile *, int const, int const, double const, double const, double const, double const)
Merges two profiles which intersect at their final (most seaward) line segments, seaward of their poi...
int m_nCoastNormalInterventionSpacing
Average spacing between coastline normals on interventions, measured in cells.
void CheckForIntersectingProfiles(void)
Checks all coastline-normal profiles for intersection, and modifies those that intersect.
double dGridCentroidYToExtCRSY(int const) const
Given the integer Y-axis ordinate of a cell in the raster grid CRS, returns the external CRS Y-axis o...
int nCreateProfile(int const, int const, int const, int const, bool const, CGeom2DIPoint const *)
Creates a single coastline-normal profile (which may be an intervention profile)
vector< CGeom2DIPoint > m_VEdgeCell
Edge cells.
static bool bCheckForIntersection(CGeomProfile *const, CGeomProfile *const, int &, int &, double &, double &, double &, double &)
Checks all line segments of a pair of coastline-normal profiles for intersection. If the lines inters...
int nGetCoastNormalEndPoint(int const, int const, int const, CGeom2DPoint const *, double const, CGeom2DPoint *, CGeom2DIPoint *, bool const)
Finds the end point of a coastline-normal line, given the start point on the vector coastline....
void TruncateProfileAndAppendNew(int const, CGeomProfile *, int const, vector< CGeom2DPoint > const *, vector< vector< pair< int, int > > > const *)
Truncate a profile at the point of intersection, and do the same for all its co-incident profiles.
double dExtCRSXToGridX(double const) const
int nInsertPointIntoProfilesIfNeededThenUpdate(int const, CGeomProfile *, double const, double const, int const, CGeomProfile *, int const, bool const)
Inserts an intersection point into the profile that is to be retained, if that point is not already p...
default_random_engine m_Rand[NUMBER_OF_RNGS]
The c++11 random number generators.
int nCreateAllProfiles(void)
Create coastline-normal profiles for all coastlines. The first profiles are created 'around' the most...
double dExtCRSYToGridY(double const) const
Transforms a Y-axis ordinate in the external CRS to the equivalent Y-axis ordinate in the raster grid...
double m_dCoastNormalRandSpacingFactor
Random factor for spacing of along-coast normals.
bool bIsWithinValidGrid(int const, int const) const
static void AppendEnsureNoGap(vector< CGeom2DIPoint > *, CGeom2DIPoint const *)
Appends a CGeom2DIPoint to a vector<CGeom2DIPoint>, making sure that the new end point touches the pr...
unsigned long m_ulIter
The number of the current iteration (time step)
double dGridCentroidXToExtCRSX(int const) const
int nLocateAndCreateGridEdgeProfile(bool const, int const, int &)
Creates a 'special' profile at each end of a coastline, at the edge of the raster grid....
void TruncateOneProfileRetainOtherProfile(int const, CGeomProfile *, CGeomProfile *, double, double, int, int, bool const)
Truncates one intersecting profile at the point of intersection, and retains the other profile.
double m_dCellSide
Length of a cell side (in external CRS units)
double m_dDepthOfClosure
Depth of closure (in m) TODO 007 can be calculated using Hallermeier, R.J. (1978) or Birkemeier (1985...
int nCheckAndMarkAllProfiles(void)
Check all coastline-normal profiles and modify the profiles if they intersect, then mark valid profil...
vector< int > m_VEdgeCellEdge
The grid edge that each edge cell belongs to.
This file contains global definitions for CoastalME.
int const RTN_ERR_COAST_CANT_FIND_EDGE_CELL
int const DIRECTION_DOWNCOAST
bool const ACCEPT_TRUNCATED_PROFILES
int const LOG_FILE_MIDDLE_DETAIL
bool bFPIsEqual(const T d1, const T d2, const T dEpsilon)
int const RTN_ERR_PROFILE_END_INSUFFICIENT_DEPTH
int const PROFILE_CHECK_DIST_FROM_COAST
int const RTN_ERR_NO_PROFILES_2
int const RTN_ERR_PROFILE_ENDPOINT_IS_INLAND
int const RTN_ERR_NO_PROFILES_1
int const RTN_ERR_CANNOT_INSERT_POINT
int const DIRECTION_UPCOAST
string strDblToStr(const T &t)
int const RTN_ERR_NO_SOLUTION_FOR_ENDPOINT
int const LF_CAT_INTERVENTION
Contains CRWCoast definitions.
Contains CSimulation definitions.
int nRound(double const d)
Version of the above that returns an int.