Commit 7b75de59a7826bb70bfb861958fc2decc1de5259

Authored by Josh Klontz
1 parent 8db62b7d

dropped BBDLDAAlignmentTransform

openbr/plugins/classification/lda.cpp
... ... @@ -15,9 +15,6 @@
15 15 * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
16 16  
17 17 #include <Eigen/Dense>
18   -#include <opencv2/imgproc/imgproc.hpp>
19   -#include <opencv2/highgui/highgui.hpp>
20   -
21 18 #include <openbr/plugins/openbr_internal.h>
22 19  
23 20 #include <openbr/core/common.h>
... ... @@ -793,439 +790,6 @@ class WCDATransform : public Transform
793 790  
794 791 BR_REGISTER(Transform, WCDATransform)
795 792  
796   -// For use in BBLDAAlignmentTransform
797   -// A single decision boundary with gaussian distributed responses
798   -struct DecisionBoundary
799   -{
800   - VectorXf function;
801   - float positiveMean, negativeMean, positiveSD, negativeSD;
802   -};
803   -
804   -QDataStream &operator<<(QDataStream &stream, const DecisionBoundary &db)
805   -{
806   - return stream << db.function << db.positiveMean << db.negativeMean << db.positiveSD << db.negativeSD;
807   -}
808   -
809   -QDataStream &operator>>(QDataStream &stream, DecisionBoundary &db)
810   -{
811   - return stream >> db.function >> db.positiveMean >> db.negativeMean >> db.positiveSD >> db.negativeSD;
812   -}
813   -
814   -/*!
815   - * \ingroup transforms
816   - * \brief Boosted Binary LDA classifier for local alignment.
817   - * \author Unknown \cite Unknown
818   - */
819   -class BBLDAAlignmentTransform : public MetaTransform
820   -{
821   - Q_OBJECT
822   - Q_PROPERTY(QList<int> anchors READ get_anchors WRITE set_anchors RESET reset_anchors STORED false) // A list of two indicies to provide initial rotation and scaling of training data
823   - Q_PROPERTY(QString detectionKey READ get_detectionKey WRITE set_detectionKey RESET reset_detectionKey STORED false) // Metadata key containing the detection bounding box
824   - Q_PROPERTY(int iad READ get_iad WRITE set_iad RESET reset_iad STORED false) // Inter-anchor distance in pixels
825   - Q_PROPERTY(float radius READ get_radius WRITE set_radius RESET reset_radius STORED false) // Kernel size as a fraction of IAD
826   - Q_PROPERTY(int resolution READ get_resolution WRITE set_resolution RESET reset_resolution STORED false) // Project resolution for sampling rotation and scale
827   - BR_PROPERTY(QList<int>, anchors, QList<int>())
828   - BR_PROPERTY(QString, detectionKey, "FrontalFace")
829   - BR_PROPERTY(int, iad, 16)
830   - BR_PROPERTY(float, radius, 2)
831   - BR_PROPERTY(int, resolution, 20)
832   -
833   - typedef Matrix<quint8, Dynamic, Dynamic> MatrixX8U;
834   -
835   - struct RegistrationMatrix
836   - {
837   - Mat src;
838   - Point2f center;
839   - double angle, scale;
840   -
841   - RegistrationMatrix(const Mat &src = Mat(), const Point2f &center = Point2f(), double angle = 0, double scale = 1)
842   - : src(src)
843   - , center(center)
844   - , angle(angle)
845   - , scale(scale)
846   - {}
847   -
848   - Mat calculateRotationMatrix2D(int size) const
849   - {
850   - Mat rotationMatrix2D = getRotationMatrix2D(center, angle, scale);
851   - rotationMatrix2D.at<double>(0, 2) -= (center.x - size / 2.0); // Adjust the center to the middle of the dst image
852   - rotationMatrix2D.at<double>(1, 2) -= (center.y - size / 2.0);
853   - return rotationMatrix2D;
854   - }
855   -
856   - Point2f invert(double x, double y, int size) const
857   - {
858   - Mat m;
859   - invertAffineTransform(calculateRotationMatrix2D(size), m);
860   -
861   - Mat src(1, 1, CV_64FC2);
862   - src.ptr<double>()[0] = y; // According to the docs these should be specified in the opposite order ...
863   - src.ptr<double>()[1] = x; // ... however running it seems to suggest this is the correct way
864   -
865   - Mat dst;
866   - transform(src, dst, m);
867   - return Point2f(dst.ptr<double>()[0], dst.ptr<double>()[1]);
868   - }
869   -
870   - MatrixX8U warp(int size) const
871   - {
872   - const Mat rotationMatrix2D = calculateRotationMatrix2D(size);
873   - Mat dst;
874   - warpAffine(src, dst, rotationMatrix2D, Size(size, size), INTER_LINEAR, BORDER_REFLECT);
875   - return Map<MatrixX8U>(dst.ptr<quint8>(), size, size);
876   - }
877   - };
878   -
879   - float rotationMax, scaleMin, scaleMax, xTranslationMax, yTranslationMin, yTranslationMax;
880   - QList<DecisionBoundary> decisionBoundaries;
881   -
882   - static void show(const char *name, const MatrixX8U &data)
883   - {
884   - const Mat mat(data.rows(), data.cols(), CV_8UC1, (void*) data.data());
885   - imshow(name, mat);
886   - waitKey(-1);
887   - }
888   -
889   - static void show(const char *name, MatrixXf data)
890   - {
891   - data.array() -= data.minCoeff();
892   - data.array() *= (255 / data.maxCoeff());
893   - show(name, MatrixX8U(data.cast<quint8>()));
894   - }
895   -
896   - static void show(const char *name, VectorXf data)
897   - {
898   - MatrixXf resized = data;
899   - const int size = int(sqrtf(float(data.rows())));
900   - resized.resize(size, size);
901   - show(name, resized);
902   - }
903   -
904   - static MatrixXf getSample(const MatrixXf &registered, int j, int k, int kernelSize)
905   - {
906   - MatrixXf sample(registered.block(j, k, kernelSize, kernelSize));
907   - const float norm = sample.norm();
908   - if (norm > 0)
909   - sample /= norm;
910   - return sample;
911   - }
912   -
913   - static MatrixXf getSample(const QList<MatrixXf> &registered, int i, int j, int k, int kernelSize)
914   - {
915   - return getSample(registered[i], j, k, kernelSize);
916   - }
917   -
918   - static float IADError(const MatrixXf &responses, const MatrixXf &mask)
919   - {
920   - const int responseSize = sqrtf(responses.cols());
921   - int numImages = 0;
922   - float totalError = 0;
923   - for (int i=0; i<responses.rows(); i++) {
924   - float maxResponse = -std::numeric_limits<float>::max();
925   - int maxX = -1;
926   - int maxY = -1;
927   - for (int j=0; j<responseSize; j++)
928   - for (int k=0; k<responseSize; k++)
929   - if (mask(i, j*responseSize + k) > 0) {
930   - const float response = responses(i, j*responseSize+k);
931   - if (response > maxResponse) {
932   - maxResponse = response;
933   - maxY = j;
934   - maxX = k;
935   - }
936   - }
937   - totalError += sqrtf(powf(maxX - responseSize/2, 2) + powf(maxY - responseSize/2, 2)) / responseSize;
938   - numImages++;
939   - }
940   - return totalError / numImages;
941   - }
942   -
943   - static bool isPositive(int j, int k, int responseSize)
944   - {
945   - return (abs(j - responseSize/2) <= 0) && (abs(k - responseSize/2) <= 0);
946   - }
947   -
948   - static MatrixXf /* output responses */ computeBoundaryRecursive(const QList<MatrixXf> &registered, int kernelSize, const MatrixXf &prior /* input responses */, QList<DecisionBoundary> &boundaries)
949   - {
950   - const int numImages = registered.size();
951   - const int responseSize = registered.first().cols() - kernelSize + 1;
952   - int positiveSamples = 0;
953   - int negativeSamples = 0;
954   -
955   - // Compute weights, a weight of zero operates as a mask
956   - MatrixXf weights(numImages, responseSize*responseSize);
957   - float totalPositiveWeight = 0, totalNegativeWeight = 0;
958   - for (int i=0; i<numImages; i++)
959   - for (int j=0; j<responseSize; j++)
960   - for (int k=0; k<responseSize; k++) {
961   - const float probability = prior(i, j*responseSize+k);
962   - const bool positive = isPositive(j, k, responseSize);
963   - // Weight by probability ...
964   - float weight = positive ? (2 - probability) // Subtracting from a number greater than 1 helps ensure that we don't overfit outliers
965   - : probability;
966   - // ... and by distance
967   - const float distance = sqrtf(powf(j - responseSize/2, 2) + powf(k - responseSize/2, 2));
968   - if (positive) weight *= 1 / (distance + 1);
969   - else weight *= distance;
970   - weights(i, j*responseSize+k) = weight;
971   - if (weight > 0) {
972   - if (positive) { totalPositiveWeight += weight; positiveSamples++; }
973   - else { totalNegativeWeight += weight; negativeSamples++; }
974   - }
975   - }
976   -
977   - // Normalize weights to preserve class sample ratio
978   - const float positiveWeightAdjustment = positiveSamples / totalPositiveWeight;
979   - const float negativeWeightAdjustment = negativeSamples / totalNegativeWeight;
980   - for (int i=0; i<numImages; i++)
981   - for (int j=0; j<responseSize; j++)
982   - for (int k=0; k<responseSize; k++)
983   - weights(i, j*responseSize+k) *= isPositive(j, k, responseSize) ? positiveWeightAdjustment : negativeWeightAdjustment;
984   -
985   - // Compute weighted class means
986   - MatrixXf positiveMean = MatrixXf::Zero(kernelSize, kernelSize);
987   - MatrixXf negativeMean = MatrixXf::Zero(kernelSize, kernelSize);
988   - for (int i=0; i<numImages; i++)
989   - for (int j=0; j<responseSize; j++)
990   - for (int k=0; k<responseSize; k++) {
991   - const MatrixXf sample(getSample(registered, i, j, k, kernelSize));
992   - const float weight = weights(i, j*responseSize+k);
993   - if (weight > 0) {
994   - if (isPositive(j, k, responseSize)) positiveMean += sample * weight;
995   - else negativeMean += sample * weight;
996   - }
997   - }
998   - positiveMean /= positiveSamples;
999   - negativeMean /= negativeSamples;
1000   -
1001   - // Compute weighted scatter matrix and decision boundary
1002   - MatrixXf scatter = MatrixXf::Zero(kernelSize*kernelSize, kernelSize*kernelSize);
1003   - for (int i=0; i<numImages; i++)
1004   - for (int j=0; j<responseSize; j++)
1005   - for (int k=0; k<responseSize; k++) {
1006   - const float weight = weights(i, j*responseSize+k);
1007   - if (weight > 0) {
1008   - const MatrixXf &centeredSampleMatrix = getSample(registered, i, j, k, kernelSize) - (isPositive(j, k, responseSize) ? positiveMean : negativeMean);
1009   - const Map<const VectorXf> centeredSample(centeredSampleMatrix.data(), kernelSize*kernelSize);
1010   - scatter.noalias() += centeredSample * centeredSample.transpose() * weights(i, j*responseSize+k);
1011   - }
1012   - }
1013   - scatter /= (positiveSamples + negativeSamples); // normalize for numerical stability
1014   - DecisionBoundary db;
1015   - db.function = scatter.inverse() * VectorXf(positiveMean - negativeMean); // fisher linear discriminant
1016   -
1017   - // Compute response values to decision boundary
1018   - MatrixXf posterior(numImages, responseSize*responseSize);
1019   - for (int i=0; i<numImages; i++)
1020   - for (int j=0; j<responseSize; j++)
1021   - for (int k=0; k<responseSize; k++)
1022   - if (weights(i, j*responseSize + k) > 0) {
1023   - const MatrixXf sample(getSample(registered, i, j, k, kernelSize));
1024   - posterior(i, j*responseSize + k) = db.function.transpose() * Map<const VectorXf>(sample.data(), kernelSize*kernelSize);
1025   - }
1026   -
1027   - // Compute class response means
1028   - db.positiveMean = 0;
1029   - db.negativeMean = 0;
1030   - for (int i=0; i<numImages; i++)
1031   - for (int j=0; j<responseSize; j++)
1032   - for (int k=0; k<responseSize; k++)
1033   - if (weights(i, j*responseSize + k) > 0) {
1034   - const float response = posterior(i, j*responseSize + k);
1035   - if (isPositive(j, k, responseSize)) db.positiveMean += response;
1036   - else db.negativeMean += response;
1037   - }
1038   - db.positiveMean /= positiveSamples;
1039   - db.negativeMean /= negativeSamples;
1040   -
1041   - // Compute class response standard deviations
1042   - db.positiveSD = 0;
1043   - db.negativeSD = 0;
1044   - for (int i=0; i<numImages; i++)
1045   - for (int j=0; j<responseSize; j++)
1046   - for (int k=0; k<responseSize; k++)
1047   - if (weights(i, j*responseSize + k) > 0) {
1048   - const float response = posterior(i, j*responseSize + k);
1049   - if (isPositive(j, k, responseSize)) db.positiveSD += powf(response-db.positiveMean, 2.f);
1050   - else db.negativeSD += powf(response-db.negativeMean, 2.f);
1051   - }
1052   - db.positiveSD = sqrtf(db.positiveSD / positiveSamples);
1053   - db.negativeSD = sqrtf(db.negativeSD / negativeSamples);
1054   -
1055   - // Normalize responses and propogating prior probabilities
1056   - for (int i=0; i<numImages; i++)
1057   - for (int j=0; j<responseSize; j++)
1058   - for (int k=0; k<responseSize; k++) {
1059   - float &response = posterior(i, j*responseSize + k);
1060   - if (weights(i, j*responseSize + k) > 0) {
1061   - const float positiveDensity = 1 / (sqrtf(2 * CV_PI) * db.positiveSD) * expf(-0.5 * powf((response - db.positiveMean) / db.positiveSD, 2.f));
1062   - const float negativeDensity = 1 / (sqrtf(2 * CV_PI) * db.negativeSD) * expf(-0.5 * powf((response - db.negativeMean) / db.negativeSD, 2.f));
1063   - response = prior(i, j*responseSize+k) * positiveDensity / (positiveDensity + negativeDensity);
1064   - } else {
1065   - response = 0;
1066   - }
1067   - }
1068   - const float previousMeanError = IADError(prior, weights);
1069   - const float meanError = IADError(posterior, weights);
1070   - qDebug() << "Samples positive & negative:" << positiveSamples << negativeSamples;
1071   - qDebug() << "Positive mean & stddev:" << db.positiveMean << db.positiveSD;
1072   - qDebug() << "Negative mean & stddev:" << db.negativeMean << db.negativeSD;
1073   - qDebug() << "Mean error before & after:" << previousMeanError << meanError;
1074   -
1075   - if (meanError < previousMeanError) {
1076   - boundaries.append(db);
1077   - return computeBoundaryRecursive(registered, kernelSize, posterior, boundaries);
1078   - } else {
1079   - return prior;
1080   - }
1081   - }
1082   -
1083   - void train(const TemplateList &data)
1084   - {
1085   - QList<RegistrationMatrix> registrationMatricies;
1086   - {
1087   - if (Globals->verbose)
1088   - qDebug("Preprocessing training data...");
1089   -
1090   - rotationMax = -std::numeric_limits<float>::max();
1091   - scaleMin = std::numeric_limits<float>::max();
1092   - scaleMax = -std::numeric_limits<float>::max();
1093   - xTranslationMax = -std::numeric_limits<float>::max();
1094   - yTranslationMin = std::numeric_limits<float>::max();
1095   - yTranslationMax = -std::numeric_limits<float>::max();
1096   -
1097   - foreach (const Template &t, data) {
1098   - const QRectF detection = t.file.get<QRectF>(detectionKey);
1099   - const QList<QPointF> points = t.file.points();
1100   - const float length = sqrt(pow(points[anchors[1]].x() - points[anchors[0]].x(), 2) +
1101   - pow(points[anchors[1]].y() - points[anchors[0]].y(), 2));
1102   - const float rotation = atan2(points[anchors[1]].y() - points[anchors[0]].y(),
1103   - points[anchors[1]].x() - points[anchors[0]].x()) * 180 / CV_PI;
1104   - const float scale = length / detection.width();
1105   - const float xCenter = (points[anchors[0]].x() + points[anchors[1]].x()) / 2;
1106   - const float yCenter = (points[anchors[0]].y() + points[anchors[1]].y()) / 2;
1107   - const float xTranslation = (xCenter - detection.x() - detection.width() /2) / detection.width();
1108   - const float yTranslation = (yCenter - detection.y() - detection.height()/2) / detection.width();
1109   -
1110   - if (detection.contains(xCenter, yCenter) /* Sometimes the face detector gets the wrong face */) {
1111   - rotationMax = max(rotationMax, fabsf(rotation));
1112   - scaleMin = min(scaleMin, scale);
1113   - scaleMax = max(scaleMax, scale);
1114   - xTranslationMax = max(xTranslationMax, fabsf(xTranslation));
1115   - yTranslationMin = min(yTranslationMin, yTranslation);
1116   - yTranslationMax = max(yTranslationMax, yTranslation);
1117   - }
1118   -
1119   - registrationMatricies.append(RegistrationMatrix(t, Point2f(xCenter, yCenter), rotation, iad / length));
1120   - }
1121   - }
1122   -
1123   - if (Globals->verbose)
1124   - qDebug("Learning affine model ...");
1125   -
1126   - // Construct the registered training data for the landmark
1127   - const int regionSize = 2 * iad * radius; // Train on a search size that is twice to the kernel size
1128   - QList<MatrixXf> registered;
1129   - foreach (const RegistrationMatrix &registrationMatrix, registrationMatricies)
1130   - registered.append(registrationMatrix.warp(regionSize).cast<float>());
1131   -
1132   - // Boosted LDA
1133   - const int numImages = registered.size();
1134   - const int kernelSize = iad * radius;
1135   - const int responseSize = regionSize - kernelSize + 1;
1136   - const MatrixXf prior = MatrixXf::Ones(numImages, responseSize*responseSize);
1137   - const MatrixXf responses = computeBoundaryRecursive(registered, kernelSize, prior, decisionBoundaries);
1138   -
1139   -// for (int i=0; i<numImages; i++) {
1140   -// float maxResponse = -std::numeric_limits<float>::max();
1141   -// int maxX = -1;
1142   -// int maxY = -1;
1143   -// for (int j=0; j<responseSize; j++)
1144   -// for (int k=0; k<responseSize; k++) {
1145   -// const float response = responses(i, j*responseSize+k);
1146   -// if (response > maxResponse) {
1147   -// maxResponse = response;
1148   -// maxY = j;
1149   -// maxX = k;
1150   -// }
1151   -// }
1152   -
1153   -// MatrixXf draw = registered[i];
1154   -// const int r = 3;
1155   -// maxX += kernelSize / 2;
1156   -// maxY += kernelSize / 2;
1157   -// for (int i=-r; i<=r; i++)
1158   -// draw(regionSize/2 + i, regionSize/2) *= 2;
1159   -// for (int i=-r; i<=r; i++)
1160   -// draw(regionSize/2, regionSize/2 + i) *= 2;
1161   -// for (int i=-r; i<=r; i++)
1162   -// draw(maxX + i, maxY) /= 2;
1163   -// for (int i=-r; i<=r; i++)
1164   -// draw(maxX, maxY + i) /= 2;
1165   -// show("draw", draw);
1166   -// show("responses", VectorXf(responses.row(i)));
1167   -// }
1168   -
1169   - if (Globals->verbose)
1170   - qDebug("Learned %d function(s) with error %g", decisionBoundaries.size(), IADError(responses, prior));
1171   - }
1172   -
1173   - void project(const Template &src, Template &dst) const
1174   - {
1175   - dst = src;
1176   - const QRectF detection = src.file.get<QRectF>(detectionKey);
1177   - RegistrationMatrix registrationMatrix(src, OpenCVUtils::toPoint(detection.center()));
1178   -
1179   - const int regionSize = iad * (1 + radius);
1180   - const int kernelSize = iad * radius;
1181   - const int responseSize = regionSize - kernelSize + 1;
1182   - const float rotationStep = (2 * rotationMax) / (resolution - 1);
1183   - const float scaleStep = (scaleMax - scaleMin) / (resolution - 1);
1184   -
1185   - float bestResponse = -std::numeric_limits<float>::max(), bestRotation = 0, bestScale = 0;
1186   - int bestX = 0, bestY =0;
1187   -
1188   - for (float rotation = -rotationMax; rotation <= rotationMax; rotation += rotationStep) {
1189   - registrationMatrix.angle = rotation;
1190   - for (float scale = scaleMin; scale <= scaleMax; scale += scaleStep) {
1191   - registrationMatrix.scale = iad / (scale * detection.width());
1192   - const MatrixXf registered = registrationMatrix.warp(regionSize).cast<float>();
1193   -
1194   - for (int j=0; j<responseSize; j++)
1195   - for (int k=0; k<responseSize; k++) {
1196   - const MatrixXf sample = getSample(registered, j, k, kernelSize);
1197   - float posterior = 1;
1198   - foreach (const DecisionBoundary &db, decisionBoundaries) {
1199   - const float response = db.function.transpose() * Map<const VectorXf>(sample.data(), kernelSize*kernelSize);
1200   - const float positiveDensity = 1 / (sqrtf(2 * CV_PI) * db.positiveSD) * expf(-0.5 * powf((response - db.positiveMean) / db.positiveSD, 2.f));
1201   - const float negativeDensity = 1 / (sqrtf(2 * CV_PI) * db.negativeSD) * expf(-0.5 * powf((response - db.negativeMean) / db.negativeSD, 2.f));
1202   - posterior *= positiveDensity / (positiveDensity + negativeDensity);
1203   - }
1204   - if (posterior > bestResponse) {
1205   - bestResponse = posterior;
1206   - bestX = k + kernelSize/2;
1207   - bestY = j + kernelSize/2;
1208   - bestRotation = rotation;
1209   - bestScale = scale;
1210   - }
1211   - }
1212   - }
1213   - }
1214   - }
1215   -
1216   - void store(QDataStream &stream) const
1217   - {
1218   - stream << rotationMax << scaleMin << scaleMax << xTranslationMax << yTranslationMin << yTranslationMax << decisionBoundaries;
1219   - }
1220   -
1221   - void load(QDataStream &stream)
1222   - {
1223   - stream >> rotationMax >> scaleMin >> scaleMax >> xTranslationMax >> yTranslationMin >> yTranslationMax >> decisionBoundaries;
1224   - }
1225   -};
1226   -
1227   -BR_REGISTER(Transform, BBLDAAlignmentTransform)
1228   -
1229 793 } // namespace br
1230 794  
1231 795 #include "classification/lda.moc"
... ...