diff --git a/.gitignore b/.gitignore index a94e449..1f678d9 100644 --- a/.gitignore +++ b/.gitignore @@ -1,2 +1,2 @@ -.DS_Store/ +.DS_Store __pycache__/ \ No newline at end of file diff --git a/app.py b/app.py new file mode 100644 index 0000000..dbdfb3f --- /dev/null +++ b/app.py @@ -0,0 +1,41 @@ +from flask import Flask, request, render_template, send_from_directory, url_for + +from flask_wtf import FlaskForm +from flask_wtf.file import FileField, FileRequired, FileAllowed +from wtforms import SubmitField +from flask_uploads import UploadSet, IMAGES, configure_uploads + +app = Flask(__name__) + +app.config['SECRET_KEY'] = 'secret' +app.config['UPLOADED_PHOTOS_DEST'] = 'uploads' + +photos = UploadSet('photos', IMAGES) +configure_uploads(app, photos) + +class PendantForm(FlaskForm): + photo = FileField('file', validators=[ + FileAllowed(photos, 'Only images are allowed'), + FileRequired('File field should not be empty'), + ]) + submit = SubmitField('Compute') + +@app.route('/upload/') +def get_file(filename): + return send_from_directory(app.config['UPLOADED_PHOTOS_DEST'], filename) + +@app.route('/', methods=['GET', 'POST']) +def index(): + form = PendantForm(); + + if form.validate_on_submit(): + filename = photos.save(form.photo.data) + file_url = url_for('get_file', filename=filename) + else: + file_url = None + return render_template('index.html', form=form, file_url=file_url) + + + +if __name__ == '__main__': + app.run() diff --git a/scripts/img_V01.py b/scripts/img_V01.py new file mode 100644 index 0000000..06234ff --- /dev/null +++ b/scripts/img_V01.py @@ -0,0 +1,799 @@ +#!/usr/bin/env python + +# Surface Tension Calculation via Pendant Drop Method (Bashforth-Adams) +# Created by : Ren Tristan A. de la Cruz +# Created on : 2016 July 26 +# Last Modified : 2025 January 27 +# Contact : rentristandelacruz@gmail.com +# Contact : radelacruz@up.edu.ph + +# ========== IMPORTS ========== + +import math +import copy +import time +import sys +import cv2 +import multiprocessing +import numpy as np +import scipy.optimize #as spo + +# ========== READ THE IMAGE & EXTRACT THE EDGE COORDINATES ========== + +def ExtractEdgeCoordinate(_ImageFilePath, _OutputFileName, _SyringeWidth): + + # Detect the edges of the image. + # All edge pixels will white(255) and rest of the pixels will be black(0). + Image = cv2.imread(_ImageFilePath, cv2.IMREAD_COLOR) + EdgeImage = cv2.Canny(Image, 100, 200) + + # Find the first edge point (white pixel). + # It is assumed that first edge point part of the syringe. + NumOfRow, NumOfCol = EdgeImage.shape + EdgePoint = (0, 0) + EdgePointFlag = 0 + for RowIndex in range(NumOfRow): + for ColIndex in range(NumOfCol): + if EdgeImage[RowIndex, ColIndex] == 255: + EdgePoint = (RowIndex, ColIndex) + EdgePointFlag = 1 + break + if EdgePointFlag == 1: + break + + + # Color the all the syringe and droplet pixels differently (70). + # Determine the left-most, right-most, and bottom droplet pixels. + EdgeImage[EdgePoint] = 70 + DropletEdgePoints = [] + DropletEdgePoints.append(EdgePoint) + LeftMostCol = 9999999999 + RightMostCol = 0 + BottomMostRow = 0 + BottomMostRow2 = 0 + CenterCol = 0 + + while (len(DropletEdgePoints) != 0): + CurrentPoint = DropletEdgePoints[0] + DropletEdgePoints.remove(CurrentPoint) + NeighborPoints = [] + for RowOffset in range(-3,4, 1): + for ColOffset in range(-3,4,1): + if(RowOffset != 0 or ColOffset != 0): + NeighborPoints.append((CurrentPoint[0]+RowOffset, CurrentPoint[1]+ColOffset)) + for PointRow, PointCol in NeighborPoints: + if (0 <= PointRow and PointRow < NumOfRow) and (0 <= PointCol and PointCol < NumOfCol): + if EdgeImage[PointRow][PointCol] == 255: + DropletEdgePoints.append((PointRow, PointCol)) + EdgeImage[PointRow, PointCol] = 70 + if PointRow > BottomMostRow: + BottomMostRow = PointRow + if PointCol > RightMostCol: + RightMostCol = PointCol + if PointCol < LeftMostCol: + LeftMostCol = PointCol + + CenterCol = (LeftMostCol + RightMostCol)//2 + for Row in range(BottomMostRow + 1): + if EdgeImage[Row][CenterCol] == 70: + BottomMostRow2 = Row + BottomMostRow = BottomMostRow2 + + + # Determine the transition row between the syringe pixels and the droplet pixels. + PreviousWidth = 0 + WidthIncCount = 5 + WidthHistory = [0]*(WidthIncCount+1) + TransitionRow = 0 + + for Row in range(BottomMostRow): + LeftEdgeCol = 0 + RightEdgeCol = 9999999 + for Col in range(LeftMostCol-1, RightMostCol+1, 1): + if EdgeImage[Row,Col] == 70: + LeftEdgeCol = Col + break + for Col in range(RightMostCol+1, LeftMostCol-1, -1): + if EdgeImage[Row,Col] == 70: + RightEdgeCol = Col + break + CurrentWidth = RightEdgeCol - LeftEdgeCol + if (CurrentWidth != PreviousWidth) and (0 < CurrentWidth) and (CurrentWidth < (9999999-NumOfCol)): + for i in range(WidthIncCount): + WidthHistory[i] = WidthHistory[i+1] + WidthHistory[WidthIncCount] = CurrentWidth + WidthIncreaseCount = 0 + for i in range(WidthIncCount): + if WidthHistory[i] < WidthHistory[i+1]: + WidthIncreaseCount += 1 + if WidthIncreaseCount == WidthIncCount: + TransitionRow = Row + break + if (0 < CurrentWidth) and (CurrentWidth < (9999999-NumOfCol)): + PreviousWidth = CurrentWidth + + + # Determine the width of the syringe in terms of pixel count + WidthSum = 0 + WidthCount = 0 + SyringeLevel = int(0.75 * TransitionRow) + for Row in range(TransitionRow): + RightEdgeCol = 9999999 + LeftEdgeCol = 0 + for Col in range(LeftMostCol, CenterCol + 1, 1): + if EdgeImage[Row, Col] == 70: + EdgeImage[Row, Col] = 255 + LeftEdgeCol = Col + break + for Col in range(RightMostCol, CenterCol - 1, -1): + if EdgeImage[Row, Col] == 70: + EdgeImage[Row, Col] = 255 + RightEdgeCol = Col + break + CurrentWidth = RightEdgeCol - LeftEdgeCol + if WidthCount <= SyringeLevel: + if (0 < CurrentWidth) and (CurrentWidth < (9999999-NumOfCol)): + WidthSum += CurrentWidth + WidthCount += 1 + AverageWidth = float(WidthSum)/float(WidthCount) + + + # Extract the pixel coordinates of all droplet pixels (color: 70) + PositiveEdgePoints = [] + NegativeEdgePoints = [] + for Row in range(BottomMostRow, TransitionRow-1 ,-1): + for Col in range(LeftMostCol, CenterCol + 1, 1): + if EdgeImage[Row, Col] == 70: + NegativeEdgePoints.append((Row, Col)) + break + for Col in range(RightMostCol, CenterCol, -1): + if EdgeImage[Row, Col] == 70: + PositiveEdgePoints.append((Row, Col)) + break + + + # Scale the point coordinates by using the syringe length + SyringeWidth = _SyringeWidth + ScalingFactor = SyringeWidth/AverageWidth + ScaledNegativeEdgePoints = [] + NegativeCoordinatesFile = open(_OutputFileName+"_Negative.dat", 'w') + for Row, Col in NegativeEdgePoints: + Row2 = (ScalingFactor*float(BottomMostRow-Row)) + 0.0001 + Col2 = ScalingFactor*float(Col - CenterCol) + NegativeCoordinatesFile.write(str(Col2) + " " + str(Row2) + "\n") + ScaledNegativeEdgePoints.append((Col2,Row2)) + NegativeCoordinatesFile.close() + + ScaledPositiveEdgePoints = [] + PositiveCoordinatesFile = open(_OutputFileName+"_Positive.dat", 'w') + for Row, Col in PositiveEdgePoints: + Row2 = (ScalingFactor*float(BottomMostRow-Row)) + 0.0001 + Col2 = ScalingFactor*float(Col - CenterCol) + PositiveCoordinatesFile.write(str(Col2) + " " + str(Row2) + "\n") + ScaledPositiveEdgePoints.append((Col2, Row2)) + PositiveCoordinatesFile.close() + + + # Draw the detected edge points of the droplet. + Image2 = copy.deepcopy(Image) + for Row in range(0, NumOfRow): + for Col in range(0, NumOfCol): + Image[(Row,Col)] = 0 + EdgeImage[(Row,Col)] = 0 + for Point in PositiveEdgePoints: + Image[Point] = (0, 0, 255) + EdgeImage[Point] = 255 + for Point in NegativeEdgePoints: + Image[Point] = (0, 0, 255) + EdgeImage[Point] = 255 + + # Draw the center column of the drop + cv2.line(Image, (CenterCol, TransitionRow), (CenterCol, BottomMostRow), (0, 255, 0), 1) # image, start, end, color(BGR), thickness + + MaxWidthRatio = 0.0 + for Row in range(BottomMostRow, TransitionRow-1 ,-1): + RightCol = 99999 + LeftCol = -1 + for Col in range(LeftMostCol, CenterCol + 1, 1): + if EdgeImage[Row, Col] == 255: + LeftCol = Col + break + for Col in range(RightMostCol, CenterCol, -1): + if EdgeImage[Row, Col] == 255: + RightCol = Col + break + WidthRatio = float(RightCol - LeftCol)/float(RightMostCol-LeftMostCol) + if WidthRatio < 1.0 and WidthRatio > MaxWidthRatio: + MaxWidthRatio = WidthRatio + + RowWithMaxWidth = [] + for Row in range(BottomMostRow, TransitionRow-1 ,-1): + RightCol = 99999 + LeftCol = -1 + for Col in range(LeftMostCol, CenterCol + 1, 1): + if EdgeImage[Row, Col] == 255: + LeftCol = Col + break + for Col in range(RightMostCol, CenterCol, -1): + if EdgeImage[Row, Col] == 255: + RightCol = Col + break + WidthRatio = float(RightCol - LeftCol)/float(RightMostCol-LeftMostCol) + if WidthRatio == MaxWidthRatio: + RowWithMaxWidth.append((Row, LeftCol, RightCol)) + + # Draw the main diameter + n = len(RowWithMaxWidth) + n = int(n/2) + Row = RowWithMaxWidth[n][0] + LeftCol = RowWithMaxWidth[n][1] + RightCol = RowWithMaxWidth[n][2] + cv2.line(Image, (LeftCol, Row), (RightCol, Row), (0, 255, 0), 1) + + # Measure and draw the ten diameters + Diameters = [] + MaxWidth = RightMostCol - LeftMostCol + for i in range(1, 13, 1): + for Col in range(LeftMostCol, CenterCol + 1, 1): + if EdgeImage[BottomMostRow-int(i*MaxWidth/10.0), Col] == 255: + LeftCol = Col + break + for Col in range(RightMostCol, CenterCol, -1): + if EdgeImage[BottomMostRow-int(i*MaxWidth/10.0), Col] == 255: + RightCol = Col + break + Diameters.append(RightCol-LeftCol) + cv2.line(Image, (LeftCol, BottomMostRow-int(i*MaxWidth/10.0)), (RightCol, BottomMostRow-int(i*MaxWidth/10.0)), (255, 0, 0), 1) + + #DELETE THIS + #cv2.imwrite("1.jpg", Image) + #sys.exit() + + + # Scale the 12 drop diameters + for i in range(12): + Diameters[i] *= ScalingFactor + + ScaledMaxWidth = ScalingFactor*float(MaxWidth) + ScaledUpWidth = ScalingFactor*float(RightCol - LeftCol) + return ScaledPositiveEdgePoints, ScaledNegativeEdgePoints, ScaledMaxWidth, Diameters + + + +# ========== RK4 for BASHFORTH-ADAMS FORMULATION ========== + +def BashforthAdamsRK4(_afRK4Table, _fApexCurvature, _fShapeFactor, _fSign): + A = _fApexCurvature + A2 = 2.0/A + B = _fShapeFactor + KX0 = 1.0 + KT0 = 1.0 + Inv6 = 1.0/6.0 + fTotalDiffSquared = 0.0 + n = 1 + while n < len(_afRK4Table): + x = _afRK4Table[n-1][3] # predicted x + y = _afRK4Table[n-1][1] # actual y + t = _afRK4Table[n-1][4] # predicted angle + dy = _afRK4Table[n][2] # actual y diff + AG = 1.0 #AG: Added + ydy = y + 0.5 * dy + KX1 = dy * (1.0 / math.tan(t)) + KT1 = dy * (1.0 / math.sin(t)) * (A2 - (B * y) - (math.sin(t) / (AG * x))) # A -> AG + KX2 = dy * (1.0 / math.tan(t + 0.5 * KT1)) + STKT1 = math.sin(t + 0.5 * KT1) + KT2 = dy * (1.0 / STKT1) * (A2 - (B * ydy) - (STKT1 / (AG * (x + 0.5 * KX1)))) # A -> AG + KX3 = dy * (1.0 / math.tan(t + 0.5 * KT2)) + STKT2 = math.sin(t + 0.5 * KT2) + KT3 = dy * (1.0 / STKT2) * (A2 - (B * ydy) - (STKT2 / (AG * (x + 0.5 * KX2)))) # A -> AG + KX4 = dy * (1.0 / math.tan(t + KT3)) + STKT3 = math.sin(t + KT3) + KT4 = dy * (1.0 / STKT3) * (A2 - (B * (y + dy)) - (STKT3 / (AG * (x + KX3)))) # A -> AG + KX = (KX1 + 2.0 * (KX2 + KX3) + KX4) * Inv6 + KT = (KT1 + 2.0 * (KT2 + KT3) + KT4) * Inv6 + _afRK4Table[n][3] = x + KX + _afRK4Table[n][4] = t + KT + fTotalDiffSquared += (_afRK4Table[n][0] - (_fSign * _afRK4Table[n][3])) ** 2 + n += 1 + return fTotalDiffSquared + +def MakeRK4table(_edgePoints): + InitialValues = [0, 0, 0, 0, 0, 0, 0, 0, 0, 0 , 0, 0, 0] + + rk4table = [] + rk4table.append(InitialValues) + yPrev = 0.0 + for x, y in _edgePoints: + EntryValues = [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0] + EntryValues[0] = x + EntryValues[1] = y + EntryValues[2] = y - yPrev + yPrev = y + rk4table.append(EntryValues) + + return rk4table + """ + RK4NegTable = [] + RK4NegTable.append(InitialValues) + yPrev = 0.0 + for x, y in _NegEdgePoints: + EntryValues = [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0] + EntryValues[0] = x + EntryValues[1] = y + EntryValues[2] = y - yPrev + yPrev = y + RK4NegTable.append(EntryValues) + + return RK4NegTable + """ + +# ========== FINDING OPTIMAL PARAMETERS for BASHFORTH-ADAMS SOLUTION ========== +def FindOptimalParameters(_PosEdgePoints, _NegEdgePoints, _ARange_radCurv, _BRange_shapeFact, _XRange, _TRange_angle, _Params, _Ticks): + + RK4PosTable = MakeRK4table(_PosEdgePoints) + RK4NegTable = MakeRK4table(_NegEdgePoints) + + MinimumRMSD = 999999.0 + BestParameters = [0, 0, 0, 0] + for A_radCurv in _ARange_radCurv: + for B_shapeFact in _BRange_shapeFact: + for X in _XRange: + for T_angle in _TRange_angle: + # x, y, dy, x*, t*, KX1, KT1, KX2, KT2, KX3, KT3, KX4, KT4 + fTotalDiffSquared = 0.0 + RK4PosTable[0][3] = X + RK4PosTable[0][4] = T_angle + RK4NegTable[0][3] = X + RK4NegTable[0][4] = T_angle + fTotalDiffSquared += BashforthAdamsRK4(RK4PosTable, A_radCurv, B_shapeFact, 1.0) + fTotalDiffSquared += BashforthAdamsRK4(RK4NegTable, A_radCurv, B_shapeFact, -1.0) + + CurrentRMSD = (fTotalDiffSquared/(len(_PosEdgePoints) + len(_NegEdgePoints)))**0.5 + if CurrentRMSD < MinimumRMSD: + MinimumRMSD = CurrentRMSD + BestParameters = [A_radCurv, B_shapeFact, X, T_angle] + _Ticks.put("#") + + A_radCurv = BestParameters[0] + B_shapeFact = BestParameters[1] + X = BestParameters[2] + T_angle = BestParameters[3] + ST = (9.8 * A_radCurv) / B_shapeFact + _Params.put((MinimumRMSD, A_radCurv, B_shapeFact, X, T_angle)) + + +# ========== FINDING OPTIMAL PARAMETERS USING MINIMZATION ========== +def OptimizeParameters(_PosEdgePoints, _NegEdgePoints, _AGuess_radCurv, _BGuess_shapeFact, _XGuess, _TGuess_angle): + + RK4PosTable = MakeRK4table(_PosEdgePoints) + RK4NegTable = MakeRK4table(_NegEdgePoints) + """ + ARG_RAD = 0 + ARG_SHAPE = 1 + ARG_X = 2 + ARG_ANGLE = 3 + """ + + def msd(args): + A_radCurv = args[0] + B_shapeFact = args[1] + X = args[2] + T_angle = args[3] + + RK4PosTable[0][3] = X + RK4PosTable[0][4] = T_angle + RK4NegTable[0][3] = X + RK4NegTable[0][4] = T_angle + fTotalDiffSquared = BashforthAdamsRK4(RK4PosTable, A_radCurv, B_shapeFact, 1.0) + BashforthAdamsRK4(RK4NegTable, A_radCurv, B_shapeFact, -1.0) + + return (fTotalDiffSquared/(len(_PosEdgePoints) + len(_NegEdgePoints)))**0.5 + + + optResult = scipy.optimize.minimize(msd, [_AGuess_radCurv, _BGuess_shapeFact, _XGuess, _TGuess_angle]) + if optResult.success: + print("minimization success") + #print(optResult.x) + #print(optResult.fun) + return [float(optResult.fun)] + optResult.x.tolist() + else: + return [9999999,0,0,0,0] + + +# ========== Juza's equation for 1/H ========== + +def HJuza(_S_diamRatio, _K): + # LowerS, UpperS, Ke, Kex, Kep, K3, K2, K1, K0 + + # 0.82 - 0.999 + Coefficient8 = [[ 0.82, 0.83, 1791452.8, 533.85850, 170.55837, -81676.15812, 202187.89690, -166836.27680, 45888.22973], + [ 0.83, 0.84, 189.07578, 272.52455, 72.730022, -19963.60152, 50019.79936, -41775.38186, 11629.85610], + [ 0.84, 0.85, 7.1549709, 165.45945, 35.285687, -6944.66848, 17609.43832, -14883.84379, 4193.34232], + [ 0.85, 0.87, 1.1496269, 95.066407, 12.600013, -2158.91585, 5571.60176, -4792.82331, 1374.26272], + [ 0.87, 0.895, 0.47873040, 50.357240, 0.089787883, -567.76534, 1503.51828, -1327.11835, 390.45562], + [ 0.895, 0.93, 0.35255000, 25.498893, -5.4176608, -165.99710, 454.58851, -414.93772, 126.23908], + [ 0.93, 0.97, 0.32002037, 6.8474560, -8.0901766, -102.84970, 293.25377, -278.69176, 88.27639], + [ 0.97, 0.99, 0.30845061, -32.343947, -10.455428, -475.69091, 1398.86173, -1371.17931, 448.00538], + [ 0.99, 0.999, 0.30110729, -333.50440, -15.711260, -11334.69334, 33822.93507, -33642.61426, 11154.37157]] + + # 0.635 - 0.997 + Coefficient9 = [[ 0.635, 0.645, 58249804, 119.64583, 89.483167, -28072.11478, 53918.28040, -34519.94034, 7366.77050], + [ 0.645, 0.66, 5100.2910, 70.920100, 46.811007, -9051.10879, 17726.62108, -15572.21470, 2518.10214], + [ 0.66, 0.685, 19.518159, 38.509198, 19.951285, -2338.00251, 4719.64936, -3175.58038, 712.17229], + [ 0.685, 0.72, 1.3823760, 20.055606, 5.9746092, -522.67397, 1102.26575, -774.75596, 181.49582], + [ 0.72, 0.77, 0.49074291, 10.484929, -0.31885445, -111.99730, 250.54286, -186.78287, 46.40581], + [ 0.77, 0.84, 0.34607161, 5.4063548, -2.9788620, -24.21100, 58.53312, -47.15244, 12.65668], + [ 0.84, 0.93, 0.31342828, 2.1140055, -4.1151543, -9.16969, 24.37544, -21.58758, 6.36954], + [ 0.93, 0.98, 0.30397966, -3.6334200, -4.9395699, -29.80626, 85.43510, -81.61766, 25.98669], + [ 0.98, 0.997, 0.30007321, -34.095653, -6.1531194, -434.11439, 1287.65238, -1273.10762, 419.56923]] + + # 0.17 - 0.983 + Coefficient10 = [[ 0.17, 0.30, 0.31081678, -0.086278355, -2.7023254, -13.95071, 10.30398, -2.49619, 0.19805], + [ 0.30, 0.45, 0.30636442, -0.094871613, -2.7246428, 0.38766, -0.44128, 0.16613, -0.02068], + [ 0.45, 0.68, 0.31156188, -0.063734909, -2.6789763, 0.41537, -0.71168, 0.40321, -0.07551], + [ 0.68, 0.90, 0.31195754, -0.092720991, -2.6863859, -0.44813, 1.06637, -0.84260, 0.22106], + [ 0.90, 0.983, 0.30712046, -1.5619311, -2.9809169, -5.74538, 16.23781, -15.29128, 4.79808]] + + # 0.06 - 0.953 + Coefficient11 = [[ 0.06, 0.08, 1.6030433, 0.043380701, -0.15208454, -1341.32950, 282.27311, -19.71010, 0.45664], + [ 0.08, 0.13, 0.85491410, -0.056761738, -0.65372414, -160.97221, 50.59168, -5.23625, 0.17842], + [ 0.13, 0.20, 0.57401485, -0.15265815, -1.0443069, -42.79943, 21.22139, -3.47357, 0.18765], + [ 0.20, 0.30, 0.41622542, -0.27871823, -1.4464670, -14.08741, 10.58863, -2.63159, 0.21621], + [ 0.30, 0.44, 0.33587066, -0.42674255, -1.8022482, -3.22540, 3.58964, -1.32194, 0.16106], + [ 0.44, 0.75, 0.31504207, -0.51197565, -1.9499037, 0.07609, -0.13681, 0.08090, -0.01572], + [ 0.75, 0.953, 0.31658938, -0.49808609, -1.9290368, -0.24018, 0.61450, -0.52258, 0.14771]] + + # 0.10 - 0.902 + Coefficient12 = [[ 0.10, 0.12, 103.56308, 1.2019976, 4.4869062, -6625.95443, 2292.46751, -264.28560, 10.15218], + [ 0.12, 0.14, 2.1585254, 0.34129627, 0.83653215, -742.15388, 289.52312, -37.59979, 1.62555], + [ 0.14, 0.19, 0.65654048, 0.031877859, -0.37693028, -66.92134, 33.18430, -5.45788, 0.29773], + [ 0.19, 0.27, 0.45076908, -0.10378254, -0.82836350, -10.66154, 7.36263, -1.68453, 0.12768], + [ 0.27, 0.40, 0.37085833, -0.21965481, -1.1287016, -2.87493, 2.88810, -0.95966, 0.10546], + [ 0.40, 0.55, 0.33605594, -0.33317467, -1.3397591, -0.70394, 1.00402, -0.47492, 0.07450], + [ 0.55, 0.902, 0.32790337, -0.39597877, -1.4180887, 0.00490, -0.01070, 0.00769, -0.00182]] + + if ( _K == 0.8): + Coeff = Coefficient8 + if ( _K == 0.9): + Coeff = Coefficient9 + if ( _K == 1.0): + Coeff = Coefficient10 + if ( _K == 1.1): + Coeff = Coefficient11 + if ( _K == 1.2): + Coeff = Coefficient12 + + # Juza 1996/1997 Equation: 8 + for Entry_k in Coeff: + if (Entry_k[0] <= _S_diamRatio) and (_S_diamRatio <= Entry_k[1]): + HInv = Entry_k[2] * _S_diamRatio**(Entry_k[3] * math.log(_S_diamRatio) + Entry_k[4]) + Entry_k[5] * (_S_diamRatio**3.0) + Entry_k[6] * (_S_diamRatio **2.0) + Entry_k[7] * _S_diamRatio + Entry_k[8] + return HInv + break + + return 0 + + + +# ========== USAGE PRINTER ========== + +def PrintCorrectUsage(): + print("\nDESCRIPTION:\n") + print("This program will take an image or a set of files (coordinate points) representing") + print("a droplet and will process the data to calculate the surface tension of the droplet.") + print("When an image file is set, the program will extract the coordinate of the points") + print("representing the edge of the droplet. When the files are set, it is assume that they") + print("will already contain the points representing the edge of the droplet. RK4 will be") + print("use to solve the Bashforth-Adams set of equations to get the parameters for surface") + print("tension calculation. parameters: A (apex radius curvature) and B (shape factor)") + print("\nUSAGE:\n") + print("Syntax 1: ./img.py -option1 value1 -option2 value2") + print("Syntax 2: python img.py -option1 value1 -option2 value2\n") + print("Options:") + print("-help : Prints this help message. ") + print("-out : Reference name for all the output files (coordinate files, RK4 dump, etc.).") + print(" : This should always be set. ") + print("-img : Image file to be processed. ") + print(" : When this is set, the file name inputs (-filep, -filen) should not be set.") + print("-filep : File containing the positive(x) point coordinates to be processed.") + print(" : If this is set, -filen should also be set and -img should not be set.") + print("-filen : File containing the negative(x) point coordinates to be processed.") + print(" : If this is set, -filep should also be set and -img should not be set.") + print("-proc : The number of processor to use.") + print(" : If not specified the default value is 2.") + print("-syrg : This is the syringe width in millimeters.") + print(" : If not specified the default value is 0.805 mm.") + print("-a : RK4 Parameter: minimum value for alpha (radius of curvature at apex).") + print(" : If not specified the default value is 0.0001") + print("-A : RK4 Parameter: maximum value for alpha (radius of curvature at apex).") + print(" : If not specified the default value is 2.0001") + print("-b : RK4 Parameter: minimum value for beta (shape factor).") + print(" : If not specified the default value is 0.0001") + print("-B : RK4 Parameter: minimum value for beta (shape factor).") + print(" : If not specified the default value is 0.5001.") + print("-x : RK4 Parameter: minimum value for x (initial guessed x coordinate).") + print(" : If not specified the default value is 0.0001.") + print("-X : RK4 Parameter: maximum value for x (initial guessed x coordinate).") + print(" : If not specified the default value is 0.5001.") + print("-t : RK4 Parameter: minimum value for t (initial guessed angle).") + print(" : If not specified the default value is 0.0001.") + print("-T : RK4 Parameter: maximum value for t (initial guessed angle).") + print(" : If not specified the default value is 0.5001.") + print("-ai : RK4 Parameter: alpha interval - the number of values to consider between min alpha and max alpha.") + print(" : If not specified the default value is 20.") + print("-bi : RK4 Parameter: beta interval - the number of values to consider between min beta and max beta.") + print(" : If not specified the default value is 20.") + print("-xi : RK4 Parameter: x interval - the number of values to consider between min x and max x.") + print(" : If not specified the default value is 5.") + print("-ti : RK4 Parameter: angle interval - the number of values to consider between min angle and max angle.") + print(" : If not specified the default value is 5.") + print("-solver : 1 for divide and conquer scanning.") + print(" : If not specified, optimization will be used to solve the parameters.") + print("\nSAMPLE USAGE:\n") + print("1. Calculate surface tension of a drop using default values.") + print("./img.py -img water-drop.jpg -out water-drop") + print("python img.py -img water-drop.jpg -out water-drop\n") + print("2. Uses 4 processors for calculation.") + print("./img.py -img water-drop.jpg -out water-drop -proc 4\n") + print("3. Setting interval sizes for A and B parameters (search for 30 possible values for A and B).") + print("./img.py -img water-drop.jpg -out water-drop -ai 30 -bi 30\n") + print("4. Setting min (0.00001) and max (1.5) values for parameter.") + print("./img.py -img water-drop.jpg -out water-drop -a 0.00001 -A 1.5\n") + + + +# ========== MAIN ========== + +if __name__ == '__main__': + + ImageFilePath = "" + OutputFileName = "" + PosPointsFileName = "" + NegPointsFileName = "" + ProcessCount = 2 + SyringeWidth = 0.70 * 1.15 + MinA_radCurv = 0.0001 + MaxA_radCurv = 2.0001 + MinB_shapeFact = 0.0001 + MaxB_shapeFact = 1.0001 + MinX = 0.0001 + MaxX = 0.5001 + MinT_angle = 0.0001 + MaxT_angle = 0.5001 + AStep = 20 + BStep = 20 + XStep = 5 + TStep = 5 + + if (len(sys.argv) <= 1): + PrintCorrectUsage() + sys.exit() + + ImageFlag = 0 + PosFileFlag = 0 + NegFileFlag = 0 + OutFlag = 0 + Counter = 1 + SolveMethod = 0 + + while (Counter <= len(sys.argv)-1): + + CurrentOption = sys.argv[Counter] + try: + CurrentValue = sys.argv[Counter+1] + except: + PrintCorrectUsage() + sys.exit() + + if CurrentOption == "-help": + PrintCorrectUsage() + sys.exit() + elif CurrentOption == "-img": + ImageFilePath = CurrentValue + ImageFlag = 1 + elif CurrentOption == "-filep": + PosPointsFileName = CurrentValue + PosFileFlag = 1 + elif CurrentOption == "-filen": + NegPointsFileName = CurrentValue + NegFileFlag = 1 + elif CurrentOption == "-out": + OutputFileName = CurrentValue + OutFlag = 1 + elif CurrentOption == "-proc": + ProcessCount = int(CurrentValue) + elif CurrentOption == "-syrg": + SyringeWidth = float(CurrentValue) + elif CurrentOption == "-a": + MinA_radCurv = float(CurrentValue) + elif CurrentOption == "-A": + MaxA_radCurv = float(CurrentValue) + elif CurrentOption == "-b": + MinB_shapeFact = float(CurrentValue) + elif CurrentOption == "-B": + MaxB_shapeFact = float(CurrentValue) + elif CurrentOption == "-x": + MinX = float(CurrentValue) + elif CurrentOption == "-X": + MaxX = float(CurrentValue) + elif CurrentOption == "-t": + MinT_angle = float(CurrentValue) + elif CurrentOption == "-T": + MaxT_angle = float(CurrentValue) + elif CurrentOption == "-ai": + AStep = int(CurrentValue) + elif CurrentOption == "-bi": + BStep = int(CurrentValue) + elif CurrentOption == "-xi": + XStep = int(CurrentValue) + elif CurrentOption == "-ti": + TStep = int(CurrentValue) + elif CurrentOption == "-solver": + SolveMethod = int(CurrentValue) + else: + PrintCorrectUsage() + sys.exit() + Counter = Counter + 2 + + if (OutFlag == 0 and ImageFlag == 1): + print("Output file name not set. (-out)") + PrintCorrectUsage() + sys.exit() + if (ImageFlag == 1 and ((PosFileFlag == 1) or (NegFileFlag == 1))): + print("Image(-img) can not be set with files option (-filep, -filen).") + print("Either set image only or files only.") + PrintCorrectUsage() + sys.exit() + if (ImageFlag == 0 and ((PosFileFlag == 1 and NegFileFlag == 0) or (PosFileFlag == 0 and NegFileFlag == 1))): + print("File inputs (-filep, -filen) should both be set if image input is not set.") + PrintCorrectUsage() + sys.exit() + + ARange_radCurv = [] + a = MinA_radCurv + ahop = float(MaxA_radCurv - MinA_radCurv)/float(AStep-1) + while (a <= MaxA_radCurv): + ARange_radCurv.append(a) + a += ahop + if ahop == 0.0: + break + if(np.abs(MaxA_radCurv - ARange_radCurv[len(ARange_radCurv)-1] - ahop) < 0.000001 and ahop != 0.0): + ARange_radCurv.append(MaxA_radCurv) + + BRange_shapeFact = [] + b = MinB_shapeFact + bhop = float(MaxB_shapeFact - MinB_shapeFact)/float(BStep-1) + while (b <= MaxB_shapeFact): + BRange_shapeFact.append(b) + b += bhop + if bhop == 0: + break + if(np.abs(MaxB_shapeFact - BRange_shapeFact[len(BRange_shapeFact)-1] - bhop) < 0.000001 and bhop != 0.0): + BRange_shapeFact.append(MaxB_shapeFact) + + XRange = [] + x = MinX + xhop = float(MaxX - MinX)/float(XStep-1) + while (x <= MaxX): + XRange.append(x) + x += xhop + if xhop == 0: + break + if(np.abs(MaxX - XRange[len(XRange)-1] - xhop) < 0.000001 and xhop != 0.0): + XRange.append(MaxX) + + TRange_angle = [] + t = MinT_angle + thop = float(MaxT_angle - MinT_angle)/float(TStep-1) + while (t <= MaxT_angle): + TRange_angle.append(t) + t += thop + if thop == 0.0: + break + if(np.abs(MaxT_angle - TRange_angle[len(TRange_angle)-1] - thop) < 0.000001 and thop != 0.0): + TRange_angle.append(MaxT_angle) + + + print("\n") + if (ImageFlag == 1): + print("Image File : " + ImageFilePath) + print("Output Name : " + OutputFileName) + print("Negative Cooridnate File : " + OutputFileName + "_Negative.dat") + print("Positive Cooridnate File : " + OutputFileName + "_Positive.dat") + print("Syringe Width (millimeter): " + str(SyringeWidth)) + else: + print("Negative Cooridnate File : " + NegPointsFileName) + print("Positive Cooridnate File : " + PosPointsFileName) + print("Radius of Curvature Range : " + str(MinA_radCurv) + "-" + str(MaxA_radCurv) + " (interval: "+str(len(ARange_radCurv))+")") + print("Shape Factor Range : " + str(MinB_shapeFact) + "-" + str(MaxB_shapeFact) + " (interval: "+str(len(BRange_shapeFact))+")") + print("Initial X Coordinate Range: " + str(MinX) + "-" + str(MaxX) + " (interval: "+str(len(XRange))+")") + print("Initial Angle Range : " + str(MinT_angle) + "-" + str(MaxT_angle) + " (interval: "+str(len(TRange_angle))+")") + print("Parameter Combinations : " + str(len(ARange_radCurv) * len(BRange_shapeFact) * len(XRange) * len(TRange_angle))) + print("Number of Process to Use : " + str(ProcessCount)) + + if (ImageFlag == 1): + Pos_edgePts, Neg_edgePts, De_scaledMaxWidth, Dms_diamters = ExtractEdgeCoordinate(ImageFilePath, OutputFileName, SyringeWidth) + HInv1 = HJuza(Dms_diamters[7]/De_scaledMaxWidth, 0.8) + HInv2 = HJuza(Dms_diamters[8]/De_scaledMaxWidth, 0.9) + HInv3 = HJuza(Dms_diamters[9]/De_scaledMaxWidth, 1.0) + HInv4 = HJuza(Dms_diamters[10]/De_scaledMaxWidth, 1.1) + HInv5 = HJuza(Dms_diamters[11]/De_scaledMaxWidth, 1.2) + HInvM = (0.20)*(HInv1 + HInv2 + HInv3 + HInv4 + HInv5) + HInvSD = math.sqrt((0.20)*( (HInv1-HInvM)**2.0 + (HInv2-HInvM)**2.0 + (HInv3-HInvM)**2.0 + (HInv4-HInvM)**2.0 + (HInv5-HInvM)**2.0)) + print("Juza 5-Plane (1/H) Std.Dev: " + str(HInvSD)) + else: + Pos_edgePts = [] + PositivePointFile = open(PosPointsFileName,'r') + for Line in PositivePointFile: + Coordinate = Line.split() + #print(Coordinate) + Pos_edgePts.append((float(Coordinate[0]), float(Coordinate[1]))) + PositivePointFile.close() + Neg_edgePts = [] + NegativePointFile = open(NegPointsFileName,'r') + for Line in NegativePointFile: + Coordinate = Line.split() + #print(Coordinate) + Neg_edgePts.append((float(Coordinate[0]), float(Coordinate[1]))) + NegativePointFile.close() + + def BruteForce(): + PosList = [] + NegList = [] + for ProcessNo in range(ProcessCount): + PosList.append(copy.deepcopy(Pos_edgePts)) + NegList.append(copy.deepcopy(Neg_edgePts)) + #print(ProcessNo) + + ARanges = [] + for ProcessNo in range(ProcessCount): + ARanges.append([]) + + i = 0 + for A in ARange_radCurv: + ARanges[i%ProcessCount].append(A) + i += 1 + + BestParams = multiprocessing.Queue() + TicksQueue = multiprocessing.Queue() + Processes = [] + + for ProcessNo in range(ProcessCount): + Processes.append(multiprocessing.Process(name='Search ' + str(ProcessNo), target=FindOptimalParameters, args=(PosList[ProcessNo], NegList[ProcessNo], ARanges[ProcessNo], BRange_shapeFact, XRange, TRange_angle, BestParams, TicksQueue))) + + for aProcess in Processes: + aProcess.start() + + sys.stdout.write("\rParameter Search Progress : 0%") + sys.stdout.flush() + Count = 0 + while(Count < len(ARange_radCurv)): + t = TicksQueue.get() + Count += 1 + sys.stdout.write("\rParameter Search Progress : " + str(int(100*float(Count)/float(len(ARange_radCurv)))) + "%") + sys.stdout.flush() + + for aProcess in Processes: + aProcess.join() + + BestOutput = (999999999, 1, 1, 1 ,1) + while not BestParams.empty(): + Output = BestParams.get() + if Output[0] < BestOutput[0]: + BestOutput = Output + return BestOutput + + def Optimize(): + return OptimizeParameters(Pos_edgePts, Neg_edgePts, 0.95, 0.3, 0.1, 0.2) + + starttime = time.time() + + BestOutput = Optimize() if 1 != SolveMethod else BruteForce() + + print("best output: " + str(BestOutput)) + + print("\nLowest RMSD : " + str(BestOutput[0])) + print("Computed Surface Tension : " + str(9.8*BestOutput[1]/BestOutput[2])) + print("Best Radius of Curvature : " + str(BestOutput[1])) + print("Best Shape Factor : " + str(BestOutput[2])) + print("Best Initial X Coordinate : " + str(BestOutput[3])) + print("Best Initial Angle : " + str(BestOutput[4])) + print("Lapsed Time in Seconds : " + str(time.time() -starttime)) diff --git a/static/css/style.css b/static/css/style.css new file mode 100644 index 0000000..79d768b --- /dev/null +++ b/static/css/style.css @@ -0,0 +1,33 @@ +body { + background-color: #7b1113; + font-family: 'Montserrat', sans-serif !important; + font-size: 15px; + color: #E0E0E0FF; +} + +.header { + font-size: 50px; +} + +.maindiv { + background-color: #E0E0E0; + color: #000000; + border-radius: 30px; +} + +.imagediv { + background-color: #9E9E9E; + border-radius: 30px; + display: flex; + justify-content: center; + align-items: center; +} + +.labels { + font-weight: bold; +} + +.btn { + border-radius: 30px; + width: 100%; +} \ No newline at end of file diff --git a/static/images/ISC Logo White.png b/static/images/ISC Logo White.png new file mode 100644 index 0000000..e0eff3b Binary files /dev/null and b/static/images/ISC Logo White.png differ diff --git a/templates/base.html b/templates/base.html new file mode 100644 index 0000000..b174bd9 --- /dev/null +++ b/templates/base.html @@ -0,0 +1,33 @@ + + + + + + + + + + + + + + UP ISC | Pendant Drop Project + + +
+ {% block content %} {% endblock %} +
+
+ © 2025 UP Intelligent Systems Center +
+
+
+ + + + + + + + + diff --git a/templates/index.html b/templates/index.html new file mode 100644 index 0000000..e709258 --- /dev/null +++ b/templates/index.html @@ -0,0 +1,71 @@ +{% extends 'base.html' %} + +{% block content %} +
+
+ + {% block title %} | PENDANT DROP TENSIOMETRY {% endblock %} +
+
+
+ +
+ {% if file_url %} +
+ +
+ {% else %} +
+ No image uploaded +
+ {% endif %} +
+
+ {{ form.hidden_tag() }} + {% for error in form.photo.errors %} + {{ error }} + {% endfor %} + +
+
+ {{ form.photo(class="form-control") }} +
+
+ {{ form.submit(class="btn btn-success") }} +
+
+
+ +
+ +
+

Results:

+
+
+ Negative Cooridnate File:
+ Positive Cooridnate File:
+ Syringe Width (millimeter):
+ Radius of Curvature Range:
+ Shape Factor Range :
+ Initial X Coordinate Range:
+ Initial Angle Range:
+ Parameter Combinations:
+ Number of Process to Use:
+ Juza 5-Plane (1/H) Std.Dev:
+ + best output:
+ Lowest RMSD:
+ Computed Surface Tension:
+ Best Radius of Curvature:
+ Best Shape Factor:
+ Best Initial X Coordinate:
+ Best Initial Angle:
+ Lapsed Time in Seconds:
+
+
+ asdsadsa +
+
+
+
+{% endblock %}